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

Comparing trunk/src/integrators/VelocityVerletIntegrator.cpp (file contents):
Revision 276 by tim, Wed Jan 26 15:26:47 2005 UTC vs.
Revision 1341 by skuang, Fri May 8 19:47:05 2009 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 52 | Line 52 | namespace oopse {
52   #include "utils/StringUtils.hpp"
53  
54   namespace oopse {
55 < VelocityVerletIntegrator::VelocityVerletIntegrator(SimInfo *info) : Integrator(info), rotAlgo(NULL) {
55 >  VelocityVerletIntegrator::VelocityVerletIntegrator(SimInfo *info) : Integrator(info), rotAlgo(NULL) {
56      dt2 = 0.5 * dt;
57      rotAlgo = new DLM();
58      rattle = new Rattle(info);
59 < }
60 <
61 < VelocityVerletIntegrator::~VelocityVerletIntegrator() {
59 >  }
60 >  
61 >  VelocityVerletIntegrator::~VelocityVerletIntegrator() {
62      delete rotAlgo;
63      delete rattle;
64 < }
65 <
66 < void VelocityVerletIntegrator::initialize(){
67 <
64 >  }
65 >  
66 >  void VelocityVerletIntegrator::initialize(){
67 >    
68      forceMan_->init();
69 <
70 <    // remove center of mass drift velocity (in case we passed in a configuration
71 <    // that was drifting
69 >    
70 >    // remove center of mass drift velocity (in case we passed in a
71 >    // configuration that was drifting)
72      velocitizer_->removeComDrift();
73 <
73 >    
74      // initialize the forces before the first step
75      calcForce(true, true);
76 <
77 <    //execute constraint algorithm to make sure at the very beginning the system is constrained  
76 >    
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();
81 <        calcForce(true, true);
82 <        rattle->constraintB();        
83 <        info_->getSnapshotManager()->advance();//copy the current snapshot to previous snapshot
80 >      rattle->constraintA();
81 >      calcForce(true, true);
82 >      rattle->constraintB();      
83 >      //copy the current snapshot to previous snapshot
84 >      info_->getSnapshotManager()->advance();
85      }
86 <
86 >    
87      if (needVelocityScaling) {
88 <        velocitizer_->velocitize(targetScalingTemp);
88 >      velocitizer_->velocitize(targetScalingTemp);
89      }
90      
91      dumpWriter = createDumpWriter();
92 <    eorWriter = createEorWriter();
92 >    
93      statWriter = createStatWriter();
94 +
95 +    dumpWriter->writeDumpAndEor();
96  
97 <    dumpWriter->writeDump();
98 <    eorWriter->writeDump();
97 >    if (simParams->getUseSolidThermInt()) {
98 >      restWriter = createRestWriter();
99 >      restWriter->writeZAngFile();
100 >    }
101      
102      //save statistics, before writeStat,  we must save statistics
103      thermo.saveStat();
104      saveConservedQuantity();
105 <    statWriter->writeStat(currentSnapshot_->statData);
105 >    if (simParams->getUseRNEMD())
106 >      rnemd_->getStatus();
107  
108 +    statWriter->writeStat(currentSnapshot_->statData);
109 +    
110      currSample = sampleTime + currentSnapshot_->getTime();
111 <    currStatus =  statusTime + currentSnapshot_->getTime();;
112 <    currThermal = thermalTime +  + currentSnapshot_->getTime();
111 >    currStatus =  statusTime + currentSnapshot_->getTime();
112 >    currThermal = thermalTime + currentSnapshot_->getTime();
113 >    if (needReset) {
114 >      currReset = resetTime + currentSnapshot_->getTime();
115 >    }
116 >    if (simParams->getUseRNEMD()){
117 >      currRNEMD = RNEMD_swapTime + currentSnapshot_->getTime();
118 >    }
119      needPotential = false;
120      needStress = false;      
121 <  
122 < }
108 <        
109 < void VelocityVerletIntegrator::doIntegrate() {
121 >    
122 >  }
123  
124 <
124 >  void VelocityVerletIntegrator::doIntegrate() {
125 >  
126 >  
127      initialize();
128 <
128 >  
129      while (currentSnapshot_->getTime() < runTime) {
130 <        
131 <        preStep();
132 <
133 <        integrateStep();
134 <        
135 <        postStep();
136 <
130 >    
131 >      preStep();
132 >    
133 >      integrateStep();
134 >    
135 >      postStep();
136 >    
137      }
138 <
138 >  
139      finalize();
140 <    
141 < }
140 >  
141 >  }
142  
143  
144 < void VelocityVerletIntegrator::preStep() {
145 <        double difference = currentSnapshot_->getTime() + dt - currStatus;
144 >  void VelocityVerletIntegrator::preStep() {
145 >    RealType difference = currentSnapshot_->getTime() + dt - currStatus;
146 >  
147 >    if (difference > 0 || fabs(difference) < oopse::epsilon) {
148 >      needPotential = true;
149 >      needStress = true;  
150 >    }
151 >  }
152  
153 <        if (difference > 0 || fabs(difference) < oopse::epsilon) {
133 <            needPotential = true;
134 <            needStress = true;  
135 <        }
153 >  void VelocityVerletIntegrator::postStep() {
154  
155 < }
155 >    //save snapshot
156 >    info_->getSnapshotManager()->advance();
157 >  
158 >    //increase time
159 >    currentSnapshot_->increaseTime(dt);        
160 >  
161 >    if (needVelocityScaling) {
162 >      if (currentSnapshot_->getTime() >= currThermal) {
163 >        velocitizer_->velocitize(targetScalingTemp);
164 >        currThermal += thermalTime;
165 >      }
166 >    }
167 >    if (useRNEMD) {
168 >      if (currentSnapshot_->getTime() >= currRNEMD) {
169 >        rnemd_->doSwap();
170 >        currRNEMD += RNEMD_swapTime;
171 >      }
172 >    }
173 >    
174 >    if (currentSnapshot_->getTime() >= currSample) {
175 >      dumpWriter->writeDumpAndEor();
176 >      
177 >      if (simParams->getUseSolidThermInt())
178 >        restWriter->writeZAngFile();
179 >      
180 >      currSample += sampleTime;
181 >    }
182 >    
183 >    if (currentSnapshot_->getTime() >= currStatus) {
184 >      //save statistics, before writeStat,  we must save statistics
185 >      thermo.saveStat();
186 >      saveConservedQuantity();
187  
188 < void VelocityVerletIntegrator::postStep() {
188 >      if (simParams->getUseRNEMD())
189 >        rnemd_->getStatus();
190  
191 <        //save snapshot
192 <        info_->getSnapshotManager()->advance();
191 >      statWriter->writeStat(currentSnapshot_->statData);
192 >      
193 >      needPotential = false;
194 >      needStress = false;
195 >      currStatus += statusTime;
196 >    }
197 >    
198 >    if (needReset && currentSnapshot_->getTime() >= currReset) {    
199 >      resetIntegrator();
200 >      currReset += resetTime;
201 >    }        
202 >  }
203  
144        //increase time
145        currentSnapshot_->increaseTime(dt);        
146        
147        if (needVelocityScaling) {
148            if (currentSnapshot_->getTime() >= currThermal) {
149                velocitizer_->velocitize(targetScalingTemp);
150                currThermal += thermalTime;
151            }
152        }
204  
205 <        if (currentSnapshot_->getTime() >= currSample) {
206 <            dumpWriter->writeDump();
207 <            eorWriter->writeDump();
208 <            currSample += sampleTime;
209 <        }
210 <
211 <        if (currentSnapshot_->getTime() >= currStatus) {
212 <            //save statistics, before writeStat,  we must save statistics
213 <            thermo.saveStat();
163 <            saveConservedQuantity();
164 <            statWriter->writeStat(currentSnapshot_->statData);
165 <
166 <            needPotential = false;
167 <            needStress = false;
168 <            currStatus += statusTime;
169 <        }
170 <
171 <        
172 < }
173 <
174 <
175 < void VelocityVerletIntegrator::finalize() {
176 <
177 <    dumpWriter->writeDump();
178 <    eorWriter->writeDump();
179 <
205 >  void VelocityVerletIntegrator::finalize() {
206 >    dumpWriter->writeEor();
207 >  
208 >    if (simParams->getUseSolidThermInt()) {
209 >      restWriter->writeZAngFile();
210 >      delete restWriter;
211 >      restWriter = NULL;
212 >    }
213 >  
214      delete dumpWriter;
181    delete eorWriter;
215      delete statWriter;
216 <
216 >  
217      dumpWriter = NULL;
185    eorWriter = NULL;
218      statWriter = NULL;
219 <    
220 < }
219 >  
220 >  }
221  
222 < void VelocityVerletIntegrator::integrateStep() {
223 <
222 >  void VelocityVerletIntegrator::integrateStep() {
223 >  
224      moveA();
225      calcForce(needPotential, needStress);
226      moveB();
227 < }
227 >  }
228  
229  
230 < void VelocityVerletIntegrator::calcForce(bool needPotential,
231 <                                         bool needStress) {
230 >  void VelocityVerletIntegrator::calcForce(bool needPotential,
231 >                                           bool needStress) {
232      forceMan_->calcForces(needPotential, needStress);
233 < }
233 >  }
234  
235 < DumpWriter* VelocityVerletIntegrator::createDumpWriter() {
236 <    return new DumpWriter(info_, info_->getDumpFileName());
237 < }
235 >  DumpWriter* VelocityVerletIntegrator::createDumpWriter() {
236 >    return new DumpWriter(info_);
237 >  }
238  
239 < DumpWriter* VelocityVerletIntegrator::createEorWriter() {
208 <    return new DumpWriter(info_, getPrefix(info_->getDumpFileName()) + ".eor");
209 < }
239 >  StatWriter* VelocityVerletIntegrator::createStatWriter() {
240  
241 < StatWriter* VelocityVerletIntegrator::createStatWriter() {
242 <    return new StatWriter(info_->getStatFileName());
243 < }
241 >    std::string statFileFormatString = simParams->getStatFileFormat();
242 >    StatsBitSet mask = parseStatFileFormat(statFileFormatString);
243 >    
244 >    // if solidThermInt is true, add extra information to the statfile
245 >    if (simParams->getUseSolidThermInt()){
246 >      mask.set(Stats::VRAW);
247 >      mask.set(Stats::VHARM);
248 >    }
249  
250 +    if (simParams->havePrintPressureTensor() &&
251 +        simParams->getPrintPressureTensor()){
252 +      mask.set(Stats::PRESSURE_TENSOR_XX);
253 +      mask.set(Stats::PRESSURE_TENSOR_XY);
254 +      mask.set(Stats::PRESSURE_TENSOR_XZ);
255 +      mask.set(Stats::PRESSURE_TENSOR_YX);
256 +      mask.set(Stats::PRESSURE_TENSOR_YY);
257 +      mask.set(Stats::PRESSURE_TENSOR_YZ);
258 +      mask.set(Stats::PRESSURE_TENSOR_ZX);
259 +      mask.set(Stats::PRESSURE_TENSOR_ZY);
260 +      mask.set(Stats::PRESSURE_TENSOR_ZZ);
261 +    }
262 +    
263 +    if (simParams->getAccumulateBoxDipole()) {
264 +      mask.set(Stats::BOX_DIPOLE_X);
265 +      mask.set(Stats::BOX_DIPOLE_Y);
266 +      mask.set(Stats::BOX_DIPOLE_Z);
267 +    }
268 +  
269 +    if (simParams->haveTaggedAtomPair() && simParams->havePrintTaggedPairDistance()) {
270 +      if (simParams->getPrintTaggedPairDistance()) {
271 +        mask.set(Stats::TAGGED_PAIR_DISTANCE);
272 +      }
273 +    }
274  
275 +    if (simParams->getUseRNEMD()) {
276 +      mask.set(Stats::RNEMD_SWAP_TOTAL);
277 +    }
278 +      
279 +
280 +     return new StatWriter(info_->getStatFileName(), mask);
281 +  }
282 +
283 +  RestWriter* VelocityVerletIntegrator::createRestWriter(){
284 +    return new RestWriter(info_);
285 +  }
286 +
287 +
288   } //end namespace oopse

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines