--- trunk/src/brains/Snapshot.cpp	2013/06/16 15:15:42	1879
+++ trunk/src/brains/Snapshot.cpp	2014/10/16 19:13:51	2024
@@ -70,6 +70,7 @@ namespace OpenMD {
     frameData.torsionPotential = 0.0;   
     frameData.inversionPotential = 0.0; 
     frameData.lrPotentials = potVec(0.0);
+    frameData.reciprocalPotential = 0.0;
     frameData.excludedPotentials = potVec(0.0); 
     frameData.restraintPotential = 0.0; 
     frameData.rawPotential = 0.0;   
@@ -103,6 +104,7 @@ namespace OpenMD {
     frameData.torsionPotential = 0.0;   
     frameData.inversionPotential = 0.0; 
     frameData.lrPotentials = potVec(0.0);
+    frameData.reciprocalPotential = 0.0;
     frameData.excludedPotentials = potVec(0.0); 
     frameData.restraintPotential = 0.0; 
     frameData.rawPotential = 0.0;       
@@ -128,7 +130,8 @@ namespace OpenMD {
     frameData.pressure = 0.0;        
     frameData.temperature = 0.0;
     frameData.pressureTensor = Mat3x3d(0.0);   
-    frameData.systemDipole = Vector3d(0.0);            
+    frameData.systemDipole = Vector3d(0.0);
+    frameData.systemQuadrupole = Mat3x3d(0.0);
     frameData.convectiveHeatFlux = Vector3d(0.0, 0.0, 0.0);
     frameData.electronicTemperature = 0.0;
     frameData.COM = V3Zero;             
@@ -152,6 +155,7 @@ namespace OpenMD {
     hasCOMw = false;
     hasPressureTensor = false;    
     hasSystemDipole = false;      
+    hasSystemQuadrupole = false;
     hasConvectiveHeatFlux = false;  
     hasInertiaTensor = false;
     hasGyrationalVolume = false;   
@@ -187,6 +191,11 @@ namespace OpenMD {
   /** Returns the number of rigid bodies */
   int Snapshot::getNumberOfCutoffGroups() {
     return cgData.getSize();
+  }
+
+    /** Returns the number of bytes in a FrameData structure */
+  int Snapshot::getFrameDataSize() {
+    return sizeof(FrameData);
   }
   
   /** Returns the H-Matrix */
@@ -297,22 +306,24 @@ namespace OpenMD {
     frameData.volume = vol;
   }
 
+
   /** Wrap a vector according to periodic boundary conditions */
   void Snapshot::wrapVector(Vector3d& pos) {
-    
-    Vector3d scaled = scaleVector(pos);
     
-    for (int i = 0; i < 3; i++) 
-      scaled[i] -= roundMe(scaled[i]);
-    
-    if( !frameData.orthoRhombic )
-      pos = frameData.hmat * scaled;    
-    else {
-      
+    if( !frameData.orthoRhombic ) {
+      Vector3d scaled = frameData.invHmat * pos;
+      for (int i = 0; i < 3; i++) {
+        scaled[i] -= roundMe( scaled[i] );        
+      }
       // calc the wrapped real coordinates from the wrapped scaled coordinates
-      for (int i=0; i<3; i++) {
-	pos[i] = scaled[i] * frameData.hmat(i, i);
-      }   
+      pos = frameData.hmat * scaled;
+    } else {
+      RealType scaled;
+      for (int i=0; i<3; i++) {      
+        scaled = pos[i] * frameData.invHmat(i,i);
+        scaled -= roundMe( scaled );
+        pos[i] = scaled * frameData.hmat(i,i);
+      }
     }
   }
 
@@ -412,6 +423,14 @@ namespace OpenMD {
     return frameData.shortRangePotential;
   }
 
+  void Snapshot::setReciprocalPotential(RealType rp){
+    frameData.reciprocalPotential = rp;
+  }
+
+  RealType Snapshot::getReciprocalPotential() {
+    return frameData.reciprocalPotential;
+  }
+
   void Snapshot::setLongRangePotential(potVec lrPot) {
     frameData.lrPotentials = lrPot;
   }
@@ -421,6 +440,7 @@ namespace OpenMD {
       for (int i = 0; i < N_INTERACTION_FAMILIES; i++) {
         frameData.longRangePotential += frameData.lrPotentials[i];
       }
+      frameData.longRangePotential += frameData.reciprocalPotential;
       hasLongRangePotential = true;
     }   
     return frameData.longRangePotential;
@@ -434,6 +454,7 @@ namespace OpenMD {
     if (!hasPotentialEnergy) {
       frameData.potentialEnergy = this->getLongRangePotential();
       frameData.potentialEnergy += this->getShortRangePotential();
+      frameData.potentialEnergy += this->getRestraintPotential();
       hasPotentialEnergy = true;
     }
     return frameData.potentialEnergy;
@@ -584,6 +605,15 @@ namespace OpenMD {
     frameData.systemDipole = bd;
   }
 
+  Mat3x3d Snapshot::getSystemQuadrupole() {
+    return frameData.systemQuadrupole;
+  }
+
+  void Snapshot::setSystemQuadrupole(const Mat3x3d& bq) {    
+    hasSystemQuadrupole = true;
+    frameData.systemQuadrupole = bq;
+  }
+
   void Snapshot::setThermostat(const pair<RealType, RealType>& thermostat) {
     frameData.thermostat = thermostat;
   }