# | Line 4 | Line 4 | |
---|---|---|
4 | #include <string> | |
5 | #include <vector> | |
6 | #include "Atom.hpp" | |
7 | + | #include "StuntDouble.hpp" |
8 | #include "Molecule.hpp" | |
9 | #include "SRI.hpp" | |
10 | #include "AbstractClasses.hpp" | |
# | Line 12 | Line 13 | |
13 | #include "Thermo.hpp" | |
14 | #include "ReadWrite.hpp" | |
15 | #include "ZConsWriter.hpp" | |
16 | + | #include "Restraints.hpp" |
17 | ||
18 | using namespace std; | |
19 | const double kB = 8.31451e-7;// boltzmann constant amu*Ang^2*fs^-2/K | |
# | Line 20 | Line 22 | const double tol = 1.0e-6; | |
22 | const int maxIteration = 300; | |
23 | const double tol = 1.0e-6; | |
24 | ||
25 | < | |
25 | > | class VelVerletConsFramework; |
26 | template<typename T = BaseIntegrator> class Integrator : public T { | |
27 | ||
28 | public: | |
# | Line 33 | Line 35 | template<typename T = BaseIntegrator> class Integrator | |
35 | protected: | |
36 | ||
37 | virtual void integrateStep( int calcPot, int calcStress ); | |
38 | < | virtual void preMove( void ); |
38 | > | //virtual void preMove( void ); |
39 | virtual void moveA( void ); | |
40 | virtual void moveB( void ); | |
41 | < | virtual void constrainA( void ); |
42 | < | virtual void constrainB( void ); |
41 | > | //virtual void constrainA( void ); |
42 | > | //virtual void constrainB( void ); |
43 | virtual int readyCheck( void ) { return 1; } | |
44 | ||
45 | virtual void resetIntegrator( void ) { } | |
# | Line 45 | Line 47 | template<typename T = BaseIntegrator> class Integrator | |
47 | virtual void calcForce( int calcPot, int calcStress ); | |
48 | virtual void thermalize(); | |
49 | ||
50 | < | virtual void rotationPropagation( DirectionalAtom* dAtom, double ji[3] ); |
50 | > | virtual bool stopIntegrator() {return false;} |
51 | ||
52 | < | void checkConstraints( void ); |
52 | > | virtual void rotationPropagation( StuntDouble* sd, double ji[3] ); |
53 | > | |
54 | > | //void checkConstraints( void ); |
55 | void rotate( int axes1, int axes2, double angle, double j[3], | |
56 | double A[3][3] ); | |
57 | ||
58 | ForceFields* myFF; | |
59 | ||
60 | SimInfo *info; // all the info we'll ever need | |
61 | + | vector<StuntDouble*> integrableObjects; |
62 | int nAtoms; /* the number of atoms */ | |
63 | int oldAtoms; | |
64 | Atom **atoms; /* array of atom pointers */ | |
65 | Molecule* molecules; | |
66 | int nMols; | |
67 | ||
68 | < | int isConstrained; // boolean to know whether the systems contains |
64 | < | // constraints. |
65 | < | int nConstrained; // counter for number of constraints |
66 | < | int *constrainedA; // the i of a constraint pair |
67 | < | int *constrainedB; // the j of a constraint pair |
68 | < | double *constrainedDsqr; // the square of the constraint distance |
68 | > | VelVerletConsFramework* consFramework; |
69 | ||
70 | < | int* moving; // tells whether we are moving atom i |
71 | < | int* moved; // tells whether we have moved atom i |
72 | < | double* oldPos; // pre constrained positions |
70 | > | //int isConstrained; // boolean to know whether the systems contains constraints. |
71 | > | //int nConstrained; // counter for number of constraints |
72 | > | //int *constrainedA; // the i of a constraint pair |
73 | > | //int *constrainedB; // the j of a constraint pair |
74 | > | //double *constrainedDsqr; // the square of the constraint distance |
75 | > | |
76 | > | //int* moving; // tells whether we are moving atom i |
77 | > | //int* moved; // tells whether we have moved atom i |
78 | > | //double* oldPos; // pre constrained positions |
79 | ||
80 | short isFirst; /*boolean for the first time integrate is called */ | |
81 | ||
# | Line 83 | Line 89 | typedef Integrator<BaseIntegrator> RealIntegrator; | |
89 | }; | |
90 | ||
91 | typedef Integrator<BaseIntegrator> RealIntegrator; | |
92 | + | |
93 | + | // ansi instantiation |
94 | + | // template class Integrator<BaseIntegrator>; |
95 | + | |
96 | ||
97 | template<typename T> class NVE : public T { | |
98 | ||
# | Line 418 | Line 428 | template<typename T> class ZConstraint : public T { (p | |
428 | ||
429 | void zeroOutVel(); | |
430 | void doZconstraintForce(); | |
431 | < | void doHarmonic(); |
431 | > | void doHarmonic(vector<double>& resPos); |
432 | bool checkZConsState(); | |
433 | ||
434 | bool haveFixedZMols(); | |
# | Line 449 | Line 459 | template<typename T> class ZConstraint : public T { (p | |
459 | ||
460 | vector<int> indexOfAllZConsMols; //index of All Z-Constraint Molecuels | |
461 | ||
462 | < | int* indexOfZConsMols; //index of local Z-Constraint Molecules |
463 | < | double* fz; |
464 | < | double* curZPos; |
455 | < | |
462 | > | vector<int> indexOfZConsMols; //index of local Z-Constraint Molecules |
463 | > | vector<double> fz; |
464 | > | vector<double> curZPos; |
465 | ||
466 | + | bool usingSMD; |
467 | + | vector<double> prevCantPos; |
468 | + | vector<double> cantPos; |
469 | + | vector<double> cantVel; |
470 | ||
471 | + | double zconsFixTime; |
472 | + | double zconsGap; |
473 | + | bool hasZConsGap; |
474 | + | vector<double> endFixTime; |
475 | + | |
476 | int whichDirection; //constraint direction | |
477 | ||
478 | private: | |
# | Line 467 | Line 485 | template<typename T> class ZConstraint : public T { (p | |
485 | double calcMovingMolsCOMVel(); | |
486 | double calcSysCOMVel(); | |
487 | double calcTotalForce(); | |
488 | < | |
488 | > | void updateZPos(); |
489 | > | void updateCantPos(); |
490 | > | |
491 | ForceSubtractionPolicy* forcePolicy; //force subtraction policy | |
492 | friend class ForceSubtractionPolicy; | |
493 | ||
494 | }; | |
495 | ||
476 | – | class OOPSEMinimizerBase : public RealIntegrator { |
477 | – | public: |
496 | ||
497 | < | OOPSEMinimizerBase ( SimInfo *theInfo, ForceFields* the_ff ); |
498 | < | virtual ~OOPSEMinimizerBase(); |
497 | > | //Sympletic quaternion Scheme Integrator |
498 | > | //Reference: |
499 | > | // T.F. Miller, M. Eleftheriou, P. Pattnaik, A. Ndirango, D. Newns and G.J. Martyna |
500 | > | //Symplectic quaternion Scheme for biophysical molecular dynamics |
501 | > | //116(20), 8649, J. Chem. Phys. (2002) |
502 | > | template<typename T> class SQSIntegrator : public T{ |
503 | > | public: |
504 | > | virtual void moveA(); |
505 | > | virtual void moveB(); |
506 | > | protected: |
507 | > | void freeRotor(); |
508 | > | void rotate(int k, double dt); |
509 | ||
482 | – | double calcGradient(vector<double>& x, vector<double>& grad); |
483 | – | void setOptCoor(vector<double>& x); |
484 | – | void getOptCoor(vector<double>& x); |
485 | – | void output(); |
486 | – | |
510 | }; | |
488 | – | |
489 | – | template<typename TMinimizer> class OOPSEMinimizer : public OOPSEMinimizerBase, TMinimizer{ |
490 | – | public: |
491 | – | void writeOutput(); |
492 | – | }; |
493 | – | |
511 | #endif |
– | Removed lines |
+ | Added lines |
< | Changed lines |
> | Changed lines |