RDKit
Open-source cheminformatics and machine learning.
FFConvenience.h
Go to the documentation of this file.
1 //
2 // Copyright (C) 2019 Paolo Tosco
3 //
4 // @@ All Rights Reserved @@
5 // This file is part of the RDKit.
6 // The contents are covered by the terms of the BSD license
7 // which is included in the file license.txt, found at the root
8 // of the RDKit source tree.
9 //
10 #include <RDGeneral/export.h>
11 #ifndef RD_FFCONVENIENCE_H
12 #define RD_FFCONVENIENCE_H
13 #include <ForceField/ForceField.h>
14 #include <RDGeneral/RDThreads.h>
15 
16 namespace RDKit {
17 class ROMol;
18 namespace ForceFieldsHelper {
19 namespace detail {
20 #ifdef RDK_THREADSAFE_SSS
21 void OptimizeMoleculeConfsHelper_(ForceFields::ForceField ff, ROMol *mol,
22  std::vector<std::pair<int, double>> *res,
23  unsigned int threadIdx,
24  unsigned int numThreads, int maxIters) {
25  PRECONDITION(mol, "mol must not be nullptr");
26  PRECONDITION(res, "res must not be nullptr");
27  PRECONDITION(res->size() >= mol->getNumConformers(), "res->size() must be >= mol->getNumConformers()");
28  unsigned int i = 0;
29  ff.positions().resize(mol->getNumAtoms());
30  for (ROMol::ConformerIterator cit = mol->beginConformers();
31  cit != mol->endConformers(); ++cit, ++i) {
32  if (i % numThreads != threadIdx) continue;
33  for (unsigned int aidx = 0; aidx < mol->getNumAtoms(); ++aidx) {
34  ff.positions()[aidx] = &(*cit)->getAtomPos(aidx);
35  }
36  ff.initialize();
37  int needsMore = ff.minimize(maxIters);
38  double e = ff.calcEnergy();
39  (*res)[i] = std::make_pair(needsMore, e);
40  }
41 }
42 
43 void OptimizeMoleculeConfsMT(ROMol &mol, const ForceFields::ForceField &ff,
44  std::vector<std::pair<int, double>> &res,
45  int numThreads, int maxIters) {
46  std::vector<std::thread> tg;
47  for (int ti = 0; ti < numThreads; ++ti) {
48  tg.emplace_back(std::thread(detail::OptimizeMoleculeConfsHelper_,
49  ff, &mol, &res, ti, numThreads, maxIters));
50  }
51  for (auto &thread : tg) {
52  if (thread.joinable()) thread.join();
53  }
54 }
55 #endif
56 
58  std::vector<std::pair<int, double>> &res,
59  int maxIters) {
60  PRECONDITION(res.size() >= mol.getNumConformers(), "res.size() must be >= mol.getNumConformers()");
61  unsigned int i = 0;
62  for (ROMol::ConformerIterator cit = mol.beginConformers();
63  cit != mol.endConformers(); ++cit, ++i) {
64  for (unsigned int aidx = 0; aidx < mol.getNumAtoms(); ++aidx) {
65  ff.positions()[aidx] = &(*cit)->getAtomPos(aidx);
66  }
67  ff.initialize();
68  int needsMore = ff.minimize(maxIters);
69  double e = ff.calcEnergy();
70  res[i] = std::make_pair(needsMore, e);
71  }
72 }
73 } // end of detail namespace
74 
75 //! Convenience function for optimizing a molecule using a pre-generated force-field
76 /*
77  \param ff the force-field
78  \param res vector of (needsMore,energy) pairs
79  \param maxIters the maximum number of force-field iterations
80 
81  \return a pair with:
82  first: -1 if parameters were missing, 0 if the optimization converged, 1 if
83  more iterations are required.
84  second: the energy
85 */
86 std::pair<int, double> OptimizeMolecule(ForceFields::ForceField &ff, int maxIters = 1000) {
87  ff.initialize();
88  int res = ff.minimize(maxIters);
89  double e = ff.calcEnergy();
90  return std::make_pair(res, e);
91 }
92 
93 //! Convenience function for optimizing all of a molecule's conformations using
94 // a pre-generated force-field
95 /*
96  \param mol the molecule to use
97  \param ff the force-field
98  \param res vector of (needsMore,energy) pairs
99  \param numThreads the number of simultaneous threads to use (only has an
100  effect if the RDKit is compiled with thread support).
101  If set to zero, the max supported by the system will be
102  used.
103  \param maxIters the maximum number of force-field iterations
104 
105 */
107  std::vector<std::pair<int, double>> &res,
108  int numThreads = 1, int maxIters = 1000) {
109  res.resize(mol.getNumConformers());
110  numThreads = getNumThreadsToUse(numThreads);
111  if (numThreads == 1) {
112  detail::OptimizeMoleculeConfsST(mol, ff, res, maxIters);
113  }
114 #ifdef RDK_THREADSAFE_SSS
115  else {
116  detail::OptimizeMoleculeConfsMT(mol, ff, res, numThreads, maxIters);
117  }
118 #endif
119 }
120 } // end of namespace ForceFieldsHelper
121 } // end of namespace RDKit
122 #endif
RDKit::ForceFieldsHelper::detail::OptimizeMoleculeConfsST
void OptimizeMoleculeConfsST(ROMol &mol, ForceFields::ForceField &ff, std::vector< std::pair< int, double >> &res, int maxIters)
Definition: FFConvenience.h:57
ForceFields::ForceField::initialize
void initialize()
does initialization
RDKit::ROMol::getNumConformers
unsigned int getNumConformers() const
Definition: ROMol.h:443
ForceFields::ForceField
A class to store forcefields and handle minimization.
Definition: ForceField.h:79
RDKit::ROMol
Definition: ROMol.h:171
ForceFields::ForceField::minimize
int minimize(unsigned int snapshotFreq, RDKit::SnapshotVect *snapshotVect, unsigned int maxIts=200, double forceTol=1e-4, double energyTol=1e-6)
minimizes the energy of the system by following gradients
RDKit::ForceFieldsHelper::OptimizeMolecule
std::pair< int, double > OptimizeMolecule(ForceFields::ForceField &ff, int maxIters=1000)
Convenience function for optimizing a molecule using a pre-generated force-field.
Definition: FFConvenience.h:86
ForceFields::ForceField::positions
RDGeom::PointPtrVect & positions()
returns a reference to our points (a PointPtrVect)
Definition: ForceField.h:182
RDThreads.h
RDKit::ROMol::beginConformers
ConformerIterator beginConformers()
Definition: ROMol.h:615
RDKit::ROMol::getNumAtoms
unsigned int getNumAtoms(bool onlyExplicit=1) const
returns our number of atoms
RDKit
Std stuff.
Definition: Atom.h:30
RDKit::ForceFieldsHelper::OptimizeMoleculeConfs
void OptimizeMoleculeConfs(ROMol &mol, ForceFields::ForceField &ff, std::vector< std::pair< int, double >> &res, int numThreads=1, int maxIters=1000)
Convenience function for optimizing all of a molecule's conformations using.
Definition: FFConvenience.h:106
ForceFields::ForceField::calcEnergy
double calcEnergy(std::vector< double > *contribs=NULL) const
calculates and returns the energy (in kcal/mol) based on existing
PRECONDITION
#define PRECONDITION(expr, mess)
Definition: Invariant.h:109
ForceField.h
RDKit::ROMol::endConformers
ConformerIterator endConformers()
Definition: ROMol.h:617
RDKit::getNumThreadsToUse
unsigned int getNumThreadsToUse(int target)
Definition: RDThreads.h:37
export.h