mirror of
https://github.com/rdkit/rdkit.git
synced 2026-06-03 21:44:30 +08:00
997 lines
34 KiB
C++
997 lines
34 KiB
C++
//
|
|
// Copyright (C) 2004-2024 Paolo Tosco and other RDKit contributors
|
|
//
|
|
// @@ All Rights Reserved @@
|
|
// This file is part of the RDKit.
|
|
// The contents are covered by the terms of the BSD license
|
|
// which is included in the file license.txt, found at the root
|
|
// of the RDKit source tree.
|
|
//
|
|
//
|
|
#include <RDGeneral/export.h>
|
|
#ifndef RD_MMFFPARAMS_H
|
|
#define RD_MMFFPARAMS_H
|
|
|
|
#include <memory>
|
|
#include <RDGeneral/Invariant.h>
|
|
#include <cmath>
|
|
#include <string>
|
|
#include <vector>
|
|
#include <algorithm>
|
|
#include <map>
|
|
#include <cstdint>
|
|
|
|
#ifndef M_PI
|
|
#define M_PI 3.14159265358979323846
|
|
#endif
|
|
|
|
// binary searches are slightly faster than std::map;
|
|
// however when I moved to binary searches I had already
|
|
// written the code for std::map, so the two methods
|
|
// can be toggled defining RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
|
|
// #define RDKIT_MMFF_PARAMS_USE_STD_MAP 1
|
|
|
|
namespace ForceFields {
|
|
namespace MMFF {
|
|
|
|
constexpr double DEG2RAD = M_PI / 180.0;
|
|
constexpr double RAD2DEG = 180.0 / M_PI;
|
|
constexpr double MDYNE_A_TO_KCAL_MOL = 143.9325;
|
|
inline bool isDoubleZero(const double x) {
|
|
return ((x < 1.0e-10) && (x > -1.0e-10));
|
|
}
|
|
inline void clipToOne(double &x) { x = std::clamp(x, -1.0, 1.0); }
|
|
|
|
//! class to store MMFF atom type equivalence levels
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFDef {
|
|
public:
|
|
std::uint8_t eqLevel[4];
|
|
};
|
|
|
|
//! class to store MMFF Properties
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFProp {
|
|
public:
|
|
std::uint8_t atno;
|
|
std::uint8_t crd;
|
|
std::uint8_t val;
|
|
std::uint8_t pilp;
|
|
std::uint8_t mltb;
|
|
std::uint8_t arom;
|
|
std::uint8_t linh;
|
|
std::uint8_t sbmb;
|
|
};
|
|
|
|
//! class to store MMFF Partial Bond Charge Increments
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFPBCI {
|
|
public:
|
|
double pbci;
|
|
double fcadj;
|
|
};
|
|
|
|
//! class to store MMFF bond-charge-increment parameters used to
|
|
//! construct MMFF partial atomic charges
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFChg {
|
|
public:
|
|
double bci;
|
|
};
|
|
|
|
//! class to store MMFF parameters for bond stretching
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFBond {
|
|
public:
|
|
double kb;
|
|
double r0;
|
|
};
|
|
|
|
//! class to store parameters for Herschbach-Laurie's version
|
|
//! of Badger's rule
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFHerschbachLaurie {
|
|
public:
|
|
double a_ij;
|
|
double d_ij;
|
|
double dp_ij;
|
|
};
|
|
|
|
//! class to store covalent radius and Pauling electronegativity
|
|
//! values for MMFF bond stretching empirical rule
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFCovRadPauEle {
|
|
public:
|
|
double r0;
|
|
double chi;
|
|
};
|
|
|
|
//! class to store MMFF parameters for angle bending
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFAngle {
|
|
public:
|
|
double ka;
|
|
double theta0;
|
|
};
|
|
|
|
//! class to store MMFF parameters for stretch-bending
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFStbn {
|
|
public:
|
|
double kbaIJK;
|
|
double kbaKJI;
|
|
};
|
|
|
|
//! class to store MMFF parameters for out-of-plane bending
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFOop {
|
|
public:
|
|
double koop;
|
|
};
|
|
|
|
//! class to store MMFF parameters for torsions
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFTor {
|
|
public:
|
|
double V1;
|
|
double V2;
|
|
double V3;
|
|
};
|
|
|
|
//! class to store MMFF parameters for non-bonded Van der Waals
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFVdW {
|
|
public:
|
|
double alpha_i;
|
|
double N_i;
|
|
double A_i;
|
|
double G_i;
|
|
double R_star;
|
|
std::uint8_t DA;
|
|
};
|
|
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFVdWRijstarEps {
|
|
public:
|
|
double R_ij_starUnscaled;
|
|
double epsilonUnscaled;
|
|
double R_ij_star;
|
|
double epsilon;
|
|
};
|
|
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFAromCollection {
|
|
public:
|
|
//! Looks up the parameters for a particular key and returns them.
|
|
/*!
|
|
\return a pointer to the MMFFArom object, NULL on failure.
|
|
*/
|
|
bool isMMFFAromatic(const unsigned int atomType) const {
|
|
return std::find(d_params.begin(), d_params.end(), atomType) !=
|
|
d_params.end();
|
|
}
|
|
|
|
MMFFAromCollection(const std::vector<std::uint8_t> *mmffArom = nullptr);
|
|
std::vector<std::uint8_t> d_params; //!< the aromatic type vector
|
|
};
|
|
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFDefCollection {
|
|
public:
|
|
//! Looks up the parameters for a particular key and returns them.
|
|
/*!
|
|
\return a pointer to the MMFFDef object, NULL on failure.
|
|
*/
|
|
const MMFFDef *operator()(const unsigned int atomType) const {
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
const auto res = d_params.find(atomType);
|
|
return ((res != d_params.end()) ? &((*res).second) : NULL);
|
|
#else
|
|
return ((atomType && (atomType <= d_params.size()))
|
|
? &d_params[atomType - 1]
|
|
: nullptr);
|
|
#endif
|
|
}
|
|
|
|
MMFFDefCollection(std::string mmffDef = "");
|
|
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
std::map<const unsigned int, MMFFDef> d_params; //!< the parameter map
|
|
#else
|
|
std::vector<MMFFDef> d_params; //!< the parameter vector
|
|
#endif
|
|
};
|
|
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFPropCollection {
|
|
public:
|
|
//! Looks up the parameters for a particular key and returns them.
|
|
/*!
|
|
\return a pointer to the MMFFProp object, NULL on failure.
|
|
*/
|
|
const MMFFProp *operator()(const unsigned int atomType) const {
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
const auto res = d_params.find(atomType);
|
|
return ((res != d_params.end()) ? &((*res).second) : NULL);
|
|
#else
|
|
auto bounds =
|
|
std::equal_range(d_iAtomType.begin(), d_iAtomType.end(), atomType);
|
|
return ((bounds.first != bounds.second)
|
|
? &d_params[bounds.first - d_iAtomType.begin()]
|
|
: nullptr);
|
|
#endif
|
|
}
|
|
|
|
MMFFPropCollection(std::string mmffProp = "");
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
std::map<const unsigned int, MMFFProp> d_params; //!< the parameter map
|
|
#else
|
|
std::vector<MMFFProp> d_params;
|
|
std::vector<std::uint8_t> d_iAtomType; //!< the parameter vector
|
|
#endif
|
|
};
|
|
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFPBCICollection {
|
|
public:
|
|
//! Looks up the parameters for a particular key and returns them.
|
|
/*!
|
|
\return a pointer to the MMFFPBCI object, NULL on failure.
|
|
*/
|
|
const MMFFPBCI *operator()(const unsigned int atomType) const {
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
const auto res = d_params.find(atomType);
|
|
return ((res != d_params.end()) ? &((*res).second) : NULL);
|
|
#else
|
|
return ((atomType && (atomType <= d_params.size()))
|
|
? &d_params[atomType - 1]
|
|
: nullptr);
|
|
#endif
|
|
}
|
|
|
|
MMFFPBCICollection(std::string mmffPBCI = "");
|
|
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
std::map<const unsigned int, MMFFPBCI> d_params; //!< the parameter map
|
|
#else
|
|
std::vector<MMFFPBCI> d_params; //!< the parameter vector
|
|
#endif
|
|
};
|
|
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFChgCollection {
|
|
public:
|
|
//! Looks up the parameters for a particular key and returns them.
|
|
/*!
|
|
\return a pointer to the MMFFChg object, NULL on failure.
|
|
*/
|
|
const std::pair<int, const MMFFChg *> getMMFFChgParams(
|
|
const unsigned int bondType, const unsigned int iAtomType,
|
|
const unsigned int jAtomType) const {
|
|
int sign = -1;
|
|
const MMFFChg *mmffChgParams = nullptr;
|
|
unsigned int canIAtomType = iAtomType;
|
|
unsigned int canJAtomType = jAtomType;
|
|
if (iAtomType > jAtomType) {
|
|
canIAtomType = jAtomType;
|
|
canJAtomType = iAtomType;
|
|
sign = 1;
|
|
}
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
const auto res1 = d_params[bondType].find(canIAtomType);
|
|
if (res1 != d_params[bondType].end()) {
|
|
const auto res2 = ((*res1).second).find(canJAtomType);
|
|
if (res2 != ((*res1).second).end()) {
|
|
mmffChgParams = &((*res2).second);
|
|
}
|
|
}
|
|
#else
|
|
auto bounds =
|
|
std::equal_range(d_iAtomType.begin(), d_iAtomType.end(), canIAtomType);
|
|
if (bounds.first != bounds.second) {
|
|
bounds = std::equal_range(
|
|
d_jAtomType.begin() + (bounds.first - d_iAtomType.begin()),
|
|
d_jAtomType.begin() + (bounds.second - d_iAtomType.begin()),
|
|
canJAtomType);
|
|
if (bounds.first != bounds.second) {
|
|
bounds = std::equal_range(
|
|
d_bondType.begin() + (bounds.first - d_jAtomType.begin()),
|
|
d_bondType.begin() + (bounds.second - d_jAtomType.begin()),
|
|
bondType);
|
|
if (bounds.first != bounds.second) {
|
|
mmffChgParams = &d_params[bounds.first - d_bondType.begin()];
|
|
}
|
|
}
|
|
}
|
|
#endif
|
|
|
|
return std::make_pair(sign, mmffChgParams);
|
|
}
|
|
|
|
MMFFChgCollection(std::string mmffChg = "");
|
|
|
|
//!< the parameter 3D-map
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
std::map<const unsigned int,
|
|
std::map<const unsigned int, std::map<const unsigned int, MMFFChg>>>
|
|
d_params; //!< the parameter 3D-map
|
|
#else
|
|
std::vector<MMFFChg> d_params; //!< the parameter vector
|
|
std::vector<std::uint8_t> d_iAtomType; //!< atom type vector for atom i
|
|
std::vector<std::uint8_t> d_jAtomType; //!< atom type vector for atom j
|
|
std::vector<std::uint8_t> d_bondType; //!< bond type vector for bond i-j
|
|
#endif
|
|
};
|
|
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFBondCollection {
|
|
public:
|
|
//! Looks up the parameters for a particular key and returns them.
|
|
/*!
|
|
\return a pointer to the MMFFBond object, NULL on failure.
|
|
*/
|
|
const MMFFBond *operator()(const unsigned int bondType,
|
|
const unsigned int atomType,
|
|
const unsigned int nbrAtomType) const {
|
|
const MMFFBond *mmffBondParams = nullptr;
|
|
unsigned int canAtomType = atomType;
|
|
unsigned int canNbrAtomType = nbrAtomType;
|
|
if (atomType > nbrAtomType) {
|
|
canAtomType = nbrAtomType;
|
|
canNbrAtomType = atomType;
|
|
}
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
const auto res1 = d_params.find(bondType);
|
|
std::map<const unsigned int,
|
|
std::map<const unsigned int, MMFFBond>>::const_iterator res2;
|
|
std::map<const unsigned int, MMFFBond>::const_iterator res3;
|
|
if (res1 != d_params.end()) {
|
|
res2 = ((*res1).second).find(canAtomType);
|
|
if (res2 != ((*res1).second).end()) {
|
|
res3 = ((*res2).second).find(canNbrAtomType);
|
|
if (res3 != ((*res2).second).end()) {
|
|
mmffBondParams = &((*res3).second);
|
|
}
|
|
}
|
|
}
|
|
#else
|
|
auto bounds =
|
|
std::equal_range(d_iAtomType.begin(), d_iAtomType.end(), canAtomType);
|
|
if (bounds.first != bounds.second) {
|
|
bounds = std::equal_range(
|
|
d_jAtomType.begin() + (bounds.first - d_iAtomType.begin()),
|
|
d_jAtomType.begin() + (bounds.second - d_iAtomType.begin()),
|
|
canNbrAtomType);
|
|
if (bounds.first != bounds.second) {
|
|
bounds = std::equal_range(
|
|
d_bondType.begin() + (bounds.first - d_jAtomType.begin()),
|
|
d_bondType.begin() + (bounds.second - d_jAtomType.begin()),
|
|
bondType);
|
|
if (bounds.first != bounds.second) {
|
|
mmffBondParams = &d_params[bounds.first - d_bondType.begin()];
|
|
}
|
|
}
|
|
}
|
|
#endif
|
|
|
|
return mmffBondParams;
|
|
}
|
|
|
|
MMFFBondCollection(std::string mmffBond = "");
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
std::map<const unsigned int,
|
|
std::map<const unsigned int, std::map<const unsigned int, MMFFBond>>>
|
|
d_params; //!< the parameter 3D-map
|
|
#else
|
|
std::vector<MMFFBond> d_params; //!< the parameter vector
|
|
std::vector<std::uint8_t> d_iAtomType; //!< atom type vector for atom i
|
|
std::vector<std::uint8_t> d_jAtomType; //!< atom type vector for atom j
|
|
std::vector<std::uint8_t> d_bondType; //!< bond type vector for bond i-j
|
|
#endif
|
|
};
|
|
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFBndkCollection {
|
|
public:
|
|
//! Looks up the parameters for a particular key and returns them.
|
|
/*!
|
|
\return a pointer to the MMFFBndk object, NULL on failure.
|
|
*/
|
|
const MMFFBond *operator()(const int atomicNum,
|
|
const int nbrAtomicNum) const {
|
|
const MMFFBond *mmffBndkParams = nullptr;
|
|
unsigned int canAtomicNum = atomicNum;
|
|
unsigned int canNbrAtomicNum = nbrAtomicNum;
|
|
if (atomicNum > nbrAtomicNum) {
|
|
canAtomicNum = nbrAtomicNum;
|
|
canNbrAtomicNum = atomicNum;
|
|
}
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
const auto res1 = d_params.find(canAtomicNum);
|
|
std::map<const unsigned int, MMFFBond>::const_iterator res2;
|
|
if (res1 != d_params.end()) {
|
|
res2 = ((*res1).second).find(canNbrAtomicNum);
|
|
if (res2 != ((*res1).second).end()) {
|
|
mmffBndkParams = &((*res2).second);
|
|
}
|
|
}
|
|
#else
|
|
auto bounds = std::equal_range(d_iAtomicNum.begin(), d_iAtomicNum.end(),
|
|
canAtomicNum);
|
|
if (bounds.first != bounds.second) {
|
|
bounds = std::equal_range(
|
|
d_jAtomicNum.begin() + (bounds.first - d_iAtomicNum.begin()),
|
|
d_jAtomicNum.begin() + (bounds.second - d_iAtomicNum.begin()),
|
|
canNbrAtomicNum);
|
|
if (bounds.first != bounds.second) {
|
|
mmffBndkParams = &d_params[bounds.first - d_jAtomicNum.begin()];
|
|
}
|
|
}
|
|
#endif
|
|
|
|
return mmffBndkParams;
|
|
}
|
|
|
|
MMFFBndkCollection(std::string mmffBndk = "");
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
std::map<const unsigned int, std::map<const unsigned int, MMFFBond>>
|
|
d_params; //!< the parameter 2D-map
|
|
#else
|
|
std::vector<MMFFBond> d_params; //!< the parameter vector
|
|
std::vector<std::uint8_t> d_iAtomicNum; //!< atomic number vector for atom i
|
|
std::vector<std::uint8_t> d_jAtomicNum; //!< atomic number vector for atom j
|
|
#endif
|
|
};
|
|
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFHerschbachLaurieCollection {
|
|
public:
|
|
//! Looks up the parameters for a particular key and returns them.
|
|
/*!
|
|
\return a pointer to the MMFFHerschbachLaurie object, NULL on failure.
|
|
*/
|
|
const MMFFHerschbachLaurie *operator()(const int iRow, const int jRow) const {
|
|
const MMFFHerschbachLaurie *mmffHerschbachLaurieParams = nullptr;
|
|
unsigned int canIRow = iRow;
|
|
unsigned int canJRow = jRow;
|
|
if (iRow > jRow) {
|
|
canIRow = jRow;
|
|
canJRow = iRow;
|
|
}
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
const auto res1 = d_params.find(canIRow);
|
|
std::map<const unsigned int, MMFFHerschbachLaurie>::const_iterator res2;
|
|
if (res1 != d_params.end()) {
|
|
res2 = ((*res1).second).find(canJRow);
|
|
if (res2 != ((*res1).second).end()) {
|
|
mmffHerschbachLaurieParams = &((*res2).second);
|
|
}
|
|
}
|
|
#else
|
|
auto bounds = std::equal_range(d_iRow.begin(), d_iRow.end(), canIRow);
|
|
if (bounds.first != bounds.second) {
|
|
bounds = std::equal_range(
|
|
d_jRow.begin() + (bounds.first - d_iRow.begin()),
|
|
d_jRow.begin() + (bounds.second - d_iRow.begin()), canJRow);
|
|
if (bounds.first != bounds.second) {
|
|
mmffHerschbachLaurieParams = &d_params[bounds.first - d_jRow.begin()];
|
|
}
|
|
}
|
|
#endif
|
|
|
|
return mmffHerschbachLaurieParams;
|
|
}
|
|
|
|
MMFFHerschbachLaurieCollection(std::string mmffHerschbachLaurie = "");
|
|
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
std::map<const unsigned int,
|
|
std::map<const unsigned int, MMFFHerschbachLaurie>>
|
|
d_params; //!< the parameter 2D-map
|
|
#else
|
|
std::vector<MMFFHerschbachLaurie> d_params; //!< the parameter vector
|
|
std::vector<std::uint8_t> d_iRow; //!< periodic row number vector for atom i
|
|
std::vector<std::uint8_t> d_jRow; //!< periodic row number vector for atom j
|
|
#endif
|
|
};
|
|
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFCovRadPauEleCollection {
|
|
public:
|
|
//! Looks up the parameters for a particular key and returns them.
|
|
/*!
|
|
\return a pointer to the MMFFCovRadPauEle object, NULL on failure.
|
|
*/
|
|
const MMFFCovRadPauEle *operator()(const unsigned int atomicNum) const {
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
const auto res = d_params.find(atomicNum);
|
|
return ((res != d_params.end()) ? &((*res).second) : NULL);
|
|
#else
|
|
auto bounds =
|
|
std::equal_range(d_atomicNum.begin(), d_atomicNum.end(), atomicNum);
|
|
return ((bounds.first != bounds.second)
|
|
? &d_params[bounds.first - d_atomicNum.begin()]
|
|
: nullptr);
|
|
#endif
|
|
}
|
|
|
|
MMFFCovRadPauEleCollection(std::string mmffCovRadPauEle = "");
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
std::map<const unsigned int, MMFFCovRadPauEle>
|
|
d_params; //!< the parameter map
|
|
#else
|
|
std::vector<MMFFCovRadPauEle> d_params; //!< the parameter vector
|
|
std::vector<std::uint8_t> d_atomicNum; //!< the atomic number vector
|
|
#endif
|
|
};
|
|
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFAngleCollection {
|
|
public:
|
|
//! Looks up the parameters for a particular key and returns them.
|
|
/*!
|
|
\return a pointer to the MMFFAngle object, NULL on failure.
|
|
*/
|
|
const MMFFAngle *operator()(const MMFFDefCollection *mmffDef,
|
|
const unsigned int angleType,
|
|
const unsigned int iAtomType,
|
|
const unsigned int jAtomType,
|
|
const unsigned int kAtomType) const {
|
|
const MMFFAngle *mmffAngleParams = nullptr;
|
|
unsigned int iter = 0;
|
|
|
|
// For bending of the i-j-k angle, a five-stage process based
|
|
// in the level combinations 1-1-1,2-2-2,3-2-3,4-2-4, and
|
|
// 5-2-5 is used. (MMFF.I, note 68, page 519)
|
|
// We skip 1-1-1 since Level 2 === Level 1
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
while ((iter < 4) && (!mmffAngleParams)) {
|
|
unsigned int canIAtomType = (*mmffDef)(iAtomType)->eqLevel[iter];
|
|
unsigned int canKAtomType = (*mmffDef)(kAtomType)->eqLevel[iter];
|
|
if (canIAtomType > canKAtomType) {
|
|
std::swap(canIAtomType, canKAtomType);
|
|
}
|
|
const auto res1 = d_params.find(angleType);
|
|
if (res1 != d_params.end()) {
|
|
const auto res2 = ((*res1).second).find(canIAtomType);
|
|
if (res2 != ((*res1).second).end()) {
|
|
const auto res3 = ((*res2).second).find(jAtomType);
|
|
if (res3 != ((*res2).second).end()) {
|
|
const auto res4 = ((*res3).second).find(canKAtomType);
|
|
if (res4 != ((*res3).second).end()) {
|
|
mmffAngleParams = &((*res4).second);
|
|
}
|
|
}
|
|
}
|
|
}
|
|
++iter;
|
|
}
|
|
#else
|
|
auto jBounds =
|
|
std::equal_range(d_jAtomType.begin(), d_jAtomType.end(), jAtomType);
|
|
if (jBounds.first != jBounds.second) {
|
|
while ((iter < 4) && (!mmffAngleParams)) {
|
|
unsigned int canIAtomType = (*mmffDef)(iAtomType)->eqLevel[iter];
|
|
unsigned int canKAtomType = (*mmffDef)(kAtomType)->eqLevel[iter];
|
|
if (canIAtomType > canKAtomType) {
|
|
std::swap(canIAtomType, canKAtomType);
|
|
}
|
|
|
|
auto bounds = std::equal_range(
|
|
d_iAtomType.begin() + (jBounds.first - d_jAtomType.begin()),
|
|
d_iAtomType.begin() + (jBounds.second - d_jAtomType.begin()),
|
|
canIAtomType);
|
|
if (bounds.first != bounds.second) {
|
|
bounds = std::equal_range(
|
|
d_kAtomType.begin() + (bounds.first - d_iAtomType.begin()),
|
|
d_kAtomType.begin() + (bounds.second - d_iAtomType.begin()),
|
|
canKAtomType);
|
|
if (bounds.first != bounds.second) {
|
|
bounds = std::equal_range(
|
|
d_angleType.begin() + (bounds.first - d_kAtomType.begin()),
|
|
d_angleType.begin() + (bounds.second - d_kAtomType.begin()),
|
|
angleType);
|
|
if (bounds.first != bounds.second) {
|
|
mmffAngleParams = &d_params[bounds.first - d_angleType.begin()];
|
|
}
|
|
}
|
|
}
|
|
++iter;
|
|
}
|
|
}
|
|
#endif
|
|
|
|
return mmffAngleParams;
|
|
}
|
|
|
|
MMFFAngleCollection(std::string mmffAngle = "");
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
std::map<const unsigned int,
|
|
std::map<const unsigned int,
|
|
std::map<const unsigned int,
|
|
std::map<const unsigned int, MMFFAngle>>>>
|
|
d_params; //!< the parameter 4D-map
|
|
#else
|
|
std::vector<MMFFAngle> d_params; //!< the parameter vector
|
|
std::vector<std::uint8_t> d_iAtomType; //!< atom type vector for atom i
|
|
std::vector<std::uint8_t> d_jAtomType; //!< atom type vector for atom j
|
|
std::vector<std::uint8_t> d_kAtomType; //!< atom type vector for atom k
|
|
std::vector<std::uint8_t> d_angleType; //!< angle type vector for angle i-j-k
|
|
#endif
|
|
};
|
|
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFStbnCollection {
|
|
public:
|
|
//! Looks up the parameters for a particular key and returns them.
|
|
/*!
|
|
\return a pointer to the MMFFStbn object, NULL on failure.
|
|
*/
|
|
const std::pair<bool, const MMFFStbn *> getMMFFStbnParams(
|
|
const unsigned int stretchBendType, const unsigned int bondType1,
|
|
const unsigned int bondType2, const unsigned int iAtomType,
|
|
const unsigned int jAtomType, const unsigned int kAtomType) const {
|
|
const MMFFStbn *mmffStbnParams = nullptr;
|
|
bool swap = false;
|
|
unsigned int canIAtomType = iAtomType;
|
|
unsigned int canKAtomType = kAtomType;
|
|
unsigned int canStretchBendType = stretchBendType;
|
|
if (iAtomType > kAtomType) {
|
|
canIAtomType = kAtomType;
|
|
canKAtomType = iAtomType;
|
|
swap = true;
|
|
} else if (iAtomType == kAtomType) {
|
|
swap = (bondType1 < bondType2);
|
|
}
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
const auto res1 = d_params.find(canStretchBendType);
|
|
if (res1 != d_params.end()) {
|
|
const auto res2 = ((*res1).second).find(canIAtomType);
|
|
if (res2 != ((*res1).second).end()) {
|
|
const auto res3 = ((*res2).second).find(jAtomType);
|
|
if (res3 != ((*res2).second).end()) {
|
|
const auto res4 = ((*res3).second).find(canKAtomType);
|
|
if (res4 != ((*res3).second).end()) {
|
|
mmffStbnParams = &((*res4).second);
|
|
}
|
|
}
|
|
}
|
|
}
|
|
#else
|
|
auto jBounds =
|
|
std::equal_range(d_jAtomType.begin(), d_jAtomType.end(), jAtomType);
|
|
if (jBounds.first != jBounds.second) {
|
|
auto bounds = std::equal_range(
|
|
d_iAtomType.begin() + (jBounds.first - d_jAtomType.begin()),
|
|
d_iAtomType.begin() + (jBounds.second - d_jAtomType.begin()),
|
|
canIAtomType);
|
|
if (bounds.first != bounds.second) {
|
|
bounds = std::equal_range(
|
|
d_kAtomType.begin() + (bounds.first - d_iAtomType.begin()),
|
|
d_kAtomType.begin() + (bounds.second - d_iAtomType.begin()),
|
|
canKAtomType);
|
|
if (bounds.first != bounds.second) {
|
|
bounds = std::equal_range(
|
|
d_stretchBendType.begin() + (bounds.first - d_kAtomType.begin()),
|
|
d_stretchBendType.begin() + (bounds.second - d_kAtomType.begin()),
|
|
canStretchBendType);
|
|
if (bounds.first != bounds.second) {
|
|
mmffStbnParams =
|
|
&d_params[bounds.first - d_stretchBendType.begin()];
|
|
}
|
|
}
|
|
}
|
|
}
|
|
#endif
|
|
|
|
return std::make_pair(swap, mmffStbnParams);
|
|
}
|
|
|
|
MMFFStbnCollection(std::string mmffStbn = "");
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
std::map<const unsigned int,
|
|
std::map<const unsigned int,
|
|
std::map<const unsigned int,
|
|
std::map<const unsigned int, MMFFStbn>>>>
|
|
d_params; //!< the parameter 4D-map
|
|
#else
|
|
std::vector<MMFFStbn> d_params; //!< the parameter vector
|
|
std::vector<std::uint8_t> d_iAtomType; //!< atom type vector for atom i
|
|
std::vector<std::uint8_t> d_jAtomType; //!< atom type vector for atom j
|
|
std::vector<std::uint8_t> d_kAtomType; //!< atom type vector for atom k
|
|
std::vector<std::uint8_t>
|
|
d_stretchBendType; //!< stretch-bend type vector for angle i-j-k
|
|
#endif
|
|
};
|
|
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFDfsbCollection {
|
|
public:
|
|
//! Looks up the parameters for a particular key and returns them.
|
|
/*!
|
|
\return a pointer to the MMFFStbn object, NULL on failure.
|
|
*/
|
|
const std::pair<bool, const MMFFStbn *> getMMFFDfsbParams(
|
|
const unsigned int periodicTableRow1,
|
|
const unsigned int periodicTableRow2,
|
|
const unsigned int periodicTableRow3) const {
|
|
const MMFFStbn *mmffDfsbParams = nullptr;
|
|
bool swap = false;
|
|
unsigned int canPeriodicTableRow1 = periodicTableRow1;
|
|
unsigned int canPeriodicTableRow3 = periodicTableRow3;
|
|
if (periodicTableRow1 > periodicTableRow3) {
|
|
canPeriodicTableRow1 = periodicTableRow3;
|
|
canPeriodicTableRow3 = periodicTableRow1;
|
|
swap = true;
|
|
}
|
|
const auto res1 = d_params.find(canPeriodicTableRow1);
|
|
if (res1 != d_params.end()) {
|
|
const auto res2 = ((*res1).second).find(periodicTableRow2);
|
|
if (res2 != ((*res1).second).end()) {
|
|
const auto res3 = ((*res2).second).find(canPeriodicTableRow3);
|
|
if (res3 != ((*res2).second).end()) {
|
|
mmffDfsbParams = &((*res3).second);
|
|
}
|
|
}
|
|
}
|
|
|
|
return std::make_pair(swap, mmffDfsbParams);
|
|
}
|
|
|
|
MMFFDfsbCollection(std::string mmffDfsb = "");
|
|
std::map<const unsigned int,
|
|
std::map<const unsigned int, std::map<const unsigned int, MMFFStbn>>>
|
|
d_params; //!< the parameter 3D-map
|
|
};
|
|
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFOopCollection {
|
|
public:
|
|
//! Looks up the parameters for a particular key and returns them.
|
|
/*!
|
|
\return a pointer to the MMFFOop object, NULL on failure.
|
|
*/
|
|
const MMFFOop *operator()(const MMFFDefCollection *mmffDef,
|
|
const unsigned int iAtomType,
|
|
const unsigned int jAtomType,
|
|
const unsigned int kAtomType,
|
|
const unsigned int lAtomType) const {
|
|
const MMFFOop *mmffOopParams = nullptr;
|
|
unsigned int iter = 0;
|
|
std::vector<unsigned int> canIKLAtomType(3);
|
|
// For out-of-plane bending ijk; I , where j is the central
|
|
// atom [cf. eq. (511, the five-stage protocol 1-1-1; 1, 2-2-2; 2,
|
|
// 3-2-3;3, 4-2-4;4, 5-2-5;5 is used. The final stage provides
|
|
// wild-card defaults for all except the central atom.
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
while ((iter < 4) && (!mmffOopParams)) {
|
|
canIKLAtomType[0] = (*mmffDef)(iAtomType)->eqLevel[iter];
|
|
unsigned int canJAtomType = jAtomType;
|
|
canIKLAtomType[1] = (*mmffDef)(kAtomType)->eqLevel[iter];
|
|
canIKLAtomType[2] = (*mmffDef)(lAtomType)->eqLevel[iter];
|
|
std::sort(canIKLAtomType.begin(), canIKLAtomType.end());
|
|
const auto res1 = d_params.find(canIKLAtomType[0]);
|
|
if (res1 != d_params.end()) {
|
|
const auto res2 = ((*res1).second).find(canJAtomType);
|
|
if (res2 != ((*res1).second).end()) {
|
|
const auto res3 = ((*res2).second).find(canIKLAtomType[1]);
|
|
if (res3 != ((*res2).second).end()) {
|
|
const auto res4 = ((*res3).second).find(canIKLAtomType[2]);
|
|
if (res4 != ((*res3).second).end()) {
|
|
mmffOopParams = &((*res4).second);
|
|
}
|
|
}
|
|
}
|
|
}
|
|
++iter;
|
|
}
|
|
#else
|
|
auto jBounds =
|
|
std::equal_range(d_jAtomType.begin(), d_jAtomType.end(), jAtomType);
|
|
if (jBounds.first != jBounds.second) {
|
|
while ((iter < 4) && (!mmffOopParams)) {
|
|
canIKLAtomType[0] = (*mmffDef)(iAtomType)->eqLevel[iter];
|
|
canIKLAtomType[1] = (*mmffDef)(kAtomType)->eqLevel[iter];
|
|
canIKLAtomType[2] = (*mmffDef)(lAtomType)->eqLevel[iter];
|
|
std::sort(canIKLAtomType.begin(), canIKLAtomType.end());
|
|
auto bounds = std::equal_range(
|
|
d_iAtomType.begin() + (jBounds.first - d_jAtomType.begin()),
|
|
d_iAtomType.begin() + (jBounds.second - d_jAtomType.begin()),
|
|
canIKLAtomType[0]);
|
|
if (bounds.first != bounds.second) {
|
|
bounds = std::equal_range(
|
|
d_kAtomType.begin() + (bounds.first - d_iAtomType.begin()),
|
|
d_kAtomType.begin() + (bounds.second - d_iAtomType.begin()),
|
|
canIKLAtomType[1]);
|
|
if (bounds.first != bounds.second) {
|
|
bounds = std::equal_range(
|
|
d_lAtomType.begin() + (bounds.first - d_kAtomType.begin()),
|
|
d_lAtomType.begin() + (bounds.second - d_kAtomType.begin()),
|
|
canIKLAtomType[2]);
|
|
if (bounds.first != bounds.second) {
|
|
mmffOopParams = &d_params[bounds.first - d_lAtomType.begin()];
|
|
}
|
|
}
|
|
}
|
|
++iter;
|
|
}
|
|
}
|
|
#endif
|
|
|
|
return mmffOopParams;
|
|
}
|
|
|
|
MMFFOopCollection(const bool isMMFFs, std::string mmffOop = "");
|
|
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
std::map<const unsigned int,
|
|
std::map<const unsigned int,
|
|
std::map<const unsigned int,
|
|
std::map<const unsigned int, MMFFOop>>>>
|
|
d_params; //!< the parameter 4D-map
|
|
#else
|
|
std::vector<MMFFOop> d_params; //!< the parameter vector
|
|
std::vector<std::uint8_t> d_iAtomType; //!< atom type vector for atom i
|
|
std::vector<std::uint8_t> d_jAtomType; //!< atom type vector for atom j
|
|
std::vector<std::uint8_t> d_kAtomType; //!< atom type vector for atom k
|
|
std::vector<std::uint8_t> d_lAtomType; //!< atom type vector for atom l
|
|
#endif
|
|
};
|
|
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFTorCollection {
|
|
public:
|
|
//! Looks up the parameters for a particular key and returns them.
|
|
/*!
|
|
\return a pointer to the MMFFTor object, NULL on failure.
|
|
*/
|
|
const std::pair<const unsigned int, const MMFFTor *> getMMFFTorParams(
|
|
const MMFFDefCollection *mmffDef,
|
|
const std::pair<unsigned int, unsigned int> torType,
|
|
const unsigned int iAtomType, const unsigned int jAtomType,
|
|
const unsigned int kAtomType, const unsigned int lAtomType) const {
|
|
const MMFFTor *mmffTorParams = nullptr;
|
|
unsigned int iter = 0;
|
|
unsigned int iWildCard = 0;
|
|
unsigned int lWildCard = 0;
|
|
unsigned int canTorType = torType.first;
|
|
unsigned int maxIter = 5;
|
|
// For i-j-k-2 torsion interactions, a five-stage
|
|
// process based on level combinations 1-1-1-1, 2-2-2-2,
|
|
// 3-2-2-5, 5-2-2-3, and 5-2-2-5 is used, where stages 3
|
|
// and 4 correspond to "half-default" or "half-wild-card" entries.
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
#else
|
|
#endif
|
|
|
|
while (((iter < maxIter) && ((!mmffTorParams) || (maxIter == 4))) ||
|
|
((iter == 4) && (torType.first == 5) && torType.second)) {
|
|
// The rule of setting the torsion type to the value it had
|
|
// before being set to 5 as a last resort in case parameters
|
|
// could not be found is not mentioned in MMFF.IV; it was
|
|
// empirically discovered due to a number of tests in the
|
|
// MMFF validation suite otherwise failing
|
|
if ((maxIter == 5) && (iter == 4)) {
|
|
maxIter = 4;
|
|
iter = 0;
|
|
canTorType = torType.second;
|
|
}
|
|
iWildCard = iter;
|
|
lWildCard = iter;
|
|
if (iter == 1) {
|
|
iWildCard = 1;
|
|
lWildCard = 3;
|
|
} else if (iter == 2) {
|
|
iWildCard = 3;
|
|
lWildCard = 1;
|
|
}
|
|
unsigned int canIAtomType = (*mmffDef)(iAtomType)->eqLevel[iWildCard];
|
|
unsigned int canJAtomType = jAtomType;
|
|
unsigned int canKAtomType = kAtomType;
|
|
unsigned int canLAtomType = (*mmffDef)(lAtomType)->eqLevel[lWildCard];
|
|
if (canJAtomType > canKAtomType) {
|
|
unsigned int temp = canKAtomType;
|
|
canKAtomType = canJAtomType;
|
|
canJAtomType = temp;
|
|
temp = canLAtomType;
|
|
canLAtomType = canIAtomType;
|
|
canIAtomType = temp;
|
|
} else if ((canJAtomType == canKAtomType) &&
|
|
(canIAtomType > canLAtomType)) {
|
|
unsigned int temp = canLAtomType;
|
|
canLAtomType = canIAtomType;
|
|
canIAtomType = temp;
|
|
}
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
const auto res1 = d_params.find(canTorType);
|
|
if (res1 != d_params.end()) {
|
|
const auto res2 = ((*res1).second).find(canIAtomType);
|
|
if (res2 != ((*res1).second).end()) {
|
|
const auto res3 = ((*res2).second).find(canJAtomType);
|
|
if (res3 != ((*res2).second).end()) {
|
|
const auto res4 = ((*res3).second).find(canKAtomType);
|
|
if (res4 != ((*res3).second).end()) {
|
|
const auto res5 = ((*res4).second).find(canLAtomType);
|
|
if (res5 != ((*res4).second).end()) {
|
|
mmffTorParams = &((*res5).second);
|
|
if (maxIter == 4) {
|
|
break;
|
|
}
|
|
}
|
|
}
|
|
}
|
|
}
|
|
}
|
|
#else
|
|
auto jBounds = std::equal_range(d_jAtomType.begin(), d_jAtomType.end(),
|
|
canJAtomType);
|
|
if (jBounds.first != jBounds.second) {
|
|
auto bounds = std::equal_range(
|
|
d_kAtomType.begin() + (jBounds.first - d_jAtomType.begin()),
|
|
d_kAtomType.begin() + (jBounds.second - d_jAtomType.begin()),
|
|
canKAtomType);
|
|
if (bounds.first != bounds.second) {
|
|
bounds = std::equal_range(
|
|
d_iAtomType.begin() + (bounds.first - d_kAtomType.begin()),
|
|
d_iAtomType.begin() + (bounds.second - d_kAtomType.begin()),
|
|
canIAtomType);
|
|
if (bounds.first != bounds.second) {
|
|
bounds = std::equal_range(
|
|
d_lAtomType.begin() + (bounds.first - d_iAtomType.begin()),
|
|
d_lAtomType.begin() + (bounds.second - d_iAtomType.begin()),
|
|
canLAtomType);
|
|
if (bounds.first != bounds.second) {
|
|
bounds = std::equal_range(
|
|
d_torType.begin() + (bounds.first - d_lAtomType.begin()),
|
|
d_torType.begin() + (bounds.second - d_lAtomType.begin()),
|
|
canTorType);
|
|
if (bounds.first != bounds.second) {
|
|
mmffTorParams = &d_params[bounds.first - d_torType.begin()];
|
|
if (maxIter == 4) {
|
|
break;
|
|
}
|
|
}
|
|
}
|
|
}
|
|
}
|
|
}
|
|
#endif
|
|
++iter;
|
|
}
|
|
|
|
return std::make_pair(canTorType, mmffTorParams);
|
|
}
|
|
|
|
MMFFTorCollection(const bool isMMFFs, std::string mmffTor = "");
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
std::map<
|
|
const unsigned int,
|
|
std::map<
|
|
const unsigned int,
|
|
std::map<const unsigned int,
|
|
std::map<const unsigned int, std::map<const unsigned int,
|
|
MMFFTor>>>>>
|
|
d_params; //!< the parameter 5D-map
|
|
#else
|
|
std::vector<MMFFTor> d_params; //!< the parameter vector
|
|
std::vector<std::uint8_t> d_iAtomType; //!< atom type vector for atom i
|
|
std::vector<std::uint8_t> d_jAtomType; //!< atom type vector for atom j
|
|
std::vector<std::uint8_t> d_kAtomType; //!< atom type vector for atom k
|
|
std::vector<std::uint8_t> d_lAtomType; //!< atom type vector for atom l
|
|
std::vector<std::uint8_t>
|
|
d_torType; //!< torsion type vector for angle i-j-k-l
|
|
#endif
|
|
};
|
|
|
|
class RDKIT_FORCEFIELD_EXPORT MMFFVdWCollection {
|
|
public:
|
|
//! gets a pointer to the singleton MMFFVdWCollection
|
|
double power;
|
|
double B;
|
|
double Beta;
|
|
double DARAD;
|
|
double DAEPS;
|
|
//! Looks up the parameters for a particular key and returns them.
|
|
/*!
|
|
\return a pointer to the MMFFVdW object, NULL on failure.
|
|
*/
|
|
const MMFFVdW *operator()(const unsigned int atomType) const {
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
const auto res = d_params.find(atomType);
|
|
return (res != d_params.end() ? &((*res).second) : NULL);
|
|
#else
|
|
auto bounds =
|
|
std::equal_range(d_atomType.begin(), d_atomType.end(), atomType);
|
|
return ((bounds.first != bounds.second)
|
|
? &d_params[bounds.first - d_atomType.begin()]
|
|
: nullptr);
|
|
#endif
|
|
}
|
|
|
|
MMFFVdWCollection(std::string mmffVdW = "");
|
|
#ifdef RDKIT_MMFF_PARAMS_USE_STD_MAP
|
|
std::map<const unsigned int, MMFFVdW> d_params; //!< the parameter map
|
|
#else
|
|
std::vector<MMFFVdW> d_params; //!< the parameter vector
|
|
std::vector<std::uint8_t> d_atomType; //!< atom type vector
|
|
#endif
|
|
};
|
|
} // namespace MMFF
|
|
} // namespace ForceFields
|
|
|
|
#endif
|