ViewVC Help

View File | Revision Log | Show Annotations | View Changeset | Root Listing

root/OpenMD/branches/development/src/integrators/VelocityVerletIntegrator.cpp

Comparing:
trunk/src/integrators/VelocityVerletIntegrator.cpp (file contents), Revision 471 by chrisfen, Tue Apr 12 14:24:29 2005 UTC vs.
branches/development/src/integrators/VelocityVerletIntegrator.cpp (file contents), Revision 1723 by gezelter, Thu May 24 20:59:54 2012 UTC

#	Line 1 | Line 1
1	<	/*
1	>	/*
2		* Copyright (c) 2005 The University of Notre Dame. All Rights Reserved.
3		*
4		* The University of Notre Dame grants you ("Licensee") a
#	Line 6 | Line 6
6		* redistribute this software in source and binary code form, provided
7		* that the following conditions are met:
8		*
9	<	* 1. Acknowledgement of the program authors must be made in any
10	<	*    publication of scientific results based in part on use of the
11	<	*    program.  An acceptable form of acknowledgement is citation of
12	<	*    the article in which the program was described (Matthew
13	<	*    A. Meineke, Charles F. Vardeman II, Teng Lin, Christopher
14	<	*    J. Fennell and J. Daniel Gezelter, "OOPSE: An Object-Oriented
15	<	*    Parallel Simulation Engine for Molecular Dynamics,"
16	<	*    J. Comput. Chem. 26, pp. 252-271 (2005))
17	<	*
18	<	* 2. Redistributions of source code must retain the above copyright
9	>	* 1. Redistributions of source code must retain the above copyright
10		*    notice, this list of conditions and the following disclaimer.
11		*
12	<	* 3. Redistributions in binary form must reproduce the above copyright
12	>	* 2. Redistributions in binary form must reproduce the above copyright
13		*    notice, this list of conditions and the following disclaimer in the
14		*    documentation and/or other materials provided with the
15		*    distribution.
#	Line 37 | Line 28
28		* arising out of the use of or inability to use software, even if the
29		* University of Notre Dame has been advised of the possibility of
30		* such damages.
31	+	*
32	+	* SUPPORT OPEN SCIENCE!  If you use OpenMD or its source code in your
33	+	* research, please cite the appropriate papers when you publish your
34	+	* work.  Good starting points are:
35	+	*
36	+	* [1]  Meineke, et al., J. Comp. Chem. 26, 252-271 (2005).
37	+	* [2]  Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006).
38	+	* [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 24107 (2008).
39	+	* [4]  Kuang & Gezelter,  J. Chem. Phys. 133, 164101 (2010).
40	+	* [5]  Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011).
41		*/
42
43		/**
#	Line 50 | Line 51
51		#include "integrators/VelocityVerletIntegrator.hpp"
52		#include "integrators/DLM.hpp"
53		#include "utils/StringUtils.hpp"
54	+	#include "utils/ProgressBar.hpp"
55
56	<	namespace oopse {
57	<	VelocityVerletIntegrator::VelocityVerletIntegrator(SimInfo *info) : Integrator(info), rotAlgo(NULL) {
56	>	namespace OpenMD {
57	>	VelocityVerletIntegrator::VelocityVerletIntegrator(SimInfo *info) : Integrator(info) {
58		dt2 = 0.5 * dt;
57	–	rotAlgo = new DLM();
58	–	rattle = new Rattle(info);
59		}
60
61		VelocityVerletIntegrator::~VelocityVerletIntegrator() {
62	–	delete rotAlgo;
63	–	delete rattle;
62		}
63
64		void VelocityVerletIntegrator::initialize(){
65
66	<	forceMan_->init();
66	>	forceMan_->initialize();
67
68	<	// remove center of mass drift velocity (in case we passed in a configuration
69	<	// that was drifting
68	>	// remove center of mass drift velocity (in case we passed in a
69	>	// configuration that was drifting)
70		velocitizer_->removeComDrift();
71	+
72	+	// find the initial fluctuating charges.
73	+	flucQ_->initialize();
74
75		// initialize the forces before the first step
76	<	calcForce(true, true);
76	>	calcForce();
77
78	<	//execute constraint algorithm to make sure at the very beginning the system is constrained
78	>	// execute the constraint algorithm to make sure that the system is
79	>	// constrained at the very beginning
80		if (info_->getNGlobalConstraints() > 0) {
81	<	rattle->constraintA();
82	<	calcForce(true, true);
83	<	rattle->constraintB();
84	<	info_->getSnapshotManager()->advance();//copy the current snapshot to previous snapshot
81	>	rattle_->constraintA();
82	>	calcForce();
83	>	rattle_->constraintB();
84	>	//copy the current snapshot to previous snapshot
85	>	info_->getSnapshotManager()->advance();
86		}
87
88		if (needVelocityScaling) {
#	Line 89 | Line 92 | namespace oopse {
92		dumpWriter = createDumpWriter();
93
94		statWriter = createStatWriter();
95	<
93	<	if (simParams->getUseSolidThermInt()) {
94	<	restWriter = createRestWriter();
95	<	restWriter->writeZangle();
96	<	}
97	<
95	>
96		dumpWriter->writeDumpAndEor();
97	<
98	<
97	>
98	>	progressBar = new ProgressBar();
99	>
100		//save statistics, before writeStat,  we must save statistics
101		thermo.saveStat();
102		saveConservedQuantity();
103	+	if (simParams->getUseRNEMD())
104	+	rnemd_->getStarted();
105	+
106		statWriter->writeStat(currentSnapshot_->statData);
107
108		currSample = sampleTime + currentSnapshot_->getTime();
109	<	currStatus =  statusTime + currentSnapshot_->getTime();;
110	<	currThermal = thermalTime +  + currentSnapshot_->getTime();
109	>	currStatus =  statusTime + currentSnapshot_->getTime();
110	>	currThermal = thermalTime + currentSnapshot_->getTime();
111	>	if (needReset) {
112	>	currReset = resetTime + currentSnapshot_->getTime();
113	>	}
114	>	if (simParams->getUseRNEMD()){
115	>	currRNEMD = RNEMD_exchangeTime + currentSnapshot_->getTime();
116	>	}
117		needPotential = false;
118		needStress = false;
119
112	–	}
113	–
114	–	void VelocityVerletIntegrator::doIntegrate() {
115	–
116	–
117	–	initialize();
118	–
119	–	while (currentSnapshot_->getTime() < runTime) {
120	–
121	–	preStep();
122	–
123	–	integrateStep();
124	–
125	–	postStep();
126	–
120		}
121	+
122	+	void VelocityVerletIntegrator::doIntegrate() {
123
124	<	finalize();
124	>	initialize();
125
126	<	}
126	>	while (currentSnapshot_->getTime() < runTime) {
127	>	preStep();
128	>	integrateStep();
129	>	postStep();
130	>	}
131	>	finalize();
132	>	}
133
134
135	<	void VelocityVerletIntegrator::preStep() {
136	<	double difference = currentSnapshot_->getTime() + dt - currStatus;
135	>	void VelocityVerletIntegrator::preStep() {
136	>	RealType difference = currentSnapshot_->getTime() + dt - currStatus;
137
138	<	if (difference > 0 || fabs(difference) < oopse::epsilon) {
139	<	needPotential = true;
140	<	needStress = true;
138	>	if (difference > 0 || fabs(difference) < OpenMD::epsilon) {
139	>	needPotential = true;
140	>	needStress = true;
141	>	}
142		}
141	–
142	–	}
143
144	<	void VelocityVerletIntegrator::postStep() {
144	>	void VelocityVerletIntegrator::postStep() {
145	>
146	>	//save snapshot
147	>	info_->getSnapshotManager()->advance();
148
149	<	//save snapshot
150	<	info_->getSnapshotManager()->advance();
151	<
152	<	//increase time
153	<	currentSnapshot_->increaseTime(dt);
154	<
155	<	if (needVelocityScaling) {
156	<	if (currentSnapshot_->getTime() >= currThermal) {
154	<	velocitizer_->velocitize(targetScalingTemp);
155	<	currThermal += thermalTime;
149	>	//increase time
150	>	currentSnapshot_->increaseTime(dt);
151	>
152	>	if (needVelocityScaling) {
153	>	if (currentSnapshot_->getTime() >= currThermal) {
154	>	velocitizer_->velocitize(targetScalingTemp);
155	>	currThermal += thermalTime;
156	>	}
157		}
158	<	}
159	<
160	<	if (currentSnapshot_->getTime() >= currSample) {
161	<	dumpWriter->writeDumpAndEor();
158	>	if (useRNEMD) {
159	>	if (currentSnapshot_->getTime() >= currRNEMD) {
160	>	rnemd_->doRNEMD();
161	>	currRNEMD += RNEMD_exchangeTime;
162	>	}
163	>	rnemd_->collectData();
164	>	}
165
166	<	if (simParams->getUseSolidThermInt())
167	<	restWriter->writeZangle();
166	>	if (currentSnapshot_->getTime() >= currSample) {
167	>	dumpWriter->writeDumpAndEor();
168	>
169	>	currSample += sampleTime;
170	>	}
171
172	<	currSample += sampleTime;
173	<	}
174	<
175	<	if (currentSnapshot_->getTime() >= currStatus) {
176	<	//save statistics, before writeStat,  we must save statistics
177	<	thermo.saveStat();
178	<	saveConservedQuantity();
179	<	statWriter->writeStat(currentSnapshot_->statData);
172	>	if (currentSnapshot_->getTime() >= currStatus) {
173	>	//save statistics, before writeStat,  we must save statistics
174	>	thermo.saveStat();
175	>	saveConservedQuantity();
176	>
177	>	if (simParams->getUseRNEMD()) {
178	>	rnemd_->getStatus();
179	>	}
180	>
181	>	statWriter->writeStat(currentSnapshot_->statData);
182	>
183	>	progressBar->setStatus(currentSnapshot_->getTime(), runTime);
184	>	progressBar->update();
185	>
186	>	needPotential = false;
187	>	needStress = false;
188	>	currStatus += statusTime;
189	>	}
190
191	<	needPotential = false;
192	<	needStress = false;
193	<	currStatus += statusTime;
191	>	if (needReset && currentSnapshot_->getTime() >= currReset) {
192	>	resetIntegrator();
193	>	currReset += resetTime;
194	>	}
195		}
178	–
179	–
180	–	}
196
197
198	<	void VelocityVerletIntegrator::finalize() {
199	<	dumpWriter->writeEor();
198	>	void VelocityVerletIntegrator::finalize() {
199	>	dumpWriter->writeEor();
200
201	<	if (simParams->getUseSolidThermInt()) {
202	<	restWriter->writeZangle();
188	<	delete restWriter;
189	<	restWriter = NULL;
190	<	}
201	>	delete dumpWriter;
202	>	delete statWriter;
203
204	<	delete dumpWriter;
205	<	delete statWriter;
206	<
195	<	dumpWriter = NULL;
196	<	statWriter = NULL;
197	<
198	<	}
204	>	dumpWriter = NULL;
205	>	statWriter = NULL;
206	>	}
207
208	<	void VelocityVerletIntegrator::integrateStep() {
208	>	void VelocityVerletIntegrator::integrateStep() {
209
210	<	moveA();
211	<	calcForce(needPotential, needStress);
212	<	moveB();
213	<	}
210	>	moveA();
211	>	calcForce();
212	>	moveB();
213	>	}
214
215
216	<	void VelocityVerletIntegrator::calcForce(bool needPotential,
217	<	bool needStress) {
218	<	forceMan_->calcForces(needPotential, needStress);
219	<	}
216	>	void VelocityVerletIntegrator::calcForce() {
217	>	forceMan_->calcForces();
218	>	flucQ_->applyConstraints();
219	>	}
220
221	<	DumpWriter* VelocityVerletIntegrator::createDumpWriter() {
222	<	return new DumpWriter(info_);
223	<	}
221	>	DumpWriter* VelocityVerletIntegrator::createDumpWriter() {
222	>	return new DumpWriter(info_);
223	>	}
224
225	<	StatWriter* VelocityVerletIntegrator::createStatWriter() {
226	<	// if solidThermInt is true, add extra information to the statfile
227	<	if (simParams->getUseSolidThermInt()){
228	<	StatsBitSet mask;
229	<	mask.set(Stats::TIME);
230	<	mask.set(Stats::TOTAL_ENERGY);
231	<	mask.set(Stats::POTENTIAL_ENERGY);
232	<	mask.set(Stats::KINETIC_ENERGY);
233	<	mask.set(Stats::TEMPERATURE);
234	<	mask.set(Stats::PRESSURE);
235	<	mask.set(Stats::VOLUME);
236	<	mask.set(Stats::CONSERVED_QUANTITY);
237	<	mask.set(Stats::VRAW);
238	<	mask.set(Stats::VHARM);
225	>	StatWriter* VelocityVerletIntegrator::createStatWriter() {
226	>
227	>	std::string statFileFormatString = simParams->getStatFileFormat();
228	>	StatsBitSet mask = parseStatFileFormat(statFileFormatString);
229	>
230	>	// if we're doing a thermodynamic integration, we'll want the raw
231	>	// potential as well as the full potential:
232	>
233	>
234	>	if (simParams->getUseThermodynamicIntegration())
235	>	mask.set(Stats::VRAW);
236	>
237	>	// if we've got restraints turned on, we'll also want a report of the
238	>	// total harmonic restraints
239	>	if (simParams->getUseRestraints()){
240	>	mask.set(Stats::VHARM);
241	>	}
242	>
243	>	if (simParams->havePrintPressureTensor() &&
244	>	simParams->getPrintPressureTensor()){
245	>	mask.set(Stats::PRESSURE_TENSOR_XX);
246	>	mask.set(Stats::PRESSURE_TENSOR_XY);
247	>	mask.set(Stats::PRESSURE_TENSOR_XZ);
248	>	mask.set(Stats::PRESSURE_TENSOR_YX);
249	>	mask.set(Stats::PRESSURE_TENSOR_YY);
250	>	mask.set(Stats::PRESSURE_TENSOR_YZ);
251	>	mask.set(Stats::PRESSURE_TENSOR_ZX);
252	>	mask.set(Stats::PRESSURE_TENSOR_ZY);
253	>	mask.set(Stats::PRESSURE_TENSOR_ZZ);
254	>	}
255	>
256	>	if (simParams->getAccumulateBoxDipole()) {
257	>	mask.set(Stats::BOX_DIPOLE_X);
258	>	mask.set(Stats::BOX_DIPOLE_Y);
259	>	mask.set(Stats::BOX_DIPOLE_Z);
260	>	}
261	>
262	>	if (simParams->getPrintHeatFlux()) {
263	>	mask.set(Stats::HEATFLUX_X);
264	>	mask.set(Stats::HEATFLUX_Y);
265	>	mask.set(Stats::HEATFLUX_Z);
266	>	}
267	>
268	>	if (simParams->haveTaggedAtomPair() && simParams->havePrintTaggedPairDistance()) {
269	>	if (simParams->getPrintTaggedPairDistance()) {
270	>	mask.set(Stats::TAGGED_PAIR_DISTANCE);
271	>	}
272	>	}
273	>
274	>	if (simParams->getUseRNEMD()) {
275	>	mask.set(Stats::RNEMD_EXCHANGE_TOTAL);
276	>	}
277	>
278	>
279		return new StatWriter(info_->getStatFileName(), mask);
280		}
233	–
234	–	return new StatWriter(info_->getStatFileName());
235	–	}
281
237	–	RestWriter* VelocityVerletIntegrator::createRestWriter(){
238	–	return new RestWriter(info_);
239	–	}
282
283	<
242	<	} //end namespace oopse
283	>	} //end namespace OpenMD

Comparing:
trunk/src/integrators/VelocityVerletIntegrator.cpp (property svn:keywords), Revision 471 by chrisfen, Tue Apr 12 14:24:29 2005 UTC vs.
branches/development/src/integrators/VelocityVerletIntegrator.cpp (property svn:keywords), Revision 1723 by gezelter, Thu May 24 20:59:54 2012 UTC

#	Line 0 | Line 1
1	+	Author Id Revision Date

Diff Legend

–	Removed lines
+	Added lines
<	Changed lines
>	Changed lines