--- trunk/OOPSE/libmdtools/SimInfo.cpp	2003/09/19 14:55:41	770
+++ trunk/OOPSE/libmdtools/SimInfo.cpp	2003/09/22 23:07:57	781
@@ -48,6 +48,7 @@ SimInfo::SimInfo(){
   haveOrigEcr = 0;
   boxIsInit = 0;
   
+  resetTime = 1e99;
   
 
   usePBC = 0;
@@ -279,6 +280,23 @@ void SimInfo::printMat9(double A[9] ){
 	    << "[ " << A[0] << ", " << A[1] << ", " << A[2] << " ]\n"
 	    << "[ " << A[3] << ", " << A[4] << ", " << A[5] << " ]\n"
 	    << "[ " << A[6] << ", " << A[7] << ", " << A[8] << " ]\n";
+}
+
+
+void SimInfo::crossProduct3(double a[3],double b[3], double out[3]){
+
+      out[0] = a[1] * b[2] - a[2] * b[1];
+      out[1] = a[2] * b[0] - a[0] * b[2] ;
+      out[2] = a[0] * b[1] - a[1] * b[0];
+      
+}
+
+double SimInfo::dotProduct3(double a[3], double b[3]){
+  return a[0]*b[0] + a[1]*b[1]+ a[2]*b[2];
+}
+
+double SimInfo::length3(double a[3]){
+  return sqrt(a[0]*a[0] + a[1]*a[1] + a[2]*a[2]);
 }
 
 void SimInfo::calcBoxL( void ){
@@ -295,27 +313,63 @@ void SimInfo::calcBoxL( void ){
   dx = Hmat[0][0]; dy = Hmat[1][0]; dz = Hmat[2][0];
   dsq = dx*dx + dy*dy + dz*dz;
   boxL[0] = sqrt( dsq );
-  maxCutoff = 0.5 * boxL[0];
+  //maxCutoff = 0.5 * boxL[0];
 
   // boxLy
   
   dx = Hmat[0][1]; dy = Hmat[1][1]; dz = Hmat[2][1];
   dsq = dx*dx + dy*dy + dz*dz;
   boxL[1] = sqrt( dsq );
-  if( (0.5 * boxL[1]) < maxCutoff ) maxCutoff = 0.5 * boxL[1];
+  //if( (0.5 * boxL[1]) < maxCutoff ) maxCutoff = 0.5 * boxL[1];
 
+
   // boxLz
   
   dx = Hmat[0][2]; dy = Hmat[1][2]; dz = Hmat[2][2];
   dsq = dx*dx + dy*dy + dz*dz;
   boxL[2] = sqrt( dsq );
-  if( (0.5 * boxL[2]) < maxCutoff ) maxCutoff = 0.5 * boxL[2];
+  //if( (0.5 * boxL[2]) < maxCutoff ) maxCutoff = 0.5 * boxL[2];
+
+  //calculate the max cutoff
+  maxCutoff =  calcMaxCutOff(); 
   
   checkCutOffs();
 
 }
 
 
+double SimInfo::calcMaxCutOff(){
+
+  double ri[3], rj[3], rk[3];
+  double rij[3], rjk[3], rki[3];
+  double minDist;
+
+  ri[0] = Hmat[0][0];
+  ri[1] = Hmat[1][0];
+  ri[2] = Hmat[2][0];
+
+  rj[0] = Hmat[0][1];
+  rj[1] = Hmat[1][1];
+  rj[2] = Hmat[2][1];
+
+  rk[0] = Hmat[0][2];
+  rk[1] = Hmat[1][2];
+  rk[2] = Hmat[2][2];
+  
+  crossProduct3(ri,rj, rij);
+  distXY = dotProduct3(rk,rij) / length3(rij);
+
+  crossProduct3(rj,rk, rjk);
+  distYZ = dotProduct3(ri,rjk) / length3(rjk);
+
+  crossProduct3(rk,ri, rki);
+  distZX = dotProduct3(rj,rki) / length3(rki);
+
+  minDist = min(min(distXY, distYZ), distZX);
+  return minDist/2;
+  
+}
+
 void SimInfo::wrapVector( double thePos[3] ){
 
   int i, j, k;
@@ -497,47 +551,49 @@ void SimInfo::checkCutOffs( void ){
   if( boxIsInit ){
     
     //we need to check cutOffs against the box
-    
+
+    //detect the change of rCut
     if(( maxCutoff > rCut )&&(usePBC)){
       if( rCut < origRcut ){
-	rCut = origRcut;
-	if (rCut > maxCutoff) rCut = maxCutoff;
-	
-	sprintf( painCave.errMsg,
-		 "New Box size is setting the long range cutoff radius "
-		 "to %lf at time %lf\n",
-		 rCut, currentTime );
-	painCave.isFatal = 0;
-	simError();
+        rCut = origRcut;
+        
+        if (rCut > maxCutoff)
+          rCut = maxCutoff;
+  
+          sprintf( painCave.errMsg,
+                    "New Box size is setting the long range cutoff radius "
+                    "to %lf at time %lf\n",
+                    rCut, currentTime );
+          painCave.isFatal = 0;
+          simError();
       }
     }
-    
-    if( maxCutoff > ecr ){
-      if( ecr < origEcr ){
-	ecr = origEcr;
-	if (ecr > maxCutoff) ecr = maxCutoff;
-	
-	sprintf( painCave.errMsg,
-		 "New Box size is setting the electrostaticCutoffRadius "
-		 "to %lf at time %lf\n",
-		 ecr, currentTime );
-	painCave.isFatal = 0;
-	simError();
-      }
-    }
-    
-    
-    if ((rCut > maxCutoff)&&(usePBC)) {
+    else if ((rCut > maxCutoff)&&(usePBC)) {
       sprintf( painCave.errMsg,
                "New Box size is setting the long range cutoff radius "
-	       "to %lf at time %lf\n",
+               "to %lf at time %lf\n",
                maxCutoff, currentTime );
       painCave.isFatal = 0;
       simError();
       rCut = maxCutoff;
     }
-    
-    if( ecr > maxCutoff){
+
+
+    //detect the change of ecr
+    if( maxCutoff > ecr ){
+      if( ecr < origEcr ){
+        ecr = origEcr;
+        if (ecr > maxCutoff) ecr = maxCutoff;
+  
+          sprintf( painCave.errMsg,
+                    "New Box size is setting the electrostaticCutoffRadius "
+                    "to %lf at time %lf\n",
+                    ecr, currentTime );
+            painCave.isFatal = 0;
+            simError();
+      }
+    }
+    else if( ecr > maxCutoff){
       sprintf( painCave.errMsg,
                "New Box size is setting the electrostaticCutoffRadius "
                "to %lf at time %lf\n",