--- trunk/OOPSE/libmdtools/SimInfo.cpp	2003/06/30 22:04:01	568
+++ trunk/OOPSE/libmdtools/SimInfo.cpp	2003/07/02 21:26:55	572
@@ -2,6 +2,8 @@
 #include <cstring>
 #include <cmath>
 
+#include <iostream>
+using namespace std;
 
 #include "SimInfo.hpp"
 #define __C
@@ -14,6 +16,11 @@
 #include "mpiSimulation.hpp"
 #endif
 
+inline double roundMe( double x ){
+  return ( x >= 0 ) ? floor( x + 0.5 ) : ceil( x - 0.5 );
+}
+	  
+
 SimInfo* currentInfo;
 
 SimInfo::SimInfo(){
@@ -54,7 +61,7 @@ void SimInfo::setBox(double newBox[3]) {
   calcHmatI();
   calcBoxL();
 
-  setFortranBoxSize(Hmat);
+  setFortranBoxSize(Hmat, HmatI, &orthoRhombic);
 
   smallestBoxL = boxLx;
   if (boxLy < smallestBoxL) smallestBoxL = boxLy;
@@ -110,7 +117,7 @@ void SimInfo::setBoxM( double theBox[9] ){
   calcHmatI();
   calcBoxL();
 
-  setFortranBoxSize(Hmat);
+  setFortranBoxSize(Hmat, HmatI, &orthoRhombic);
  
   smallestBoxL = boxLx;
   if (boxLy < smallestBoxL) smallestBoxL = boxLy;
@@ -158,7 +165,7 @@ void SimInfo::setBoxM( double theBox[9] ){
 }
  
 
-void SimInfo::getBox(double theBox[9]) {
+void SimInfo::getBoxM (double theBox[9]) {
 
   int i;
   for(i=0; i<9; i++) theBox[i] = Hmat[i];
@@ -170,6 +177,9 @@ void SimInfo::calcHmatI( void ) {
   double C[3][3];
   double detHmat;
   int i, j, k;
+  double smallDiag;
+  double tol;
+  double sanity[3][3];
 
   // calculate the adjunct of Hmat;
 
@@ -201,6 +211,41 @@ void SimInfo::calcHmatI( void ) {
       i++;
     }
   }
+
+  // sanity check
+
+  for(i=0; i<3; i++){
+    for(j=0; j<3; j++){
+      
+      sanity[i][j] = 0.0;
+      for(k=0; k<3; k++){
+	sanity[i][j] += Hmat[3*k+i] * HmatI[3*j+k];
+      }
+    }
+  }
+
+  cerr << "sanity => \n" 
+       << sanity[0][0] << "\t" << sanity[0][1] << "\t" << sanity [0][2] << "\n"
+       << sanity[1][0] << "\t" << sanity[1][1] << "\t" << sanity [1][2] << "\n"
+       << sanity[2][0] << "\t" << sanity[2][1] << "\t" << sanity [2][2] 
+       << "\n";
+    
+
+  // check to see if Hmat is orthorhombic
+  
+  smallDiag = Hmat[0];
+  if(smallDiag > Hmat[4]) smallDiag = Hmat[4];
+  if(smallDiag > Hmat[8]) smallDiag = Hmat[8];
+  tol = smallDiag * 1E-6;
+
+  orthoRhombic = 1;
+  for(i=0; (i<9) && orthoRhombic; i++){
+    
+    if( (i%4) ){ // ignore the diagonals (0, 4, and 8)
+      orthoRhombic = (Hmat[i] <= tol);
+    }
+  }
+    
 }
 
 void SimInfo::calcBoxL( void ){
@@ -241,17 +286,42 @@ void SimInfo::wrapVector( double thePos[3] ){
   int i, j, k;
   double scaled[3];
 
-  // calc the scaled coordinates.
-  
-  for(i=0; i<3; i++)
-    scaled[i] = thePos[0]*Hmat[i] + thePos[1]*Hat[i+3] + thePos[3]*Hmat[i+6];
-
-  // wrap the scaled coordinates
-
-  for(i=0; i<3; i++)
-    scaled[i] -= (copysign(1,scaled[i]) * (int)(fabs(scaled[i]) + 0.5));
-  
-
+  if( !orthoRhombic ){
+    // calc the scaled coordinates.
+    
+    for(i=0; i<3; i++)
+      scaled[i] = 
+	thePos[0]*HmatI[i] + thePos[1]*HmatI[i+3] + thePos[3]*HmatI[i+6];
+    
+    // wrap the scaled coordinates
+    
+    for(i=0; i<3; i++)
+      scaled[i] -= roundMe(scaled[i]);
+    
+    // calc the wrapped real coordinates from the wrapped scaled coordinates
+    
+    for(i=0; i<3; i++)
+      thePos[i] = 
+	scaled[0]*Hmat[i] + scaled[1]*Hmat[i+3] + scaled[2]*Hmat[i+6];
+  }
+  else{
+    // calc the scaled coordinates.
+    
+    for(i=0; i<3; i++)
+      scaled[i] = thePos[i]*HmatI[i*4];
+    
+    // wrap the scaled coordinates
+    
+    for(i=0; i<3; i++)
+      scaled[i] -= roundMe(scaled[i]);
+    
+    // calc the wrapped real coordinates from the wrapped scaled coordinates
+    
+    for(i=0; i<3; i++)
+      thePos[i] = scaled[i]*Hmat[i*4];
+  }
+    
+    
 }
 
 
@@ -297,10 +367,6 @@ void SimInfo::refreshSim(){
   fInfo.rt = 0.0;
   fInfo.dielect = 0.0;
 
-  fInfo.box[0] = box_x;
-  fInfo.box[1] = box_y;
-  fInfo.box[2] = box_z;
-
   fInfo.rlist = rList;
   fInfo.rcut = rCut;