--- trunk/OOPSE/libmdtools/NPTi.cpp	2003/07/08 21:06:14	575
+++ trunk/OOPSE/libmdtools/NPTi.cpp	2003/07/09 22:14:06	586
@@ -1,3 +1,4 @@
+#include <cmath>
 #include "Atom.hpp"
 #include "SRI.hpp"
 #include "AbstractClasses.hpp"
@@ -42,6 +43,7 @@ void NPTi::moveA() {
   double tt2, tb2;
   double angle;
 
+
   tt2 = tauThermostat * tauThermostat;
   tb2 = tauBarostat * tauBarostat;
 
@@ -49,10 +51,11 @@ void NPTi::moveA() {
   instaPress = tStats->getPressure();
   instaVol = tStats->getVolume();
    
-  // first evolve chi a half step
+   // first evolve chi a half step
   
   chi += dt2 * ( instaTemp / targetTemp - 1.0) / tt2;
-  eta += dt2 * ( instaVol * (instaPress - targetPressure) / (NkBT*tb2));
+  eta += dt2 * ( instaVol * (instaPress - targetPressure) / 
+		 (p_convert*NkBT*tb2));
 
   for( i=0; i<nAtoms; i++ ){
     atomIndex = i * 3;
@@ -65,21 +68,15 @@ void NPTi::moveA() {
 
     // position whole step    
 
-    for( j=atomIndex; j<(atomIndex+3); j=j+3 ) {
-      rj[0] = pos[j];
-      rj[1] = pos[j+1];
-      rj[2] = pos[j+2];
+    rj[0] = pos[atomIndex];
+    rj[1] = pos[atomIndex+1];
+    rj[2] = pos[atomIndex+2];
+    
+    info->wrapVector(rj);
 
-      info->wrapVector(rj);
-
-      pos[j] += dt * (vel[j] + eta*rj[0]);
-      pos[j+1] += dt * (vel[j+1] + eta*rj[1]);
-      pos[j+2] += dt * (vel[j+2] + eta*rj[2]);
-    }
-
-    // Scale the box after all the positions have been moved:
-
-    info->scaleBox(exp(dt*eta));
+    pos[atomIndex] += dt * (vel[atomIndex] + eta*rj[0]);
+    pos[atomIndex+1] += dt * (vel[atomIndex+1] + eta*rj[1]);
+    pos[atomIndex+2] += dt * (vel[atomIndex+2] + eta*rj[2]);
    
     if( atoms[i]->isDirectional() ){
 
@@ -132,6 +129,13 @@ void NPTi::moveA() {
     }
     
   }
+  // Scale the box after all the positions have been moved:
+
+  cerr << "eta = " << eta 
+       << "; exp(dt*eta) = " << exp(eta*dt) << "\n";
+
+  info->scaleBox(exp(dt*eta));
+
 }
 
 void NPTi::moveB( void ){
@@ -151,7 +155,8 @@ void NPTi::moveB( void ){
   instaVol = tStats->getVolume();
 
   chi += dt2 * ( instaTemp / targetTemp - 1.0) / tt2;
-  eta += dt2 * ( instaVol * (instaPress - targetPressure) / (NkBT*tb2));
+  eta += dt2 * ( instaVol * (instaPress - targetPressure) / 
+		 (p_convert*NkBT*tb2));
   
   for( i=0; i<nAtoms; i++ ){
     atomIndex = i * 3;