OOPSE/libmdtools/Integrator.hpp

#ifndef _INTEGRATOR_H_
#define _INTEGRATOR_H_

#include "Atom.hpp"
#include "SRI.hpp"
#include "AbstractClasses.hpp"
#include "SimInfo.hpp"
#include "ForceFields.hpp"
#include "Thermo.hpp"
#include "ReadWrite.hpp"

const double kB = 8.31451e-7;// boltzmann constant amu*Ang^2*fs^-2/K
const double eConvert = 4.184e-4; // converts kcal/mol -> amu*A^2/fs^2
const double p_convert = 1.63882576e8; //converts amu*fs^-2*Ang^-1 -> atm
const int maxIteration = 300;
const double tol = 1.0e-6;

class Integrator : public BaseIntegrator {

public:
  Integrator( SimInfo *theInfo, ForceFields* the_ff );
  virtual ~Integrator();
  void integrate( void );


protected:
  
  virtual void integrateStep( int calcPot, int calcStress );
  virtual void preMove( void );
  virtual void moveA( void );
  virtual void moveB( void );
  virtual void constrainA( void );
  virtual void constrainB( void );
  virtual int  readyCheck( void ) { return 1; }
  
  void checkConstraints( void );
  void rotate( int axes1, int axes2, double angle, double j[3], 
               double A[3][3] );


  ForceFields* myFF;

  SimInfo *info; // all the info we'll ever need
  int nAtoms;  /* the number of atoms */
  int oldAtoms;
  Atom **atoms; /* array of atom pointers */
  Molecule* molecules;
  int nMols;

  int isConstrained; // boolean to know whether the systems contains
                     // constraints.
  int nConstrained;  // counter for number of constraints 
  int *constrainedA; // the i of a constraint pair 
  int *constrainedB; // the j of a constraint pair 
  double *constrainedDsqr; // the square of the constraint distance 
  
  int* moving; // tells whether we are moving atom i
  int* moved;  // tells whether we have moved atom i
  double* oldPos; // pre constrained positions 

  short isFirst; /*boolean for the first time integrate is called */
  
  double dt;
  double dt2;

  Thermo *tStats;
  StatWriter*  statOut;
  DumpWriter*  dumpOut;
  
};

class NVE : public Integrator{

public:
  NVE ( SimInfo *theInfo, ForceFields* the_ff ):
    Integrator( theInfo, the_ff ){}
  virtual ~NVE(){}

  
};

class NVT : public Integrator{

public:

  NVT ( SimInfo *theInfo, ForceFields* the_ff);
  virtual ~NVT() {}

  void setTauThermostat(double tt) {tauThermostat = tt; have_tau_thermostat=1;}
  void setTargetTemp(double tt) {targetTemp = tt; have_target_temp = 1;}

protected:

  virtual void moveA( void );
  virtual void moveB( void );

  virtual int readyCheck();

  // chi is a propagated degree of freedom.

  double chi;

  // targetTemp must be set.  tauThermostat must also be set;

  double targetTemp;
  double tauThermostat;
  
  short int have_tau_thermostat, have_target_temp;

};


class NPTi : public Integrator{

public:

  NPTi ( SimInfo *theInfo, ForceFields* the_ff);
  virtual ~NPTi() {};

  virtual void integrateStep( int calcPot, int calcStress ){
    calcStress = 1;
    Integrator::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;}

protected:

  virtual void  moveA( void );
  virtual void moveB( void );

  virtual int readyCheck();

  // chi and eta are the propagated degrees of freedom

  double chi;
  double eta;
  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;

};

class NPTim : public Integrator{

public:

  NPTim ( SimInfo *theInfo, ForceFields* the_ff);
  virtual ~NPTim() {};

  virtual void integrateStep( int calcPot, int calcStress ){
    calcStress = 1;
    Integrator::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;}

protected:

  virtual void  moveA( void );
  virtual void moveB( void );

