ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/Verlet.cpp
(Generate patch)

Comparing trunk/OOPSE/libmdtools/Verlet.cpp (file contents):
Revision 468 by gezelter, Mon Apr 7 16:56:38 2003 UTC vs.
Revision 478 by gezelter, Tue Apr 8 14:39:40 2003 UTC

# Line 202 | Line 202 | void Verlet::integrate( void ){
202    DumpWriter*  dump_out = new DumpWriter( entry_plug );
203  
204    // the first time integrate is called, the forces need to be initialized
205 –
205    
206    myFF->doForces(1,1);
207    
# Line 215 | Line 214 | void Verlet::integrate( void ){
214    e_out->writeStat( 0.0 );
215  
216    calcPot = 0;
218 –  calcStress = 0;
217  
218 +  if (!strcasecmp( entry_plug->ensemble, "NPT")) {
219 +    calcStress = 1;
220 +  } else {
221 +    calcStress = 0;
222 +  }
223 +
224    if( c_is_constrained ){
225      for(i = 0; i < n_loops; i++){
226        
227 +      if (!strcasecmp( entry_plug->ensemble, "NVT"))
228 +        myES->NoseHooverNVT( dt / 2.0 , tStats->getKinetic() );
229 +      
230        // fill R, V, and F arrays and RATTLE in fortran
231 <
231 >      
232        for( j=0; j<c_natoms; j++ ){
233 <
233 >        
234          Rx[j] = c_atoms[j]->getX();
235          Ry[j] = c_atoms[j]->getY();
236          Rz[j] = c_atoms[j]->getZ();
# Line 276 | Line 283 | void Verlet::integrate( void ){
283          Fz[j] = c_atoms[j]->getFz();
284        }
285          
286 +
287        v_constrain_b_( dt, c_natoms, c_mass, Rx, Ry, Rz, Vx, Vy, Vz,
288                        Fx, Fy, Fz,
289                        kE, c_n_constrained, c_constrained_dsqr,
# Line 293 | Line 301 | void Verlet::integrate( void ){
301          c_atoms[j]->set_vz(Vz[j]);
302        }
303        
304 +      if (!strcasecmp( entry_plug->ensemble, "NVT"))
305 +        myES->NoseHooverNVT( dt / 2.0, tStats->getKinetic() );
306 +      
307 +      if (!strcasecmp( entry_plug->ensemble, "NPT") )
308 +        myES->NoseHooverAndersonNPT( dt,
309 +                                     tStats->getKinetic(),
310 +                                     tStats->getPressure());
311 +
312        time = i + 1;
313        
314        if( entry_plug->setTemp ){
# Line 301 | Line 317 | void Verlet::integrate( void ){
317        if( !(time % sample_n) ) dump_out->writeDump( time * dt );
318        if( !((time+1) % status_n) ) {
319          calcPot = 1;
320 <        calcStress = 1;
320 >        // bitwise masking in case we need it for NPT
321 >        calcStress = (!strcasecmp(entry_plug->ensemble,"NPT")) && 1;
322        }
323        if( !(time % status_n) ){
324          e_out->writeStat( time * dt );
325          calcPot = 0;
326 <        calcStress = 0;
326 >        // bitwise masking in case we need it for NPT
327 >        calcStress = (!strcasecmp(entry_plug->ensemble,"NPT")) && 0;
328        }
329      }
330    }
331    else{
332      for(i = 0; i < n_loops; i++){
333 <      
333 >
334 >      if (!strcasecmp( entry_plug->ensemble, "NVT"))
335 >        myES->NoseHooverNVT( dt / 2.0, tStats->getKinetic() );
336 >    
337        move_a( dt );
338        
339        // calculate the forces
# Line 323 | Line 344 | void Verlet::integrate( void ){
344        
345        move_b( dt );
346  
347 +      if (!strcasecmp( entry_plug->ensemble, "NVT"))
348 +        myES->NoseHooverNVT( dt / 2.0 , tStats->getKinetic() );
349 +
350 +      if (!strcasecmp( entry_plug->ensemble, "NPT") )
351 +        myES->NoseHooverAndersonNPT( dt,
352 +                                     tStats->getKinetic(),
353 +                                     tStats->getPressure());
354 +
355        time = i + 1;
356        
357        if( entry_plug->setTemp ){
# Line 331 | Line 360 | void Verlet::integrate( void ){
360        if( !(time % sample_n) )  dump_out->writeDump( time * dt );
361        if( !((time+1) % status_n) ) {
362          calcPot = 1;
363 <        calcStress = 1;
363 >        // bitwise masking in case we need it for NPT
364 >        calcStress = (!strcasecmp(entry_plug->ensemble,"NPT")) && 1;
365        }
366        if( !(time % status_n) ){
367          e_out->writeStat( time * dt );
368          calcPot = 0;
369 <        calcStress = 0;
369 >        // bitwise masking in case we need it for NPT
370 >        calcStress = (!strcasecmp(entry_plug->ensemble,"NPT")) && 0;
371        }
372      }
373    }

Diff Legend

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