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 376 by tim, Thu Feb 24 20:55:07 2005 UTC vs.
Revision 1329 by gezelter, Sat Mar 14 01:23:53 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 <
92 >    
93      statWriter = createStatWriter();
94 <
94 >
95      dumpWriter->writeDumpAndEor();
96 +
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);
106 <
106 >    
107      currSample = sampleTime + currentSnapshot_->getTime();
108 <    currStatus =  statusTime + currentSnapshot_->getTime();;
109 <    currThermal = thermalTime +  + currentSnapshot_->getTime();
108 >    currStatus =  statusTime + currentSnapshot_->getTime();
109 >    currThermal = thermalTime + currentSnapshot_->getTime();
110 >    if (needReset) {
111 >      currReset = resetTime + currentSnapshot_->getTime();
112 >    }
113 >    if (simParams->getUseRNEMD()){
114 >      currRNEMD = RNEMD_swapTime + currentSnapshot_->getTime();
115 >    }
116      needPotential = false;
117      needStress = false;      
118 <  
119 < }
107 <        
108 < void VelocityVerletIntegrator::doIntegrate() {
118 >    
119 >  }
120  
121 <
121 >  void VelocityVerletIntegrator::doIntegrate() {
122 >  
123 >  
124      initialize();
125 <
125 >  
126      while (currentSnapshot_->getTime() < runTime) {
127 <        
128 <        preStep();
129 <
130 <        integrateStep();
131 <        
132 <        postStep();
133 <
127 >    
128 >      preStep();
129 >    
130 >      integrateStep();
131 >    
132 >      postStep();
133 >    
134      }
135 <
135 >  
136      finalize();
137 <    
138 < }
137 >  
138 >  }
139  
140  
141 < void VelocityVerletIntegrator::preStep() {
142 <        double difference = currentSnapshot_->getTime() + dt - currStatus;
141 >  void VelocityVerletIntegrator::preStep() {
142 >    RealType difference = currentSnapshot_->getTime() + dt - currStatus;
143 >  
144 >    if (difference > 0 || fabs(difference) < oopse::epsilon) {
145 >      needPotential = true;
146 >      needStress = true;  
147 >    }
148 >  }
149  
150 <        if (difference > 0 || fabs(difference) < oopse::epsilon) {
132 <            needPotential = true;
133 <            needStress = true;  
134 <        }
150 >  void VelocityVerletIntegrator::postStep() {
151  
152 < }
152 >    //save snapshot
153 >    info_->getSnapshotManager()->advance();
154 >  
155 >    //increase time
156 >    currentSnapshot_->increaseTime(dt);        
157 >  
158 >    if (needVelocityScaling) {
159 >      if (currentSnapshot_->getTime() >= currThermal) {
160 >        velocitizer_->velocitize(targetScalingTemp);
161 >        currThermal += thermalTime;
162 >      }
163 >    }
164 >    if (useRNEMD) {
165 >      if (currentSnapshot_->getTime() >= currRNEMD) {
166 >        rnemd_->doSwap();
167 >        currRNEMD += RNEMD_swapTime;
168 >      }
169 >    }
170 >    
171 >    if (currentSnapshot_->getTime() >= currSample) {
172 >      dumpWriter->writeDumpAndEor();
173 >      
174 >      if (simParams->getUseSolidThermInt())
175 >        restWriter->writeZAngFile();
176 >      
177 >      currSample += sampleTime;
178 >    }
179 >    
180 >    if (currentSnapshot_->getTime() >= currStatus) {
181 >      //save statistics, before writeStat,  we must save statistics
182 >      thermo.saveStat();
183 >      saveConservedQuantity();
184 >      statWriter->writeStat(currentSnapshot_->statData);
185 >      
186 >      needPotential = false;
187 >      needStress = false;
188 >      currStatus += statusTime;
189 >    }
190 >    
191 >    if (needReset && currentSnapshot_->getTime() >= currReset) {    
192 >      resetIntegrator();
193 >      currReset += resetTime;
194 >    }        
195 >  }
196  
138 – void VelocityVerletIntegrator::postStep() {
197  
198 <        //save snapshot
199 <        info_->getSnapshotManager()->advance();
200 <
201 <        //increase time
202 <        currentSnapshot_->increaseTime(dt);        
203 <        
204 <        if (needVelocityScaling) {
205 <            if (currentSnapshot_->getTime() >= currThermal) {
206 <                velocitizer_->velocitize(targetScalingTemp);
149 <                currThermal += thermalTime;
150 <            }
151 <        }
152 <
153 <        if (currentSnapshot_->getTime() >= currSample) {
154 <            dumpWriter->writeDumpAndEor();
155 <            currSample += sampleTime;
156 <        }
157 <
158 <        if (currentSnapshot_->getTime() >= currStatus) {
159 <            //save statistics, before writeStat,  we must save statistics
160 <            thermo.saveStat();
161 <            saveConservedQuantity();
162 <            statWriter->writeStat(currentSnapshot_->statData);
163 <
164 <            needPotential = false;
165 <            needStress = false;
166 <            currStatus += statusTime;
167 <        }
168 <
169 <        
170 < }
171 <
172 <
173 < void VelocityVerletIntegrator::finalize() {
174 <     dumpWriter->writeEor();
175 <
198 >  void VelocityVerletIntegrator::finalize() {
199 >    dumpWriter->writeEor();
200 >  
201 >    if (simParams->getUseSolidThermInt()) {
202 >      restWriter->writeZAngFile();
203 >      delete restWriter;
204 >      restWriter = NULL;
205 >    }
206 >  
207      delete dumpWriter;
208      delete statWriter;
209 <
209 >  
210      dumpWriter = NULL;
211      statWriter = NULL;
212 <    
213 < }
212 >  
213 >  }
214  
215 < void VelocityVerletIntegrator::integrateStep() {
216 <
215 >  void VelocityVerletIntegrator::integrateStep() {
216 >  
217      moveA();
218      calcForce(needPotential, needStress);
219      moveB();
220 < }
220 >  }
221  
222  
223 < void VelocityVerletIntegrator::calcForce(bool needPotential,
224 <                                         bool needStress) {
223 >  void VelocityVerletIntegrator::calcForce(bool needPotential,
224 >                                           bool needStress) {
225      forceMan_->calcForces(needPotential, needStress);
226 < }
226 >  }
227  
228 < DumpWriter* VelocityVerletIntegrator::createDumpWriter() {
228 >  DumpWriter* VelocityVerletIntegrator::createDumpWriter() {
229      return new DumpWriter(info_);
230 < }
230 >  }
231  
232 < StatWriter* VelocityVerletIntegrator::createStatWriter() {
202 <    return new StatWriter(info_->getStatFileName());
203 < }
232 >  StatWriter* VelocityVerletIntegrator::createStatWriter() {
233  
234 +    std::string statFileFormatString = simParams->getStatFileFormat();
235 +    StatsBitSet mask = parseStatFileFormat(statFileFormatString);
236 +    
237 +    // if solidThermInt is true, add extra information to the statfile
238 +    if (simParams->getUseSolidThermInt()){
239 +      mask.set(Stats::VRAW);
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->haveTaggedAtomPair() && simParams->havePrintTaggedPairDistance()) {
263 +      if (simParams->getPrintTaggedPairDistance()) {
264 +        mask.set(Stats::TAGGED_PAIR_DISTANCE);
265 +      }
266 +    }
267 +
268 +    if (simParams->getUseRNEMD()) {
269 +      mask.set(Stats::RNEMD_SWAP_TOTAL);
270 +    }
271 +      
272 +
273 +     return new StatWriter(info_->getStatFileName(), mask);
274 +  }
275 +
276 +  RestWriter* VelocityVerletIntegrator::createRestWriter(){
277 +    return new RestWriter(info_);
278 +  }
279 +
280 +
281   } //end namespace oopse

Diff Legend

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