--- trunk/OOPSE/libmdtools/NPTf.cpp	2003/07/14 22:38:13	600
+++ trunk/OOPSE/libmdtools/NPTf.cpp	2003/07/15 19:56:08	617
@@ -1,3 +1,4 @@
+#include <cmath>
 #include "Atom.hpp"
 #include "SRI.hpp"
 #include "AbstractClasses.hpp"
@@ -50,6 +51,7 @@ void NPTf::moveA() {
   double sc[3];
   double eta2ij;
   double press[3][3], vScale[3][3], hm[3][3], hmnew[3][3], scaleMat[3][3];
+  double bigScale, smallScale, offDiagMax;
 
   tt2 = tauThermostat * tauThermostat;
   tb2 = tauBarostat * tauBarostat;
@@ -108,6 +110,8 @@ void NPTf::moveA() {
 
     for (j = 0; j < 3; j++ ) 
       pos[j] += dt * (vel[j] + sc[j]);
+
+    atoms[i]->setPos( pos );
    
     if( atoms[i]->isDirectional() ){
 
@@ -161,6 +165,9 @@ void NPTf::moveA() {
   // Use a taylor expansion for eta products:  Hmat = Hmat . exp(dt * etaMat)
   //  Hmat = Hmat . ( Ident + dt * etaMat  + dt^2 * etaMat*etaMat / 2)
   
+  bigScale = 1.0;
+  smallScale = 1.0;
+  offDiagMax = 0.0;
   
   for(i=0; i<3; i++){
     for(j=0; j<3; j++){
@@ -178,13 +185,46 @@ void NPTf::moveA() {
       if (i == j) scaleMat[i][j] = 1.0;
       // Taylor expansion for the exponential truncated at second order:
       scaleMat[i][j] += dt*eta[i][j]  + 0.5*dt*dt*eta2ij;
+
+      if (i != j)
+        if (fabs(scaleMat[i][j]) > offDiagMax) 
+          offDiagMax = fabs(scaleMat[i][j]);
       
     }
+
+    if (scaleMat[i][i] > bigScale) bigScale = scaleMat[i][i];
+    if (scaleMat[i][i] < smallScale) smallScale = scaleMat[i][i];
   }
   
-  info->getBoxM(hm);
-  info->matMul3(hm, scaleMat, hmnew);
-  info->setBoxM(hmnew);
+  if ((bigScale > 1.1) || (smallScale < 0.9)) {
+    sprintf( painCave.errMsg,
+             "NPTf error: Attempting a Box scaling of more than 10 percent.\n"
+             " Check your tauBarostat, as it is probably too small!\n\n"
+             " scaleMat = [%lf\t%lf\t%lf]\n"
+             "            [%lf\t%lf\t%lf]\n"
+             "            [%lf\t%lf\t%lf]\n",
+             scaleMat[0][0],scaleMat[0][1],scaleMat[0][2],
+             scaleMat[1][0],scaleMat[1][1],scaleMat[1][2],
+             scaleMat[2][0],scaleMat[2][1],scaleMat[2][2]);
+    painCave.isFatal = 1;
+    simError();
+  } else if (offDiagMax > 0.1) {
+    sprintf( painCave.errMsg,
+             "NPTf error: Attempting an off-diagonal Box scaling of more than 10 percent.\n"
+             " Check your tauBarostat, as it is probably too small!\n\n"
+             " scaleMat = [%lf\t%lf\t%lf]\n"
+             "            [%lf\t%lf\t%lf]\n"
+             "            [%lf\t%lf\t%lf]\n",
+             scaleMat[0][0],scaleMat[0][1],scaleMat[0][2],
+             scaleMat[1][0],scaleMat[1][1],scaleMat[1][2],
+             scaleMat[2][0],scaleMat[2][1],scaleMat[2][2]);
+    painCave.isFatal = 1;
+    simError();
+  } else {
+    info->getBoxM(hm);
+    info->matMul3(hm, scaleMat, hmnew);
+    info->setBoxM(hmnew);
+  }
   
 }