  virtual int readyCheck();

  // chi and eta are the propagated degrees of freedom

  double chi;
  double eta;
  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;

};

class NPTf : public Integrator{

public:

  NPTf ( SimInfo *theInfo, ForceFields* the_ff);
  virtual ~NPTf() {};

  virtual void integrateStep( int calcPot, int calcStress ){
    calcStress = 1;
    Integrator::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;}

protected:

  virtual void  moveA( void );
  virtual void moveB( void );

  virtual int readyCheck();

  // 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;

};

class NPTfm : public Integrator{

public:

  NPTfm ( SimInfo *theInfo, ForceFields* the_ff);
  virtual ~NPTfm() {};

  virtual void integrateStep( int calcPot, int calcStress ){
    calcStress = 1;
    Integrator::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;}

protected:

  virtual void  moveA( void );
  virtual void moveB( void );

  virtual int readyCheck();

  // 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;

};

#endif
Revision:	600
Committed:	Mon Jul 14 22:38:13 2003 UTC (20 years, 11 months ago) by gezelter
File size:	7038 byte(s)
Log Message:	Fixes for get and set routines in Atom and DirectionalAtom
#	User	Rev	Content
1	mmeineke	377	#ifndef _INTEGRATOR_H_
2			#define _INTEGRATOR_H_
3
4			#include "Atom.hpp"
5			#include "SRI.hpp"
6			#include "AbstractClasses.hpp"
7			#include "SimInfo.hpp"
8			#include "ForceFields.hpp"
9	mmeineke	540	#include "Thermo.hpp"
10			#include "ReadWrite.hpp"
11	mmeineke	377
12	mmeineke	561	const double kB = 8.31451e-7;// boltzmann constant amuAng^2fs^-2/K
13			const double eConvert = 4.184e-4; // converts kcal/mol -> amu*A^2/fs^2
14	mmeineke	586	const double p_convert = 1.63882576e8; //converts amufs^-2Ang^-1 -> atm
15	mmeineke	561	const int maxIteration = 300;
16			const double tol = 1.0e-6;
17
18	mmeineke	540	class Integrator : public BaseIntegrator {
19	mmeineke	377
20			public:
21	mmeineke	561	Integrator( SimInfo theInfo, ForceFields the_ff );
22	mmeineke	548	virtual ~Integrator();
23	mmeineke	377	void integrate( void );
24
25
26	mmeineke	540	protected:
27	mmeineke	377
28	mmeineke	540	virtual void integrateStep( int calcPot, int calcStress );
29	mmeineke	548	virtual void preMove( void );
30	mmeineke	540	virtual void moveA( void );
31			virtual void moveB( void );
32			virtual void constrainA( void );
33			virtual void constrainB( void );
34	mmeineke	559	virtual int readyCheck( void ) { return 1; }
35	mmeineke	377
36	mmeineke	540	void checkConstraints( void );
37	mmeineke	377	void rotate( int axes1, int axes2, double angle, double j[3],
38	gezelter	600	double A[3][3] );
39	mmeineke	377
40	mmeineke	540
41	mmeineke	377	ForceFields* myFF;
42
43	mmeineke	540	SimInfo *info; // all the info we'll ever need
44			int nAtoms; /* the number of atoms */
45	mmeineke	548	int oldAtoms;
46	mmeineke	540	Atom *atoms; / array of atom pointers */
47	mmeineke	423	Molecule* molecules;
48			int nMols;
49
50	mmeineke	548	int isConstrained; // boolean to know whether the systems contains
51			// constraints.
52			int nConstrained; // counter for number of constraints
53			int *constrainedA; // the i of a constraint pair
54			int *constrainedB; // the j of a constraint pair
55			double *constrainedDsqr; // the square of the constraint distance
56
57			int* moving; // tells whether we are moving atom i
58			int* moved; // tells whether we have moved atom i
59	mmeineke	561	double* oldPos; // pre constrained positions
60	mmeineke	548
61	mmeineke	540	short isFirst; /boolean for the first time integrate is called /
62
63			double dt;
64	mmeineke	541	double dt2;
65	gezelter	560
66	mmeineke	540	Thermo *tStats;
67			StatWriter* statOut;
68			DumpWriter* dumpOut;
69	mmeineke	377
70			};
71
72	mmeineke	555	class NVE : public Integrator{
73	mmeineke	540
74	mmeineke	561	public:
75			NVE ( SimInfo theInfo, ForceFields the_ff ):
76	mmeineke	555	Integrator( theInfo, the_ff ){}
77			virtual ~NVE(){}
78
79
80
81			};
82
83	mmeineke	540	class NVT : public Integrator{
84
85	gezelter	560	public:
86
87	mmeineke	561	NVT ( SimInfo theInfo, ForceFields the_ff);
88			virtual ~NVT() {}
89	mmeineke	540
90	gezelter	560	void setTauThermostat(double tt) {tauThermostat = tt; have_tau_thermostat=1;}
91			void setTargetTemp(double tt) {targetTemp = tt; have_target_temp = 1;}
92
93	mmeineke	540	protected:
94	gezelter	560
95	mmeineke	561	virtual void moveA( void );
96			virtual void moveB( void );
97	mmeineke	540
98	mmeineke	561	virtual int readyCheck();
99	gezelter	560
100	gezelter	565	// chi is a propagated degree of freedom.
101	gezelter	560
102	gezelter	565	double chi;
103	gezelter	560
104	gezelter	565	// targetTemp must be set. tauThermostat must also be set;
105	gezelter	560
106			double targetTemp;
107			double tauThermostat;
108	mmeineke	561
109	gezelter	565	short int have_tau_thermostat, have_target_temp;
110	gezelter	560
111	mmeineke	540	};
112
113
114	gezelter	574	class NPTi : public Integrator{
115	mmeineke	540
116	gezelter	560	public:
117
118	gezelter	574	NPTi ( SimInfo theInfo, ForceFields the_ff);
119			virtual ~NPTi() {};
120	gezelter	560
121	mmeineke	594	virtual void integrateStep( int calcPot, int calcStress ){
122			calcStress = 1;
123			Integrator::integrateStep( calcPot, calcStress );
124			}
125
126	gezelter	560	void setTauThermostat(double tt) {tauThermostat = tt; have_tau_thermostat=1;}
127			void setTauBarostat(double tb) {tauBarostat = tb; have_tau_barostat=1;}
128			void setTargetTemp(double tt) {targetTemp = tt; have_target_temp = 1;}
129			void setTargetPressure(double tp) {targetPressure = tp; have_target_pressure = 1;}
130
131			protected:
132
133	mmeineke	561	virtual void moveA( void );
134			virtual void moveB( void );
135	gezelter	560
136	mmeineke	561	virtual int readyCheck();
137	gezelter	560
138	gezelter	565	// chi and eta are the propagated degrees of freedom
139	gezelter	560
140	gezelter	565	double chi;
141			double eta;
142	gezelter	574	double NkBT;
143	gezelter	560
144			// targetTemp, targetPressure, and tauBarostat must be set.
145			// One of qmass or tauThermostat must be set;
146
147			double targetTemp;
148			double targetPressure;
149			double tauThermostat;
150			double tauBarostat;
151
152			short int have_tau_thermostat, have_tau_barostat, have_target_temp;
153	gezelter	565	short int have_target_pressure;
154	gezelter	560
155			};
156
157	gezelter	596	class NPTim : public Integrator{
158
159			public:
160
161			NPTim ( SimInfo theInfo, ForceFields the_ff);
162			virtual ~NPTim() {};
163
164			virtual void integrateStep( int calcPot, int calcStress ){
165			calcStress = 1;
166			Integrator::integrateStep( calcPot, calcStress );
167			}
168
169			void setTauThermostat(double tt) {tauThermostat = tt; have_tau_thermostat=1;}
170			void setTauBarostat(double tb) {tauBarostat = tb; have_tau_barostat=1;}
171			void setTargetTemp(double tt) {targetTemp = tt; have_target_temp = 1;}
172			void setTargetPressure(double tp) {targetPressure = tp; have_target_pressure = 1;}
173
174			protected:
175
176			virtual void moveA( void );
177			virtual void moveB( void );
178
179			virtual int readyCheck();
180
181			// chi and eta are the propagated degrees of freedom
182
183			double chi;
184			double eta;
185			double NkBT;
186
187			// targetTemp, targetPressure, and tauBarostat must be set.
188			// One of qmass or tauThermostat must be set;
189
190			double targetTemp;
191			double targetPressure;
192			double tauThermostat;
193			double tauBarostat;
194
195			short int have_tau_thermostat, have_tau_barostat, have_target_temp;
196			short int have_target_pressure;
197
198			};
199
200	gezelter	576	class NPTf : public Integrator{
201
202			public:
203
204			NPTf ( SimInfo theInfo, ForceFields the_ff);
205			virtual ~NPTf() {};
206
207	mmeineke	594	virtual void integrateStep( int calcPot, int calcStress ){
208			calcStress = 1;
209			Integrator::integrateStep( calcPot, calcStress );
210			}
211
212	gezelter	576	void setTauThermostat(double tt) {tauThermostat = tt; have_tau_thermostat=1;}
213			void setTauBarostat(double tb) {tauBarostat = tb; have_tau_barostat=1;}
214			void setTargetTemp(double tt) {targetTemp = tt; have_target_temp = 1;}
215			void setTargetPressure(double tp) {targetPressure = tp; have_target_pressure = 1;}
216
217			protected:
218
219			virtual void moveA( void );
220			virtual void moveB( void );
221
222			virtual int readyCheck();
223
224			// chi and eta are the propagated degrees of freedom
225
226			double chi;
227	gezelter	588	double eta[3][3];
228	gezelter	576	double NkBT;
229
230			// targetTemp, targetPressure, and tauBarostat must be set.
231			// One of qmass or tauThermostat must be set;
232
233			double targetTemp;
234			double targetPressure;
235			double tauThermostat;
236			double tauBarostat;
237
238			short int have_tau_thermostat, have_tau_barostat, have_target_temp;
239			short int have_target_pressure;
240
241			};
242
243	gezelter	596	class NPTfm : public Integrator{
244
245			public:
246
247			NPTfm ( SimInfo theInfo, ForceFields the_ff);
248			virtual ~NPTfm() {};
249
250			virtual void integrateStep( int calcPot, int calcStress ){
251			calcStress = 1;
252			Integrator::integrateStep( calcPot, calcStress );
253			}
254
255			void setTauThermostat(double tt) {tauThermostat = tt; have_tau_thermostat=1;}
256			void setTauBarostat(double tb) {tauBarostat = tb; have_tau_barostat=1;}
257			void setTargetTemp(double tt) {targetTemp = tt; have_target_temp = 1;}
258			void setTargetPressure(double tp) {targetPressure = tp; have_target_pressure = 1;}
259
260			protected:
261
262			virtual void moveA( void );
263			virtual void moveB( void );
264
265			virtual int readyCheck();
266
267			// chi and eta are the propagated degrees of freedom
268
269			double chi;
270			double eta[3][3];
271			double NkBT;
272
273			// targetTemp, targetPressure, and tauBarostat must be set.
274			// One of qmass or tauThermostat must be set;
275
276			double targetTemp;
277			double targetPressure;
278			double tauThermostat;
279			double tauBarostat;
280
281			short int have_tau_thermostat, have_tau_barostat, have_target_temp;
282			short int have_target_pressure;
283
284			};
285
286	mmeineke	377	#endif