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

Comparing trunk/OOPSE/libmdtools/Symplectic.cpp (file contents):
Revision 468 by gezelter, Mon Apr 7 16:56:38 2003 UTC vs.
Revision 486 by mmeineke, Thu Apr 10 16:22:00 2003 UTC

# Line 162 | Line 162 | void Symplectic::integrate( void ){
162    double dt2;                       // half the dt
163  
164    double vx, vy, vz;    // the velocities
165 < //  double vx2, vy2, vz2; // the square of the velocities
165 >  double vx2, vy2, vz2; // the square of the velocities
166    double rx, ry, rz;    // the postitions
167    
168    double ji[3];   // the body frame angular momentum
# Line 170 | Line 170 | void Symplectic::integrate( void ){
170    double Tb[3];   // torque in the body frame
171    double angle;   // the angle through which to rotate the rotation matrix
172    double A[3][3]; // the rotation matrix
173 +  double press[9];
174  
175    int time;
176  
# Line 186 | Line 187 | void Symplectic::integrate( void ){
187  
188    int calcPot, calcStress;
189  
190 <   Thermo *tStats = new Thermo( entry_plug );
190 >  Thermo *tStats;
191 >  StatWriter*  e_out;
192 >  DumpWriter*  dump_out;
193  
194 <  StatWriter*  e_out    = new StatWriter( entry_plug );
195 <  DumpWriter*  dump_out = new DumpWriter( entry_plug );
194 >  tStats   = new Thermo( entry_plug );
195 >  e_out    = new StatWriter( entry_plug );
196 >  dump_out = new DumpWriter( entry_plug );
197  
198    Atom** atoms = entry_plug->atoms;
199    DirectionalAtom* dAtom;
# Line 208 | Line 212 | void Symplectic::integrate( void ){
212    e_out->writeStat( 0.0 );
213    
214    calcPot = 0;
215 +
216 +  if (!strcasecmp( entry_plug->ensemble, "NPT")) {
217 +    calcStress = 1;
218 +  } else {
219 +    calcStress = 0;
220 +  }
221  
222    if( n_constrained ){
223  
# Line 225 | Line 235 | void Symplectic::integrate( void ){
235      
236  
237      for( tl=0; tl < n_loops; tl++ ){
238 +
239 +      if (!strcasecmp( entry_plug->ensemble, "NVT"))
240 +        myES->NoseHooverNVT( dt / 2.0 , tStats->getKinetic() );
241        
242        for( j=0; j<nAtoms; j++ ){
243  
# Line 390 | Line 403 | void Symplectic::integrate( void ){
403          }
404        }
405      
406 +
407 +      if (!strcasecmp( entry_plug->ensemble, "NVT"))
408 +        myES->NoseHooverNVT( dt / 2.0, tStats->getKinetic() );
409 +
410 +      if (!strcasecmp( entry_plug->ensemble, "NPT") ) {
411 +        tStats->getPressureTensor(press);
412 +        myES->NoseHooverAndersonNPT( dt,
413 +                                     tStats->getKinetic(),
414 +                                     press);
415 +      }
416 +
417        time = tl + 1;
418        
419        if( entry_plug->setTemp ){
# Line 398 | Line 422 | void Symplectic::integrate( void ){
422        if( !(time % sample_n) ) dump_out->writeDump( time * dt );
423        if( !((time+1) % status_n) ) {
424          calcPot = 1;
425 <        calcStress = 1;
425 >        calcStress = 1;
426        }
427        if( !(time % status_n) ){
428          e_out->writeStat( time * dt );
429          calcPot = 0;
430 <        calcStress = 0;
430 >        if (!strcasecmp(entry_plug->ensemble, "NPT")) calcStress = 1;
431 >        else calcStress = 0;
432        }
433      }
434    }
# Line 414 | Line 439 | void Symplectic::integrate( void ){
439        kE = 0.0;
440        rot_kE= 0.0;
441        trans_kE = 0.0;
442 +
443 +      if (!strcasecmp( entry_plug->ensemble, "NVT"))
444 +        myES->NoseHooverNVT( dt / 2.0, tStats->getKinetic() );
445        
446        for( i=0; i<nAtoms; i++ ){
447          
# Line 519 | Line 547 | void Symplectic::integrate( void ){
547          atoms[i]->set_vy( vy );
548          atoms[i]->set_vz( vz );
549          
550 < //      vx2 = vx * vx;
551 < //      vy2 = vy * vy;
552 < //      vz2 = vz * vz;
550 >        vx2 = vx * vx;
551 >        vy2 = vy * vy;
552 >        vz2 = vz * vz;
553          
554          if( atoms[i]->isDirectional() ){
555  
# Line 553 | Line 581 | void Symplectic::integrate( void ){
581            dAtom->setJy( ji[1] );
582            dAtom->setJz( ji[2] );
583          }
584 +
585        }
586 <      
586 >    
587 >      if (!strcasecmp( entry_plug->ensemble, "NVT"))
588 >        myES->NoseHooverNVT( dt / 2.0, tStats->getKinetic() );
589 >
590 >      if (!strcasecmp( entry_plug->ensemble, "NPT") ) {
591 >        tStats->getPressureTensor(press);
592 >        myES->NoseHooverAndersonNPT( dt,
593 >                                     tStats->getKinetic(),
594 >                                     press);
595 >      }
596 >  
597        time = tl + 1;
598        
599        if( entry_plug->setTemp ){
# Line 563 | Line 602 | void Symplectic::integrate( void ){
602        if( !(time % sample_n) ) dump_out->writeDump( time * dt );
603        if( !((time+1) % status_n) ) {
604          calcPot = 1;
605 <        calcStress = 1;
605 >        calcStress = 1;
606        }
607        if( !(time % status_n) ){
608          e_out->writeStat( time * dt );
609          calcPot = 0;
610 <        calcStress = 0;
610 >        if (!strcasecmp(entry_plug->ensemble, "NPT")) calcStress = 1;
611 >        else calcStress = 0;
612        }      
613      }
614    }

Diff Legend

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