--- trunk/OOPSE/libmdtools/Integrator.hpp 2003/09/16 20:02:11 767 +++ trunk/OOPSE/libmdtools/Integrator.hpp 2003/09/29 21:16:11 790 @@ -30,7 +30,7 @@ template class Integrator virtual double getConservedQuantity(void); protected: - + virtual void integrateStep( int calcPot, int calcStress ); virtual void preMove( void ); virtual void moveA( void ); @@ -44,6 +44,8 @@ template class Integrator virtual void calcForce( int calcPot, int calcStress ); virtual void thermalize(); + virtual void rotationPropagation( DirectionalAtom* dAtom, double ji[3] ); + void checkConstraints( void ); void rotate( int axes1, int axes2, double angle, double j[3], double A[3][3] ); @@ -135,20 +137,19 @@ template class NPTi : public T{ -template class NPTi : public T{ +template class NPT : public T{ public: - NPTi ( SimInfo *theInfo, ForceFields* the_ff); - virtual ~NPTi(); + NPT ( SimInfo *theInfo, ForceFields* the_ff); + virtual ~NPT(); virtual void integrateStep( int calcPot, int calcStress ){ calcStress = 1; T::integrateStep( calcPot, calcStress ); - /* accIntegralOfChidt(); */ } - virtual double getConservedQuantity(void); + virtual double getConservedQuantity(void) = 0; void setTauThermostat(double tt) {tauThermostat = tt; have_tau_thermostat=1;} void setTauBarostat(double tb) {tauBarostat = tb; have_tau_barostat=1;} @@ -167,15 +168,35 @@ template class NPTi : public T{ (protected virtual void resetIntegrator( void ); + virtual void getVelScaleA( double sc[3], double vel[3] ) = 0; + virtual void getVelScaleB( double sc[3], int index ) = 0; + virtual void getPosScale(double pos[3], double COM[3], + int index, double sc[3]) = 0; + + virtual bool chiConverged( void ); + virtual bool etaConverged( void ) = 0; + + virtual void evolveChiA( void ); + virtual void evolveEtaA( void ) = 0; + virtual void evolveChiB( void ); + virtual void evolveEtaB( void ) = 0; + + virtual void scaleSimBox( void ) = 0; + void accIntegralOfChidt(void) { integralOfChidt += dt * chi;} // chi and eta are the propagated degrees of freedom + double oldChi; + double prevChi; double chi; - double eta; double NkBT; double fkBT; + double tt2, tb2; + double instaTemp, instaPress, instaVol; + double press[3][3]; + int Nparticles; double integralOfChidt; @@ -204,64 +225,32 @@ template class NPTim : public T{ }; -template class NPTim : public T{ - +template class NPTi : public T{ + public: + NPTi( SimInfo *theInfo, ForceFields* the_ff); + ~NPTi(); - NPTim ( SimInfo *theInfo, ForceFields* the_ff); - virtual ~NPTim() {} - - virtual void integrateStep( int calcPot, int calcStress ){ - calcStress = 1; - T::integrateStep( calcPot, calcStress ); - accIntegralOfChidt(); - } - virtual double getConservedQuantity(void); + virtual void resetIntegrator(void); - void setTauThermostat(double tt) {tauThermostat = tt; have_tau_thermostat=1;} - void setTauBarostat(double tb) {tauBarostat = tb; have_tau_barostat=1;} - void setTargetTemp(double tt) {targetTemp = tt; have_target_temp = 1;} - void setTargetPressure(double tp) {targetPressure = tp; have_target_pressure = 1;} - void setChiTolerance(double tol) {chiTolerance = tol;} - void setPosIterTolerance(double tol) {posIterTolerance = tol;} - protected: - virtual void moveA( void ); - virtual void moveB( void ); + - virtual int readyCheck(); + virtual void evolveEtaA(void); + virtual void evolveEtaB(void); - virtual void resetIntegrator( void ); + virtual bool etaConverged( void ); - void accIntegralOfChidt(void) { integralOfChidt += dt * chi;} - - Molecule* myMolecules; - Atom** myAtoms; + virtual void scaleSimBox( void ); - // chi and eta are the propagated degrees of freedom - - double chi; - double eta; - double NkBT; - double integralOfChidt; + virtual void getVelScaleA( double sc[3], double vel[3] ); + virtual void getVelScaleB( double sc[3], int index ); + virtual void getPosScale(double pos[3], double COM[3], + int index, double sc[3]); - // targetTemp, targetPressure, and tauBarostat must be set. - // One of qmass or tauThermostat must be set; - - double targetTemp; - double targetPressure; - double tauThermostat; - double tauBarostat; - - short int have_tau_thermostat, have_tau_barostat, have_target_temp; - short int have_target_pressure; - double chiTolerance; - short int have_chi_tolerance; - double posIterTolerance; - short int have_pos_iter_tolerance; - + double eta, oldEta, prevEta; }; template class NPTzm : public T{ @@ -320,61 +309,26 @@ template class NPTf : public T{ (public) NPTf ( SimInfo *theInfo, ForceFields* the_ff); virtual ~NPTf(); - virtual void integrateStep( int calcPot, int calcStress ){ - calcStress = 1; - T::integrateStep( calcPot, calcStress ); - } - virtual double getConservedQuantity(void); + virtual void resetIntegrator(void); - void setTauThermostat(double tt) {tauThermostat = tt; have_tau_thermostat=1;} - void setTauBarostat(double tb) {tauBarostat = tb; have_tau_barostat=1;} - void setTargetTemp(double tt) {targetTemp = tt; have_target_temp = 1;} - void setTargetPressure(double tp) {targetPressure = tp; have_target_pressure = 1;} - void setChiTolerance(double tol) {chiTolerance = tol;} - void setPosIterTolerance(double tol) {posIterTolerance = tol;} - protected: - virtual void moveA( void ); - virtual void moveB( void ); + virtual void evolveEtaA(void); + virtual void evolveEtaB(void); - virtual void resetIntegrator( void ); + virtual bool etaConverged( void ); - virtual int readyCheck(); + virtual void scaleSimBox( void ); + virtual void getVelScaleA( double sc[3], double vel[3] ); + virtual void getVelScaleB( double sc[3], int index ); + virtual void getPosScale(double pos[3], double COM[3], + int index, double sc[3]); - // chi and eta are the propagated degrees of freedom - - double chi; double eta[3][3]; - double NkBT; - double fkBT; - - int Nparticles; - - double *oldPos; - double *oldVel; - double *oldJi; - - double integralOfChidt; - - // targetTemp, targetPressure, and tauBarostat must be set. - // One of qmass or tauThermostat must be set; - - double targetTemp; - double targetPressure; - double tauThermostat; - double tauBarostat; - - short int have_tau_thermostat, have_tau_barostat, have_target_temp; - short int have_target_pressure; - double chiTolerance; - short int have_chi_tolerance; - double posIterTolerance; - short int have_pos_iter_tolerance; - double etaTolerance; - short int have_eta_tolerance; + double oldEta[3][3]; + double prevEta[3][3]; }; template class NPTxym : public T{ @@ -427,120 +381,6 @@ template class NPTxym : public T{ (protect }; - -template class NPTfm : public T{ - -public: - - NPTfm ( SimInfo *theInfo, ForceFields* the_ff); - virtual ~NPTfm() {}; - - virtual void integrateStep( int calcPot, int calcStress ){ - calcStress = 1; - T::integrateStep( calcPot, calcStress ); - accIntegralOfChidt(); - } - - virtual double getConservedQuantity(void); - - void setTauThermostat(double tt) {tauThermostat = tt; have_tau_thermostat=1;} - void setTauBarostat(double tb) {tauBarostat = tb; have_tau_barostat=1;} - void setTargetTemp(double tt) {targetTemp = tt; have_target_temp = 1;} - void setTargetPressure(double tp) {targetPressure = tp; have_target_pressure = 1;} - void setChiTolerance(double tol) {chiTolerance = tol;} - void setPosIterTolerance(double tol) {posIterTolerance = tol;} - -protected: - - virtual void moveA( void ); - virtual void moveB( void ); - - virtual void resetIntegrator( void ); - - virtual int readyCheck(); - - void accIntegralOfChidt(void) { integralOfChidt += dt * chi;} - - Molecule* myMolecules; - Atom** myAtoms; - - // chi and eta are the propagated degrees of freedom - - double chi; - double eta[3][3]; - double NkBT; - double integralOfChidt; - - // targetTemp, targetPressure, and tauBarostat must be set. - // One of qmass or tauThermostat must be set; - - double targetTemp; - double targetPressure; - double tauThermostat; - double tauBarostat; - - short int have_tau_thermostat, have_tau_barostat, have_target_temp; - short int have_target_pressure; - double chiTolerance; - short int have_chi_tolerance; - double posIterTolerance; - short int have_pos_iter_tolerance; - -}; - - -template class NPTpr : public T{ - -public: - - NPTpr ( SimInfo *theInfo, ForceFields* the_ff); - virtual ~NPTpr() {}; - - virtual void integrateStep( int calcPot, int calcStress ){ - calcStress = 1; - T::integrateStep( calcPot, calcStress ); - } - - void setTauThermostat(double tt) {tauThermostat = tt; have_tau_thermostat=1;} - void setTauBarostat(double tb) {tauBarostat = tb; have_tau_barostat=1;} - void setTargetTemp(double tt) {targetTemp = tt; have_target_temp = 1;} - void setTargetPressure(double tp) {targetPressure = tp; have_target_pressure = 1;} - void setChiTolerance(double tol) {chiTolerance = tol;} - void setPosIterTolerance(double tol) {posIterTolerance = tol;} - -protected: - - virtual void moveA( void ); - virtual void moveB( void ); - - virtual int readyCheck(); - - virtual void resetIntegrator( void ); - - // chi and eta are the propagated degrees of freedom - - double chi; - double eta[3][3]; - double NkBT; - - // targetTemp, targetPressure, and tauBarostat must be set. - // One of qmass or tauThermostat must be set; - - double targetTemp; - double targetPressure; - double tauThermostat; - double tauBarostat; - - short int have_tau_thermostat, have_tau_barostat, have_target_temp; - short int have_target_pressure; - double chiTolerance; - short int have_chi_tolerance; - double posIterTolerance; - short int have_pos_iter_tolerance; - -}; - - template class ZConstraint : public T { public: @@ -555,7 +395,7 @@ template class ZConstraint : public T { virtual double getHFOfUnconsMols(Atom* atom, double totalForce) = 0; protected: - ZConstraint* zconsIntegrator;; + ZConstraint* zconsIntegrator; }; class PolicyByNumber : public ForceSubtractionPolicy{ @@ -608,9 +448,20 @@ template class ZConstraint : public T { virtual void update(); //which is called to indicate the molecules' migration #endif + enum ZConsState {zcsMoving, zcsFixed}; + + vector zconsMols; //z-constraint molecules array + vector states; //state of z-constraint molecules + + + + int totNumOfUnconsAtoms; //total number of uncontraint atoms + double totalMassOfUncons; //total mas of unconstraint molecules + + protected: - enum ZConsState {zcsMoving, zcsFixed}; + virtual void calcForce( int calcPot, int calcStress ); virtual void thermalize(void); @@ -633,17 +484,17 @@ template class ZConstraint : public T { double zForceConst; //base force constant term //which is estimate by OOPSE - vector zconsMols; //z-constraint molecules array + vector massOfZConsMols; //mass of z-constraint molecule vector kz; //force constant array - vector states; //state of z-constraint molecules + vector zPos; // vector unconsMols; //unconstraint molecules array vector massOfUnconsMols; //mass array of unconstraint molecules - double totalMassOfUncons; //total mas of unconstraint molecules + vector* parameters; // vector indexOfAllZConsMols; //index of All Z-Constraint Molecuels @@ -652,8 +503,8 @@ template class ZConstraint : public T { double* fz; double* curZPos; - int totNumOfUnconsAtoms; //total number of uncontraint atoms + int whichDirection; //constraint direction private: