ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/branches/new_design/OOPSE-3.0/src/integrators/VelocityVerletIntegrator.cpp
(Generate patch)

Comparing branches/new_design/OOPSE-3.0/src/integrators/VelocityVerletIntegrator.cpp (file contents):
Revision 1756 by tim, Thu Nov 18 23:26:27 2004 UTC vs.
Revision 1820 by tim, Thu Dec 2 00:09:35 2004 UTC

# Line 32 | Line 32
32   */
33  
34   #include "integrators/VelocityVerletIntegrator.hpp"
35 + #include "integrators/DLM.hpp"
36  
37   namespace oopse {
38 < VelocityVerletIntegrator::VelocityVerletIntegrator(SimInfo *info) : Integrator(info) { }
38 > VelocityVerletIntegrator::VelocityVerletIntegrator(SimInfo *info) : Integrator(info), rotAlgo(NULL) {
39 >    rotAlgo = new DLM();
40 >    dt2 = 0.5 * dt;
41 > }
42  
43 < VelocityVerletIntegrator::~VelocityVerletIntegrator() { }
43 > VelocityVerletIntegrator::~VelocityVerletIntegrator() {
44 >    delete rotAlgo;
45 > }
46  
47 < void VelocityVerletIntegrator::integrate() {
42 <    double runTime = info_->run_time;
47 > void VelocityVerletIntegrator::initialize(){
48  
49 <    int calcPot;
50 <    int calcStress;
49 >    currSample = 0.0;
50 >    currStatus = 0.0;
51 >    currThermal = 0.0;
52 >    needPotential = false;
53 >    needStress = false;      
54  
47
48    fullStep_ = info_->dt;
49    halfStep_ = 0.5 * fullStep_;
50
51    readyCheck();
52
55      // remove center of mass drift velocity (in case we passed in a configuration
56      // that was drifting
57 <    tStats->removeCOMdrift();
57 >    velocitizer_->removeComDrift();
58  
59      // initialize the forces before the first step
60 +    calcForce(true, true);
61  
59    calcForce(1, 1);
60
62      //execute constraint algorithm to make sure at the very beginning the system is constrained  
63 <    if (nConstrained) {
64 <        constrainA();
65 <        calcForce(1, 1);
66 <        constrainB();
67 <    }
63 >    //if (nConstrained) {
64 >    //    constrainA();
65 >    //    calcForce(true, true);
66 >    //    constrainB();
67 >    //}
68  
69 <    if (info_->setTemp) {
70 <        thermalize();
69 >    if (needVelocityScaling) {
70 >        velocitizer_->velocitize(targetScalingTemp);
71      }
72 +    
73 +    dumpWriter = createDumpWriter();
74 +    statWriter = createStatWriter();
75  
76 <    calcPot = 0;
77 <    calcStress = 0;
76 >    dumpWriter->writeDump();
77 >    statWriter->writeStat(currentSnapshot_->statData);
78 >  
79 > }
80 >        
81 > void VelocityVerletIntegrator::doIntegrate() {
82  
75    tStats = new Thermo(info_);
76    statOut = new StatWriter(info_);
77    dumpOut = new DumpWriter(info_);
78    
79    dumpOut->writeDump(info_->getTime());
80    statOut->writeStat(info_->getTime());
83  
84 < #ifdef IS_MPI
84 >    initialize();
85  
86 <    strcpy(checkPointMsg, "The integrator is ready to go.");
87 <    MPIcheckPoint();
86 >    while (currentSnapshot_->getTime() < runTime) {
87 >        
88 >        preStep();
89  
90 < #endif // is_mpi
90 >        integrateStep();
91 >        
92 >        postStep();
93  
94 <    while (info_->getTime() < runTime) {
90 <        difference = info_->getTime() + fullStep_ - currStatus;
94 >    }
95  
96 <        if (difference > 0 || fabs(difference) < 1e - 4) {
97 <            calcPot = 1;
98 <            calcStress = 1;
96 >    finalize();
97 >    
98 > }
99 >
100 >
101 > void VelocityVerletIntegrator::preStep() {
102 >        double difference = currentSnapshot_->getTime() + dt - currStatus;
103 >
104 >        if (difference > 0 || fabs(difference) < oopse::epsilon) {
105 >            needPotential = true;
106 >            needStress = true;  
107          }
108  
109 <        //notify before integrateStep
110 <        notify()
111 <        integrateStep(calcPot, calcStress);
109 > }
110 >
111 > void VelocityVerletIntegrator::postStep() {
112 >
113 >        //save snapshot
114 >        info_->getSnapshotManager()->advance();
115 >
116 >        //increase time
117 >        currentSnapshot_->increaseTime(dt);        
118 >
119 >        //save statistics
120 >        thermo.saveStat();
121          
122 <        info_->incrTime(fullStep_);
123 <        
124 <        //notify after integratreStep        
104 <        notify();
105 <        
106 <        if (info_->setTemp) {
107 <            if (info_->getTime() >= currThermal) {
108 <                thermalize();
122 >        if (needVelocityScaling) {
123 >            if (currentSnapshot_->getTime() >= currThermal) {
124 >                velocitizer_->velocitize(targetScalingTemp);
125                  currThermal += thermalTime;
126              }
127          }
128  
129 <        if (info_->getTime() >= currSample) {
130 <            dumpOut->writeDump(info_->getTime());
129 >        if (currentSnapshot_->getTime() >= currSample) {
130 >            dumpWriter->writeDump();
131              currSample += sampleTime;
132          }
133  
134 <        if (info_->getTime() >= currStatus) {
135 <            statOut->writeStat(info_->getTime());
136 <            calcPot = 0;
137 <            calcStress = 0;
134 >        if (currentSnapshot_->getTime() >= currStatus) {
135 >            statWriter->writeStat(currentSnapshot_->statData);
136 >            needPotential = false;
137 >            needStress = false;
138              currStatus += statusTime;
139          }
140  
141 <        if (info_->resetIntegrator) {
142 <            if (info_->getTime() >= currReset) {
127 <                this->resetIntegrator();
128 <                currReset += resetTime;
129 <            }
130 <        }
141 >        
142 > }
143  
132 #ifdef IS_MPI
144  
145 <        strcpy(checkPointMsg, "successfully took a time step.");
135 <        MPIcheckPoint();
145 > void VelocityVerletIntegrator::finalize() {
146  
147 < #endif // is_mpi
147 >    dumpWriter->writeDump();
148  
149 <    }
149 >    delete dumpWriter;
150 >    delete statWriter;
151  
152 <    dumpOut->writeFinal(info_->getTime());
153 <
154 <    delete dumpOut;
144 <    delete statOut;
152 >    dumpWriter = NULL;
153 >    statWriter = NULL;
154 >    
155   }
156  
157 < void VelocityVerletIntegrator::integrateStep() { }
157 > void VelocityVerletIntegrator::integrateStep() {
158  
159 +    moveA();
160 +    calcForce(needPotential, needStress);
161 +    moveB();
162 + }
163  
150 void VelocityVerletIntegrator::thermalize() {
164  
165 <    if (!info__->have_target_temp) {
166 <        sprintf(painCave.errMsg,
167 <                "You can't resample the velocities without a targetTemp!\n");
168 <        painCave.isFatal = 1;
156 <        painCave.severity = OOPSE_ERROR;
157 <        simError();
158 <        return;
159 <    }
165 > void VelocityVerletIntegrator::calcForce(bool needPotential,
166 >                                         bool needStress) {
167 >    forceMan_->calcForces(needPotential, needStress);
168 > }
169  
170 + DumpWriter* VelocityVerletIntegrator::createDumpWriter() {
171 +    return new DumpWriter(info_, info_->getDumpFileName());
172   }
173  
174 < void VelocityVerletIntegrator::calcForce(bool needPotential,
175 <                                         bool needStress) { }
174 > StatWriter* VelocityVerletIntegrator::createStatWriter() {
175 >    return new StatWriter(info_->getStatFileName());
176 > }
177  
178  
179   } //end namespace oopse

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines