ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/branches/development/src/integrators/VelocityVerletIntegrator.cpp
(Generate patch)

Comparing branches/development/src/integrators/VelocityVerletIntegrator.cpp (file contents):
Revision 1665 by gezelter, Tue Nov 22 20:38:56 2011 UTC vs.
Revision 1808 by gezelter, Mon Oct 22 20:42:10 2012 UTC

# Line 44 | Line 44
44   * @file VelocityVerletIntegrator.cpp
45   * @author tlin
46   * @date 11/09/2004
47 – * @time 16:16am
47   * @version 1.0
48   */
49  
# Line 54 | Line 53 | namespace OpenMD {
53   #include "utils/ProgressBar.hpp"
54  
55   namespace OpenMD {
56 <  VelocityVerletIntegrator::VelocityVerletIntegrator(SimInfo *info) : Integrator(info), rotAlgo(NULL) {
56 >  VelocityVerletIntegrator::VelocityVerletIntegrator(SimInfo *info) : Integrator(info) {
57      dt2 = 0.5 * dt;
59 –    rotAlgo = new DLM();
60 –    rattle = new Rattle(info);
58    }
59    
60    VelocityVerletIntegrator::~VelocityVerletIntegrator() {
64 –    delete rotAlgo;
65 –    delete rattle;
61    }
62    
63    void VelocityVerletIntegrator::initialize(){
# Line 72 | Line 67 | namespace OpenMD {
67      // remove center of mass drift velocity (in case we passed in a
68      // configuration that was drifting)
69      velocitizer_->removeComDrift();
70 +
71 +    // find the initial fluctuating charges.
72 +    flucQ_->initialize();
73      
74      // initialize the forces before the first step
75      calcForce();
# Line 79 | Line 77 | namespace OpenMD {
77      // execute the constraint algorithm to make sure that the system is
78      // constrained at the very beginning  
79      if (info_->getNGlobalConstraints() > 0) {
80 <      rattle->constraintA();
80 >      rattle_->constraintA();
81        calcForce();
82 <      rattle->constraintB();      
82 >      rattle_->constraintB();      
83        //copy the current snapshot to previous snapshot
84        info_->getSnapshotManager()->advance();
85      }
# Line 99 | Line 97 | namespace OpenMD {
97      progressBar = new ProgressBar();
98  
99      //save statistics, before writeStat,  we must save statistics
102 –    thermo.saveStat();
100      saveConservedQuantity();
101 <    if (simParams->getUseRNEMD())
101 >    stats->collectStats();
102 >
103 >    if (simParams->getRNEMDParameters()->getUseRNEMD())
104        rnemd_->getStarted();
105  
106 <    statWriter->writeStat(currentSnapshot_->statData);
106 >    statWriter->writeStat();
107      
108 <    currSample = sampleTime + currentSnapshot_->getTime();
109 <    currStatus =  statusTime + currentSnapshot_->getTime();
110 <    currThermal = thermalTime + currentSnapshot_->getTime();
108 >    currSample = sampleTime + snap->getTime();
109 >    currStatus =  statusTime + snap->getTime();
110 >    currThermal = thermalTime + snap->getTime();
111      if (needReset) {
112 <      currReset = resetTime + currentSnapshot_->getTime();
112 >      currReset = resetTime + snap->getTime();
113      }
114 <    if (simParams->getUseRNEMD()){
115 <      currRNEMD = RNEMD_exchangeTime + currentSnapshot_->getTime();
114 >    if (simParams->getRNEMDParameters()->getUseRNEMD()){
115 >      currRNEMD = RNEMD_exchangeTime + snap->getTime();
116      }
117      needPotential = false;
118      needStress = false;      
# Line 122 | Line 121 | namespace OpenMD {
121  
122    void VelocityVerletIntegrator::doIntegrate() {
123    
125 –  
124      initialize();
125    
126 <    while (currentSnapshot_->getTime() < runTime) {
127 <    
128 <      preStep();
129 <    
130 <      integrateStep();
133 <    
134 <      postStep();
135 <      
136 <    }
137 <  
126 >    while (snap->getTime() < runTime) {    
127 >      preStep();    
128 >      integrateStep();    
129 >      postStep();      
130 >    }  
131      finalize();
139 –  
132    }
133  
134  
135    void VelocityVerletIntegrator::preStep() {
136 <    RealType difference = currentSnapshot_->getTime() + dt - currStatus;
136 >    RealType difference = snap->getTime() + dt - currStatus;
137    
138      if (difference > 0 || fabs(difference) < OpenMD::epsilon) {
139        needPotential = true;
# Line 155 | Line 147 | namespace OpenMD {
147      info_->getSnapshotManager()->advance();
148    
149      //increase time
150 <    currentSnapshot_->increaseTime(dt);        
150 >    snap->increaseTime(dt);        
151    
152      if (needVelocityScaling) {
153 <      if (currentSnapshot_->getTime() >= currThermal) {
153 >      if (snap->getTime() >= currThermal) {
154          velocitizer_->velocitize(targetScalingTemp);
155          currThermal += thermalTime;
156        }
157      }
158      if (useRNEMD) {
159 <      if (currentSnapshot_->getTime() >= currRNEMD) {
159 >      if (snap->getTime() >= currRNEMD) {
160          rnemd_->doRNEMD();
161          currRNEMD += RNEMD_exchangeTime;
162        }
163        rnemd_->collectData();
164      }
165      
166 <    if (currentSnapshot_->getTime() >= currSample) {
166 >    if (snap->getTime() >= currSample) {
167        dumpWriter->writeDumpAndEor();
168        
169        currSample += sampleTime;
170      }
171      
172 <    if (currentSnapshot_->getTime() >= currStatus) {
172 >    if (snap->getTime() >= currStatus) {
173        //save statistics, before writeStat,  we must save statistics
174 <      thermo.saveStat();
174 >      stats->collectStats();
175        saveConservedQuantity();
176  
177 <      if (simParams->getUseRNEMD()) {
178 <        rnemd_->getStatus();
177 >      if (simParams->getRNEMDParameters()->getUseRNEMD()) {
178 >        rnemd_->writeOutputFile();
179        }
180  
181 <      statWriter->writeStat(currentSnapshot_->statData);
181 >      statWriter->writeStat();
182  
183 <      progressBar->setStatus(currentSnapshot_->getTime(), runTime);
183 >      progressBar->setStatus(snap->getTime(), runTime);
184        progressBar->update();
185  
186        needPotential = false;
# Line 196 | Line 188 | namespace OpenMD {
188        currStatus += statusTime;
189      }
190      
191 <    if (needReset && currentSnapshot_->getTime() >= currReset) {    
191 >    if (needReset && snap->getTime() >= currReset) {    
192        resetIntegrator();
193        currReset += resetTime;
194      }        
# Line 205 | Line 197 | namespace OpenMD {
197  
198    void VelocityVerletIntegrator::finalize() {
199      dumpWriter->writeEor();
200 <  
200 >    if (simParams->getRNEMDParameters()->getUseRNEMD()) {
201 >      rnemd_->writeOutputFile();
202 >    }
203 >
204      delete dumpWriter;
205      delete statWriter;
206    
207      dumpWriter = NULL;
208      statWriter = NULL;
214 –  
209    }
210  
211    void VelocityVerletIntegrator::integrateStep() {
# Line 224 | Line 218 | namespace OpenMD {
218  
219    void VelocityVerletIntegrator::calcForce() {
220      forceMan_->calcForces();
221 +    flucQ_->applyConstraints();
222    }
223  
224    DumpWriter* VelocityVerletIntegrator::createDumpWriter() {
# Line 231 | Line 226 | namespace OpenMD {
226    }
227  
228    StatWriter* VelocityVerletIntegrator::createStatWriter() {
234 –
235 –    std::string statFileFormatString = simParams->getStatFileFormat();
236 –    StatsBitSet mask = parseStatFileFormat(statFileFormatString);
237 –  
238 –    // if we're doing a thermodynamic integration, we'll want the raw
239 –    // potential as well as the full potential:
240 –
241 –
242 –    if (simParams->getUseThermodynamicIntegration())
243 –      mask.set(Stats::VRAW);
244 –
245 –    // if we've got restraints turned on, we'll also want a report of the
246 –    // total harmonic restraints
247 –    if (simParams->getUseRestraints()){
248 –      mask.set(Stats::VHARM);
249 –    }
250 –
251 –    if (simParams->havePrintPressureTensor() &&
252 –        simParams->getPrintPressureTensor()){
253 –      mask.set(Stats::PRESSURE_TENSOR_XX);
254 –      mask.set(Stats::PRESSURE_TENSOR_XY);
255 –      mask.set(Stats::PRESSURE_TENSOR_XZ);
256 –      mask.set(Stats::PRESSURE_TENSOR_YX);
257 –      mask.set(Stats::PRESSURE_TENSOR_YY);
258 –      mask.set(Stats::PRESSURE_TENSOR_YZ);
259 –      mask.set(Stats::PRESSURE_TENSOR_ZX);
260 –      mask.set(Stats::PRESSURE_TENSOR_ZY);
261 –      mask.set(Stats::PRESSURE_TENSOR_ZZ);
262 –    }
229      
230 <    if (simParams->getAccumulateBoxDipole()) {
231 <      mask.set(Stats::BOX_DIPOLE_X);
266 <      mask.set(Stats::BOX_DIPOLE_Y);
267 <      mask.set(Stats::BOX_DIPOLE_Z);
268 <    }
269 <  
270 <    if (simParams->haveTaggedAtomPair() && simParams->havePrintTaggedPairDistance()) {
271 <      if (simParams->getPrintTaggedPairDistance()) {
272 <        mask.set(Stats::TAGGED_PAIR_DISTANCE);
273 <      }
274 <    }
275 <
276 <    if (simParams->getUseRNEMD()) {
277 <      mask.set(Stats::RNEMD_EXCHANGE_TOTAL);
278 <    }
230 >    stats = new Stats(info_);
231 >    statWriter = new StatWriter(info_->getStatFileName(), stats);
232      
233 <
281 <    return new StatWriter(info_->getStatFileName(), mask);
233 >    return statWriter;
234    }
235  
284 –
236   } //end namespace OpenMD

Diff Legend

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