00001 /********************************************************************** 00002 forcefield.h - Handle OBForceField class. 00003 00004 Copyright (C) 2006-2007 by Tim Vandermeersch <[email protected]> 00005 00006 This file is part of the Open Babel project. 00007 For more information, see <http://openbabel.sourceforge.net/> 00008 00009 This program is free software; you can redistribute it and/or modify 00010 it under the terms of the GNU General Public License as published by 00011 the Free Software Foundation version 2 of the License. 00012 00013 This program is distributed in the hope that it will be useful, 00014 but WITHOUT ANY WARRANTY; without even the implied warranty of 00015 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00016 GNU General Public License for more details. 00017 ***********************************************************************/ 00018 00019 #ifndef OB_FORCEFIELD_H 00020 #define OB_FORCEFIELD_H 00021 00022 #include <vector> 00023 #include <string> 00024 #include <map> 00025 00026 #include <list> 00027 #include <set> 00028 #include <openbabel/babelconfig.h> 00029 #include <openbabel/base.h> 00030 #include <openbabel/mol.h> 00031 #include <openbabel/plugin.h> 00032 #include <openbabel/grid.h> 00033 #include <openbabel/griddata.h> 00034 #include <float.h> 00035 00036 namespace OpenBabel 00037 { 00038 // log levels 00039 #define OBFF_LOGLVL_NONE 0 00040 #define OBFF_LOGLVL_LOW 1 00041 #define OBFF_LOGLVL_MEDIUM 2 00042 #define OBFF_LOGLVL_HIGH 3 00043 00044 // terms 00045 #define OBFF_ENERGY (1 << 0) 00046 #define OBFF_EBOND (1 << 1) 00047 #define OBFF_EANGLE (1 << 2) 00048 #define OBFF_ESTRBND (1 << 3) 00049 #define OBFF_ETORSION (1 << 4) 00050 #define OBFF_EOOP (1 << 5) 00051 #define OBFF_EVDW (1 << 6) 00052 #define OBFF_EELECTROSTATIC (1 << 7) 00053 00054 // constraint types 00055 #define OBFF_CONST_IGNORE (1 << 0) 00056 #define OBFF_CONST_ATOM (1 << 1) 00057 #define OBFF_CONST_ATOM_X (1 << 2) 00058 #define OBFF_CONST_ATOM_Y (1 << 3) 00059 #define OBFF_CONST_ATOM_Z (1 << 4) 00060 #define OBFF_CONST_DISTANCE (1 << 5) 00061 #define OBFF_CONST_ANGLE (1 << 6) 00062 #define OBFF_CONST_TORSION (1 << 7) 00063 #define OBFF_CONST_CHIRAL (1 << 8) 00064 00065 // mode arguments for SteepestDescent, ConjugateGradients, ... 00066 #define OBFF_NUMERICAL_GRADIENT (1 << 0) 00067 #define OBFF_ANALYTICAL_GRADIENT (1 << 1) 00068 00069 #define KCAL_TO_KJ 4.1868 00070 00071 // inline if statements for logging. 00072 #define IF_OBFF_LOGLVL_LOW if(_loglvl >= OBFF_LOGLVL_LOW) 00073 #define IF_OBFF_LOGLVL_MEDIUM if(_loglvl >= OBFF_LOGLVL_MEDIUM) 00074 #define IF_OBFF_LOGLVL_HIGH if(_loglvl >= OBFF_LOGLVL_HIGH) 00075 00077 struct LineSearchType 00078 { 00079 enum { 00080 Simple, Newton2Num 00081 }; 00082 }; 00083 /* 00084 struct ConstraintType 00085 { 00086 enum { 00087 Ignore, Atom, AtomX, AtomY, AtomZ, Distance, Angle, Torsion, Chiral 00088 }; 00089 }; 00090 */ 00091 00094 class OBFPRT OBFFParameter { 00095 public: 00097 int a, b, c, d; 00099 std::string _a, _b, _c, _d; 00100 00102 std::vector<int> _ipar; 00104 std::vector<double> _dpar; 00105 00107 OBFFParameter& operator=(const OBFFParameter &ai) 00108 { 00109 if (this != &ai) { 00110 a = ai.a; 00111 b = ai.b; 00112 c = ai.c; 00113 d = ai.d; 00114 _a = ai._a; 00115 _b = ai._b; 00116 _c = ai._c; 00117 _d = ai._d; 00118 _ipar = ai._ipar; 00119 _dpar = ai._dpar; 00120 } 00121 00122 return *this; 00123 } 00124 00126 void clear () 00127 { 00128 a = b = c = d = 0; 00129 _ipar.clear(); 00130 _dpar.clear(); 00131 } 00132 }; // class OBFFParameter 00133 00134 // specific class introductions in forcefieldYYYY.cpp (for YYYY calculations) 00135 00138 class OBFPRT OBFFCalculation2 00139 { 00140 public: 00142 double energy; 00144 OBAtom *a, *b; 00146 int idx_a, idx_b; 00148 double *pos_a, *pos_b; 00150 double force_a[3], force_b[3]; 00152 virtual ~OBFFCalculation2() 00153 { 00154 } 00157 virtual void SetupPointers() 00158 { 00159 if (!a || !b) return; 00160 pos_a = a->GetCoordinate(); 00161 idx_a = a->GetIdx(); 00162 pos_b = b->GetCoordinate(); 00163 idx_b = b->GetIdx(); 00164 } 00165 }; 00166 00169 class OBFPRT OBFFCalculation3: public OBFFCalculation2 00170 { 00171 public: 00173 OBAtom *c; 00175 int idx_c; 00177 double *pos_c; 00179 double force_c[3]; 00181 virtual ~OBFFCalculation3() 00182 { 00183 } 00186 virtual void SetupPointers() 00187 { 00188 if (!a || !b || !c) return; 00189 pos_a = a->GetCoordinate(); 00190 idx_a = a->GetIdx(); 00191 pos_b = b->GetCoordinate(); 00192 idx_b = b->GetIdx(); 00193 pos_c = c->GetCoordinate(); 00194 idx_c = c->GetIdx(); 00195 } 00196 }; 00197 00200 class OBFPRT OBFFCalculation4: public OBFFCalculation3 00201 { 00202 public: 00204 OBAtom *d; 00206 int idx_d; 00208 double *pos_d; 00210 double force_d[3]; 00212 virtual ~OBFFCalculation4() 00213 { 00214 } 00217 void SetupPointers() 00218 { 00219 if (!a || !b || !c || !d) return; 00220 pos_a = a->GetCoordinate(); 00221 idx_a = a->GetIdx(); 00222 pos_b = b->GetCoordinate(); 00223 idx_b = b->GetIdx(); 00224 pos_c = c->GetCoordinate(); 00225 idx_c = c->GetIdx(); 00226 pos_d = d->GetCoordinate(); 00227 idx_d = d->GetIdx(); 00228 } 00229 }; 00230 00234 class OBFPRT OBFFConstraint 00235 { 00236 public: 00238 double factor, constraint_value; 00239 double rab0, rbc0; 00241 int type, ia, ib, ic, id; 00243 OBAtom *a, *b, *c, *d; 00245 vector3 grada, gradb, gradc, gradd; 00246 00248 OBFFConstraint() 00249 { 00250 a = b = c = d = NULL; 00251 ia = ib = ic = id = 0; 00252 constraint_value = 0.0; 00253 factor = 0.0; 00254 } 00256 ~OBFFConstraint() 00257 { 00258 } 00259 00260 vector3 GetGradient(int a) 00261 { 00262 if (a == ia) 00263 return grada; 00264 else if (a == ib) 00265 return gradb; 00266 else if (a == ic) 00267 return gradc; 00268 else if (a == id) 00269 return gradd; 00270 else 00271 return VZero; 00272 } 00273 }; 00274 00278 class OBFPRT OBFFConstraints 00279 { 00280 public: 00282 OBFFConstraints(); 00284 ~OBFFConstraints() 00285 { 00286 _constraints.clear(); 00287 _ignored.Clear(); 00288 _fixed.Clear(); 00289 _Xfixed.Clear(); 00290 _Yfixed.Clear(); 00291 _Zfixed.Clear(); 00292 } 00294 void Clear(); 00296 double GetConstraintEnergy(); 00298 vector3 GetGradient(int a); 00300 OBFFConstraints& operator=(const OBFFConstraints &ai) 00301 { 00302 if (this != &ai) { 00303 _constraints = ai._constraints; 00304 _ignored = ai._ignored; 00305 _fixed = ai._fixed; 00306 _Xfixed = ai._Xfixed; 00307 _Yfixed = ai._Yfixed; 00308 _Zfixed = ai._Zfixed; 00309 } 00310 return *this; 00311 } 00312 00316 void Setup(OBMol &mol); 00317 00319 // Set Constraints // 00322 00323 00324 void SetFactor(double factor); 00326 void AddIgnore(int a); 00328 void AddAtomConstraint(int a); 00330 void AddAtomXConstraint(int a); 00332 void AddAtomYConstraint(int a); 00334 void AddAtomZConstraint(int a); 00336 void AddDistanceConstraint(int a, int b, double length); 00338 void AddAngleConstraint(int a, int b, int c, double angle); 00340 void AddTorsionConstraint(int a, int b, int c, int d, double torsion); 00343 void DeleteConstraint(int index); 00345 00346 // Get Constraints // 00349 00350 00351 double GetFactor(); 00353 int Size() const; 00363 int GetConstraintType(int index) const; 00367 double GetConstraintValue(int index) const; 00370 int GetConstraintAtomA(int index) const; 00373 int GetConstraintAtomB(int index) const; 00376 int GetConstraintAtomC(int index) const; 00379 int GetConstraintAtomD(int index) const; 00382 bool IsIgnored(int a); 00385 bool IsFixed(int a); 00388 bool IsXFixed(int a); 00391 bool IsYFixed(int a); 00394 bool IsZFixed(int a); 00398 OBBitVec GetIgnoredBitVec() { return _ignored; } 00400 00401 private: 00402 std::vector<OBFFConstraint> _constraints; 00403 OBBitVec _ignored; 00404 OBBitVec _fixed; 00405 OBBitVec _Xfixed; 00406 OBBitVec _Yfixed; 00407 OBBitVec _Zfixed; 00408 double _factor; 00409 }; 00410 00411 // Class OBForceField 00412 // class introduction in forcefield.cpp 00413 class OBFPRT OBForceField : public OBPlugin 00414 { 00415 MAKE_PLUGIN(OBForceField) 00416 00417 protected: 00418 00460 OBFFParameter* GetParameter(int a, int b, int c, int d, std::vector<OBFFParameter> ¶meter); 00462 OBFFParameter* GetParameter(const char* a, const char* b, const char* c, const char* d, 00463 std::vector<OBFFParameter> ¶meter); 00465 int GetParameterIdx(int a, int b, int c, int d, std::vector<OBFFParameter> ¶meter); 00466 00475 vector3 NumericalDerivative(OBAtom *a, int terms = OBFF_ENERGY); 00477 vector3 NumericalSecondDerivative(OBAtom *a, int terms = OBFF_ENERGY); 00478 00479 /* 00480 * NEW gradients functions 00481 */ 00482 00485 void SetGradient(double *grad, int idx) 00486 { 00487 const int coordIdx = (idx - 1) * 3; 00488 for (unsigned int i = 0; i < 3; ++i) { 00489 _gradientPtr[coordIdx + i] = grad[i]; 00490 } 00491 } 00492 00495 void AddGradient(double *grad, int idx) 00496 { 00497 const int coordIdx = (idx - 1) * 3; 00498 for (unsigned int i = 0; i < 3; ++i) { 00499 _gradientPtr[coordIdx + i] += grad[i]; 00500 } 00501 } 00502 00505 virtual vector3 GetGradient(OBAtom *a, int /*terms*/ = OBFF_ENERGY) 00506 { 00507 const int coordIdx = (a->GetIdx() - 1) * 3; 00508 return _gradientPtr + coordIdx; 00509 } 00510 00513 double* GetGradientPtr() 00514 { 00515 return _gradientPtr; 00516 } 00517 00520 virtual void ClearGradients() 00521 { 00522 // We cannot use memset because IEEE floating point representations 00523 // are not guaranteed by C/C++ standard, but this loop can be 00524 // unrolled or vectorized by compilers 00525 for (unsigned int i = 0; i < _ncoords; ++i) 00526 _gradientPtr[i] = 0.0; 00527 // memset(_gradientPtr, '\0', sizeof(double)*_ncoords); 00528 } 00529 00537 bool IsInSameRing(OBAtom* a, OBAtom* b); 00538 00539 // general variables 00540 OBMol _mol; 00541 bool _init; 00542 std::string _parFile; 00543 bool _validSetup; 00544 double *_gradientPtr; 00545 // logging variables 00546 std::ostream* _logos; 00547 char _logbuf[BUFF_SIZE+1]; 00548 int _loglvl; 00549 int _origLogLevel; 00550 // conformer genereation (rotor search) variables 00551 int _current_conformer; 00552 std::vector<double> _energies; 00553 // minimization variables 00554 double _econv, _e_n1; 00555 int _cstep, _nsteps; 00556 double *_grad1; 00557 unsigned int _ncoords; 00558 int _linesearch; 00559 // molecular dynamics variables 00560 double _timestep; 00561 double _temp; 00562 double *_velocityPtr; 00563 // contraint varibles 00564 static OBFFConstraints _constraints; 00565 static int _fixAtom; 00566 static int _ignoreAtom; 00567 // cut-off variables 00568 bool _cutoff; 00569 double _rvdw; 00570 double _rele; 00571 OBBitVec _vdwpairs; 00572 OBBitVec _elepairs; 00573 int _pairfreq; 00574 // group variables 00575 std::vector<OBBitVec> _intraGroup; 00576 std::vector<OBBitVec> _interGroup; 00577 std::vector<std::pair<OBBitVec, OBBitVec> > _interGroups; 00578 00579 public: 00583 virtual OBForceField* MakeNewInstance()=0; 00584 00586 virtual ~OBForceField() 00587 { 00588 if (_grad1 != NULL) { 00589 delete [] _grad1; 00590 _grad1 = NULL; 00591 } 00592 } 00593 00595 const char* TypeID() 00596 { 00597 return "forcefields"; 00598 } 00599 00603 static OBForceField* FindForceField(const std::string& ID) 00604 { 00605 return FindType(ID.c_str()); 00606 } 00610 static OBForceField* FindForceField(const char *ID) 00611 { 00612 return FindType(ID); 00613 } 00614 /* 00615 * 00616 */ 00617 void SetParameterFile(const std::string &filename) 00618 { 00619 _parFile = filename; 00620 _init = false; 00621 } 00624 virtual std::string GetUnit() { return std::string("au"); } 00625 /* Does this force field have analytical gradients defined for all 00626 * calculation components (bonds, angles, non-bonded, etc.) 00627 * If this is true, code should default to using OBFF_ANALYTICAL_GRADIENT 00628 * for SteepestDescent() or ConjugateGradients(). 00629 * \return True if all analytical gradients are implemented. 00630 */ 00631 virtual bool HasAnalyticalGradients() { return false; } 00636 bool Setup(OBMol &mol); 00642 bool Setup(OBMol &mol, OBFFConstraints &constraints); 00646 // move to protected in future version 00647 virtual bool ParseParamFile() { return false; } 00651 // move to protected in future version 00652 virtual bool SetTypes() { return false; } 00656 // move to protected in future version 00657 virtual bool SetFormalCharges() { return false; } 00661 // move to protected in future version 00662 virtual bool SetPartialCharges() { return false; } 00666 // move to protected in future version 00667 virtual bool SetupCalculations() { return false; } 00672 // move to protected in future version 00673 virtual bool SetupPointers() { return false; } 00679 bool IsSetupNeeded(OBMol &mol); 00695 bool GetAtomTypes(OBMol &mol); 00711 bool GetPartialCharges(OBMol &mol); 00712 00713 00714 00719 bool GetCoordinates(OBMol &mol); 00721 bool UpdateCoordinates(OBMol &mol) {return GetCoordinates(mol); } 00726 bool GetConformers(OBMol &mol); 00728 bool UpdateConformers(OBMol &mol) { return GetConformers(mol); } 00733 bool SetCoordinates(OBMol &mol); 00738 bool SetConformers(OBMol &mol); 00748 OBGridData *GetGrid(double step, double padding, const char *type, double pchg); 00749 00751 // Interacting groups // 00753 00755 00756 00760 void AddIntraGroup(OBBitVec &group); 00765 void AddInterGroup(OBBitVec &group); 00773 void AddInterGroups(OBBitVec &group1, OBBitVec &group2); 00776 void ClearGroups(); 00779 bool HasGroups(); 00781 00783 // Cut-off // 00785 00787 00788 00791 void EnableCutOff(bool enable) 00792 { 00793 _cutoff = enable; 00794 } 00797 bool IsCutOffEnabled() 00798 { 00799 return _cutoff; 00800 } 00804 void SetVDWCutOff(double r) 00805 { 00806 _rvdw = r; 00807 } 00811 double GetVDWCutOff() 00812 { 00813 return _rvdw; 00814 } 00819 void SetElectrostaticCutOff(double r) 00820 { 00821 _rele = r; 00822 } 00826 double GetElectrostaticCutOff() 00827 { 00828 return _rele; 00829 } 00835 void SetUpdateFrequency(int f) 00836 { 00837 _pairfreq = f; 00838 } 00842 int GetUpdateFrequency() 00843 { 00844 return _pairfreq; 00845 } 00850 void UpdatePairsSimple(); 00851 00852 //void UpdatePairsGroup(); TODO 00853 00857 unsigned int GetNumPairs(); 00861 void EnableAllPairs() 00862 { 00863 // TODO: OBBitVec doesn't seem to be allocating it's memory correctly 00864 //_vdwpairs.SetRangeOn(0, _numpairs-1); 00865 //_elepairs.SetRangeOn(0, _numpairs-1); 00866 } 00868 00870 // Energy Evaluation // 00872 00874 00875 00884 virtual double Energy(bool = true) { return 0.0f; } 00891 virtual double E_Bond(bool = true) { return 0.0f; } 00898 virtual double E_Angle(bool = true) { return 0.0f; } 00905 virtual double E_StrBnd(bool = true) { return 0.0f; } 00912 virtual double E_Torsion(bool = true) { return 0.0f; } 00919 virtual double E_OOP(bool = true) { return 0.0f; } 00926 virtual double E_VDW(bool = true) { return 0.0f; } 00933 virtual double E_Electrostatic(bool = true) { return 0.0f; } 00935 00937 // Logging // 00939 00941 00942 00944 void PrintTypes(); 00948 void PrintFormalCharges(); 00951 void PrintPartialCharges(); 00954 void PrintVelocities(); 00959 bool SetLogFile(std::ostream *pos); 00984 bool SetLogLevel(int level); 00987 int GetLogLevel() { return _loglvl; } 00991 void OBFFLog(std::string msg) 00992 { 00993 if (!_logos) 00994 return; 00995 00996 *_logos << msg; 00997 } 01001 void OBFFLog(const char *msg) 01002 { 01003 if (!_logos) 01004 return; 01005 01006 *_logos << msg; 01007 } 01009 01011 // Structure Generation // 01013 01015 01016 01017 void DistanceGeometry(); 01037 void SystematicRotorSearch(unsigned int geomSteps = 2500); 01055 int SystematicRotorSearchInitialize(unsigned int geomSteps = 2500); 01060 bool SystematicRotorSearchNextConformer(unsigned int geomSteps = 2500); 01081 void RandomRotorSearch(unsigned int conformers, unsigned int geomSteps = 2500); 01099 void RandomRotorSearchInitialize(unsigned int conformers, unsigned int geomSteps = 2500); 01104 bool RandomRotorSearchNextConformer(unsigned int geomSteps = 2500); 01126 void WeightedRotorSearch(unsigned int conformers, unsigned int geomSteps); 01127 01129 // Energy Minimization // 01131 01133 01134 01137 void SetLineSearchType(int type) 01138 { 01139 _linesearch = type; 01140 } 01144 int GetLineSearchType() 01145 { 01146 return _linesearch; 01147 } 01151 vector3 LineSearch(OBAtom *atom, vector3 &direction); 01165 double LineSearch(double *currentCoords, double *direction); 01178 double Newton2NumLineSearch(double *direction); 01184 void LineSearchTakeStep(double *origCoords, double *direction, double step); 01200 void SteepestDescent(int steps, double econv = 1e-6f, int method = OBFF_ANALYTICAL_GRADIENT); 01227 void SteepestDescentInitialize(int steps = 1000, double econv = 1e-6f, int method = OBFF_ANALYTICAL_GRADIENT); 01242 bool SteepestDescentTakeNSteps(int n); 01258 void ConjugateGradients(int steps, double econv = 1e-6f, int method = OBFF_ANALYTICAL_GRADIENT); 01286 void ConjugateGradientsInitialize(int steps = 1000, double econv = 1e-6f, int method = OBFF_ANALYTICAL_GRADIENT); 01302 bool ConjugateGradientsTakeNSteps(int n); 01304 01306 // Molecular Dynamics // 01308 01310 01311 01313 void GenerateVelocities(); 01333 void CorrectVelocities(); 01349 void MolecularDynamicsTakeNSteps(int n, double T, double timestep = 0.001, int method = OBFF_ANALYTICAL_GRADIENT); 01351 01353 // Constraints // 01355 01357 01358 01361 OBFFConstraints& GetConstraints(); 01365 void SetConstraints(OBFFConstraints& constraints); 01373 void SetFixAtom(int index); 01377 void UnsetFixAtom(); 01386 void SetIgnoreAtom(int index); 01390 void UnsetIgnoreAtom(); 01391 01393 static bool IgnoreCalculation(int a, int b); 01395 static bool IgnoreCalculation(int a, int b, int c); 01397 static bool IgnoreCalculation(int a, int b, int c, int d); 01399 01400 01402 // Validation // 01404 01406 01407 01408 bool DetectExplosion(); 01410 vector3 ValidateLineSearch(OBAtom *atom, vector3 &direction); 01412 void ValidateSteepestDescent(int steps); 01414 void ValidateConjugateGradients(int steps); 01416 virtual bool Validate() { return false; } 01421 virtual bool ValidateGradients() { return false; } 01426 vector3 ValidateGradientError(vector3 &numgrad, vector3 &anagrad); 01428 01430 // Vector Analysis // 01432 01434 01435 01443 static double VectorBondDerivative(double *pos_a, double *pos_b, 01444 double *force_a, double *force_b); 01448 static double VectorDistanceDerivative(const double* const pos_i, const double* const pos_j, 01449 double *force_i, double *force_j); 01451 static double VectorLengthDerivative(vector3 &a, vector3 &b); 01452 01463 static double VectorAngleDerivative(double *pos_a, double *pos_b, double *pos_c, 01464 double *force_a, double *force_b, double *force_c); 01466 static double VectorAngleDerivative(vector3 &a, vector3 &b, vector3 &c); 01479 static double VectorOOPDerivative(double *pos_a, double *pos_b, double *pos_c, double *pos_d, 01480 double *force_a, double *force_b, double *force_c, double *force_d); 01482 static double VectorOOPDerivative(vector3 &a, vector3 &b, vector3 &c, vector3 &d); 01494 static double VectorTorsionDerivative(double *pos_a, double *pos_b, double *pos_c, double *pos_d, 01495 double *force_a, double *force_b, double *force_c, double *force_d); 01497 static double VectorTorsionDerivative(vector3 &a, vector3 &b, vector3 &c, vector3 &d); 01498 01504 static void VectorSubtract(double *i, double *j, double *result) 01505 { 01506 for (unsigned int c = 0; c < 3; ++c) 01507 result[c] = i[c] - j[c]; 01508 } 01509 01510 static void VectorSubtract(const double* const i, const double* const j, double *result) 01511 { 01512 for (unsigned int c = 0; c < 3; ++c) 01513 result[c] = i[c] - j[c]; 01514 } 01515 01521 static void VectorAdd(double *i, double *j, double *result) 01522 { 01523 for (unsigned int c = 0; c < 3; ++c) 01524 result[c] = i[c] + j[c]; 01525 } 01526 01532 static void VectorDivide(double *i, double n, double *result) 01533 { 01534 for (unsigned int c = 0; c < 3; ++c) 01535 result[c] = i[c] / n; 01536 } 01537 01543 static void VectorMultiply(double *i, double n, double *result) 01544 { 01545 for (unsigned int c = 0; c < 3; ++c) 01546 result[c] = i[c] * n; 01547 } 01548 01549 static void VectorMultiply(const double* const i, const double n, double *result) 01550 { 01551 for (unsigned int c = 0; c < 3; ++c) 01552 result[c] = i[c] * n; 01553 } 01554 01559 static void VectorSelfMultiply(double *i, double n) 01560 { 01561 for (unsigned int c = 0; c < 3; ++c) 01562 i[c] *= n; 01563 } 01564 01568 static void VectorNormalize(double *i) 01569 { 01570 double length = VectorLength(i); 01571 for (unsigned int c = 0; c < 3; ++c) 01572 i[c] /= length; 01573 } 01574 01579 static void VectorCopy(double *from, double *to) 01580 { 01581 for (unsigned int c = 0; c < 3; ++c) 01582 to[c] = from[c]; 01583 } 01584 01589 static double VectorLength(double *i) 01590 { 01591 return sqrt( i[0]*i[0] + i[1]*i[1] + i[2]*i[2] ); 01592 } 01593 01594 static double VectorDistance(double *pos_i, double *pos_j) 01595 { 01596 double ij[3]; 01597 VectorSubtract(pos_i, pos_j, ij); 01598 const double rij = VectorLength(ij); 01599 return rij; 01600 } 01601 01608 static double VectorAngle(double *i, double *j, double *k); 01609 01617 static double VectorTorsion(double *i, double *j, double *k, double *l); 01618 01626 static double VectorOOP(double *i, double *j, double *k, double *l); 01627 01631 static void VectorClear(double *i) 01632 { 01633 for (unsigned int c = 0; c < 3; ++c) 01634 i[c] = 0.0; 01635 } 01636 01642 static double VectorDot(double *i, double *j) 01643 { 01644 double result = 0.0; 01645 // Written as a loop for vectorization 01646 // Loop will be unrolled by compiler otherwise 01647 for (unsigned int c = 0; c < 3; ++c) 01648 result += i[c]*j[c]; 01649 return result; 01650 } 01651 01657 static void VectorCross(double *i, double *j, double *result) 01658 { 01659 result[0] = i[1]*j[2] - i[2]*j[1]; 01660 result[1] = - i[0]*j[2] + i[2]*j[0]; 01661 result[2] = i[0]*j[1] - i[1]*j[0]; 01662 } 01663 01664 static void PrintVector(double *i) 01665 { 01666 std::cout << "<" << i[0] << ", " << i[1] << ", " << i[2] << ">" << std::endl; 01667 } 01669 01670 }; // class OBForceField 01671 01672 }// namespace OpenBabel 01673 01674 #endif // OB_FORCEFIELD_H 01675
This file is part of the documentation for Open Babel, version 2.2.0.
1.5.6.