--- trunk/src/primitives/StuntDouble.hpp	2005/08/10 19:17:39	573
+++ trunk/src/primitives/StuntDouble.hpp	2013/07/19 21:25:45	1908
@@ -6,19 +6,10 @@
  * redistribute this software in source and binary code form, provided
  * that the following conditions are met:
  *
- * 1. Acknowledgement of the program authors must be made in any
- *    publication of scientific results based in part on use of the
- *    program.  An acceptable form of acknowledgement is citation of
- *    the article in which the program was described (Matthew
- *    A. Meineke, Charles F. Vardeman II, Teng Lin, Christopher
- *    J. Fennell and J. Daniel Gezelter, "OOPSE: An Object-Oriented
- *    Parallel Simulation Engine for Molecular Dynamics,"
- *    J. Comput. Chem. 26, pp. 252-271 (2005))
- *
- * 2. Redistributions of source code must retain the above copyright
+ * 1. Redistributions of source code must retain the above copyright
  *    notice, this list of conditions and the following disclaimer.
  *
- * 3. Redistributions in binary form must reproduce the above copyright
+ * 2. Redistributions in binary form must reproduce the above copyright
  *    notice, this list of conditions and the following disclaimer in the
  *    documentation and/or other materials provided with the
  *    distribution.
@@ -37,6 +28,16 @@
  * arising out of the use of or inability to use software, even if the
  * University of Notre Dame has been advised of the possibility of
  * such damages.
+ *
+ * SUPPORT OPEN SCIENCE!  If you use OpenMD or its source code in your
+ * research, please cite the appropriate papers when you publish your
+ * work.  Good starting points are:
+ *                                                                      
+ * [1]  Meineke, et al., J. Comp. Chem. 26, 252-271 (2005).             
+ * [2]  Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006).          
+ * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008).          
+ * [4]  Kuang & Gezelter,  J. Chem. Phys. 133, 164101 (2010).
+ * [5]  Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011).
  */
  
 /**
@@ -58,7 +59,7 @@
 #include "utils/PropertyMap.hpp"
 #include "brains/Snapshot.hpp"
 #include "brains/SnapshotManager.hpp"
-namespace oopse{
+namespace OpenMD{
 
 
    
@@ -83,7 +84,7 @@ namespace oopse{
    * doesn't try to do something stupid like torque an Atom (The
    * quotes above are from Spaceballs...)
    *
-   * @note the dynamic data of stuntdouble will be stored outside of the class
+   * @note the dynamic data of stuntDouble will be stored outside of the class
    */
   class StuntDouble{
   public:    
@@ -97,95 +98,125 @@ namespace oopse{
     virtual ~StuntDouble();
         
     /**
-     * Returns the global index of this stuntdouble.
-     * @return  the global index of this stuntdouble 
+     * Returns the global index of this stuntDouble.
+     * @return  the global index of this stuntDouble 
      */
     int getGlobalIndex() {
       return globalIndex_;
     }
 
     /**
-     * Sets the global index of this stuntdouble.
-     * @param new global index to be set
+     * Sets the global index of this stuntDouble.
+     * @param index new global index to be set
      */
     void setGlobalIndex(int index) {
       globalIndex_ = index;
     }
-        
+    
     /** 
-     * Returns the local index of this stuntdouble 
-     * @return the local index of this stuntdouble
+     * Returns the local index of this stuntDouble 
+     * @return the local index of this stuntDouble
      */
     int getLocalIndex() {
       return localIndex_;
     }
 
     /**
-     * Sets the local index of this stuntdouble
+     * Sets the local index of this stuntDouble
      * @param index new index to be set
      */        
     void setLocalIndex(int index) {
       localIndex_ = index;
     }
+    
+    int getGlobalIntegrableObjectIndex(){
+      return globalIntegrableObjectIndex_; 
+    }
+    void setGlobalIntegrableObjectIndex(int index) {
+      globalIntegrableObjectIndex_ = index;
+    }
 
     /**
-     * Sets the Snapshot Manager of this stuntdouble
+     * Sets the Snapshot Manager of this stuntDouble
      */
     void setSnapshotManager(SnapshotManager* sman) {
       snapshotMan_ = sman;
     }
 
     /**
-     * Tests if this stuntdouble is an atom 
-     * @return true is this stuntdouble is an atom(or a directional atom), return false otherwise
+     * Tests if this stuntDouble is an atom 
+     * @return true is this stuntDouble is an atom(or a directional atom), return false otherwise
      */
     bool isAtom(){
       return objType_ == otAtom || objType_ == otDAtom;
     }
 
     /** 
-     * Tests if this stuntdouble is an directional atom 
-     * @return true if this stuntdouble is an directional atom, return false otherwise
+     * Tests if this stuntDouble is an directional atom 
+     * @return true if this stuntDouble is an directional atom, return false otherwise
      */
     bool isDirectionalAtom(){
       return objType_ == otDAtom;
     }
 
     /**
-     * Tests if this stuntdouble is a rigid body. 
-     * @return true if this stuntdouble is a rigid body, otherwise return false
+     * Tests if this stuntDouble is a rigid body. 
+     * @return true if this stuntDouble is a rigid body, otherwise return false
      */
     bool isRigidBody(){
       return objType_ == otRigidBody;
     }
 
     /**
-     * Tests if this stuntdouble is a directional one. 
-     * @return true is this stuntdouble is a directional atom or a rigid body, return false otherwise
+     * Tests if this stuntDouble is a directional one. 
+     * @return true is this stuntDouble is a directional atom or a rigid body, return false otherwise
      */
     bool isDirectional(){
       return isDirectionalAtom() || isRigidBody();
     }
 
     /**
-     * Returns the previous position of this stuntdouble
-     * @return the position of this stuntdouble
+     * Freezes out all velocity, angular velocity, forces and torques
+     * on this StuntDouble.  Also computes the number of frozen degrees
+     * of freedom.
+     * @return the total number of frozen degrees of freedom
+     */   
+    int freeze() {
+      
+      int fdf = 3;
+
+      setVel(V3Zero);
+      setFrc(V3Zero);
+      if (isDirectional()){
+        setJ(V3Zero);
+        setTrq(V3Zero);
+        if (isLinear()) 
+          fdf +=2;
+        else 
+          fdf +=3;        
+      }      
+      return fdf;
+    }
+
+    /**
+     * Returns the previous position of this stuntDouble
+     * @return the position of this stuntDouble
      */    
     Vector3d getPrevPos() {
       return ((snapshotMan_->getPrevSnapshot())->*storage_).position[localIndex_];
     }
        
     /**
-     * Returns the current position of this stuntdouble
-     * @return the position of this stuntdouble
+     * Returns the current position of this stuntDouble
+     * @return the position of this stuntDouble
      */    
     Vector3d getPos() {
       return ((snapshotMan_->getCurrentSnapshot())->*storage_).position[localIndex_];
     }
 
     /**
-     * Returns the position of this stuntdouble in specified snapshot 
-     * @return the position of this stuntdouble
+     * Returns the position of this stuntDouble in specified snapshot 
+     * @return the position of this stuntDouble
      * @param snapshotNo
      */    
     Vector3d getPos(int snapshotNo) {
@@ -193,7 +224,7 @@ namespace oopse{
     }
 
     /**
-     * Sets  the previous position of this stuntdouble
+     * Sets  the previous position of this stuntDouble
      * @param pos  new position 
      * @see #getPos
      */         
@@ -202,7 +233,7 @@ namespace oopse{
     }
        
     /**
-     * Sets  the current position of this stuntdouble
+     * Sets  the current position of this stuntDouble
      * @param pos  new position 
      */         
     void setPos(const Vector3d& pos) {
@@ -212,7 +243,7 @@ namespace oopse{
     }
 
     /**
-     * Sets  the position of this stuntdouble in specified snapshot
+     * Sets  the position of this stuntDouble in specified snapshot
      * @param pos position to be set 
      * @param snapshotNo 
      * @see #getPos
@@ -224,24 +255,24 @@ namespace oopse{
     }
        
     /**
-     * Returns the previous velocity of this stuntdouble
-     * @return the velocity of this stuntdouble
+     * Returns the previous velocity of this stuntDouble
+     * @return the velocity of this stuntDouble
      */    
     Vector3d getPrevVel() {
       return ((snapshotMan_->getPrevSnapshot())->*storage_).velocity[localIndex_];
     }
        
     /**
-     * Returns the current velocity of this stuntdouble
-     * @return the velocity of this stuntdouble
+     * Returns the current velocity of this stuntDouble
+     * @return the velocity of this stuntDouble
      */    
     Vector3d getVel() {
       return ((snapshotMan_->getCurrentSnapshot())->*storage_).velocity[localIndex_];
     }
 
     /**
-     * Returns the velocity of this stuntdouble in specified snapshot 
-     * @return the velocity of this stuntdouble
+     * Returns the velocity of this stuntDouble in specified snapshot 
+     * @return the velocity of this stuntDouble
      * @param snapshotNo
      */    
     Vector3d getVel(int snapshotNo) {
@@ -249,7 +280,7 @@ namespace oopse{
     }
 
     /**
-     * Sets  the previous velocity of this stuntdouble
+     * Sets  the previous velocity of this stuntDouble
      * @param vel  new velocity 
      * @see #getVel
      */         
@@ -258,7 +289,7 @@ namespace oopse{
     }
        
     /**
-     * Sets  the current velocity of this stuntdouble
+     * Sets  the current velocity of this stuntDouble
      * @param vel  new velocity 
      */         
     void setVel(const Vector3d& vel) {
@@ -266,7 +297,7 @@ namespace oopse{
     }
 
     /**
-     * Sets  the velocity of this stuntdouble in specified snapshot
+     * Sets  the velocity of this stuntDouble in specified snapshot
      * @param vel velocity to be set 
      * @param snapshotNo 
      * @see #getVel
@@ -276,25 +307,25 @@ namespace oopse{
     }
 
     /**
-     * Returns the previous rotation matrix of this stuntdouble
-     * @return the rotation matrix of this stuntdouble
+     * Returns the previous rotation matrix of this stuntDouble
+     * @return the rotation matrix of this stuntDouble
      */    
     RotMat3x3d getPrevA() {
       return ((snapshotMan_->getPrevSnapshot())->*storage_).aMat[localIndex_];
     }
        
     /**
-     * Returns the current rotation matrix of this stuntdouble
-     * @return the rotation matrix of this stuntdouble
+     * Returns the current rotation matrix of this stuntDouble
+     * @return the rotation matrix of this stuntDouble
      */    
     RotMat3x3d getA() {
       return ((snapshotMan_->getCurrentSnapshot())->*storage_).aMat[localIndex_];
     }
 
     /**
-     * Returns the rotation matrix of this stuntdouble in specified snapshot 
+     * Returns the rotation matrix of this stuntDouble in specified snapshot 
      *
-     * @return the rotation matrix of this stuntdouble
+     * @return the rotation matrix of this stuntDouble
      * @param snapshotNo
      */    
     RotMat3x3d getA(int snapshotNo) {
@@ -302,7 +333,7 @@ namespace oopse{
     }
 
     /**
-     * Sets  the previous rotation matrix of this stuntdouble
+     * Sets  the previous rotation matrix of this stuntDouble
      * @param a  new rotation matrix 
      * @see #getA
      */         
@@ -311,7 +342,7 @@ namespace oopse{
     }
        
     /**
-     * Sets  the current rotation matrix of this stuntdouble
+     * Sets  the current rotation matrix of this stuntDouble
      * @param a  new rotation matrix 
      */         
     virtual void setA(const RotMat3x3d& a) {
@@ -319,7 +350,7 @@ namespace oopse{
     }
 
     /**
-     * Sets  the rotation matrix of this stuntdouble in specified snapshot
+     * Sets  the rotation matrix of this stuntDouble in specified snapshot
      * @param a rotation matrix to be set 
      * @param snapshotNo 
      * @see #getA
@@ -329,24 +360,24 @@ namespace oopse{
     }       
 
     /**
-     * Returns the previous angular momentum of this stuntdouble (body-fixed).
-     * @return the angular momentum of this stuntdouble
+     * Returns the previous angular momentum of this stuntDouble (body-fixed).
+     * @return the angular momentum of this stuntDouble
      */    
     Vector3d getPrevJ() {
       return ((snapshotMan_->getPrevSnapshot())->*storage_).angularMomentum[localIndex_];
     }
        
     /**
-     * Returns the current angular momentum of this stuntdouble (body -fixed).
-     * @return the angular momentum of this stuntdouble
+     * Returns the current angular momentum of this stuntDouble (body -fixed).
+     * @return the angular momentum of this stuntDouble
      */    
     Vector3d getJ() {
       return ((snapshotMan_->getCurrentSnapshot())->*storage_).angularMomentum[localIndex_];
     }
 
     /**
-     * Returns the angular momentum of this stuntdouble in specified snapshot (body-fixed).
-     * @return the angular momentum of this stuntdouble
+     * Returns the angular momentum of this stuntDouble in specified snapshot (body-fixed).
+     * @return the angular momentum of this stuntDouble
      * @param snapshotNo
      */    
     Vector3d getJ(int snapshotNo) {
@@ -354,7 +385,7 @@ namespace oopse{
     }
 
     /**
-     * Sets  the previous angular momentum of this stuntdouble (body-fixed).
+     * Sets  the previous angular momentum of this stuntDouble (body-fixed).
      * @param angMom  new angular momentum 
      * @see #getJ
      */         
@@ -363,7 +394,7 @@ namespace oopse{
     }
        
     /**
-     * Sets  the current angular momentum of this stuntdouble (body-fixed).
+     * Sets  the current angular momentum of this stuntDouble (body-fixed).
      * @param angMom  new angular momentum 
      */         
     void setJ(const Vector3d& angMom) {
@@ -371,7 +402,7 @@ namespace oopse{
     }
 
     /**
-     * Sets the angular momentum of this stuntdouble in specified snapshot(body-fixed).
+     * Sets the angular momentum of this stuntDouble in specified snapshot(body-fixed).
      * @param angMom angular momentum to be set 
      * @param snapshotNo 
      * @see #getJ
@@ -379,26 +410,76 @@ namespace oopse{
     void setJ(const Vector3d& angMom, int snapshotNo) {
       ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).angularMomentum[localIndex_] = angMom;
     }
-         
+    
     /**
-     * Returns the previous quaternion of this stuntdouble
-     * @return the quaternion of this stuntdouble
+     * Returns system Center of Mass for stuntDouble frame from snapshot
+     *
+     */ 
+    Vector3d getCOM(){
+      return (snapshotMan_->getCurrentSnapshot())->getCOM();
+    }
+    
+    /**
+      * Returns system Center of Mass velocity for stuntDouble frame from snapshot
+     *
+     */ 
+    
+    Vector3d getCOMvel(){
+      return (snapshotMan_->getCurrentSnapshot())->getCOMvel();
+    }
+    
+    /**
+      * Returns system Center of Mass angular momentum for stuntDouble frame from snapshot
+     *
+     */ 
+    Vector3d getCOMw(){
+      return (snapshotMan_->getCurrentSnapshot())->getCOMw();
+    }
+    
+/**
+     * Returns system Center of Mass for stuntDouble frame from snapshot
+     *
+     */ 
+    Vector3d getCOM(int snapshotNo){
+      return (snapshotMan_->getSnapshot(snapshotNo))->getCOM();
+    }
+    
+    /**
+      * Returns system Center of Mass velocity for stuntDouble frame from snapshot
+     *
+     */ 
+    
+    Vector3d getCOMvel(int snapshotNo){
+      return (snapshotMan_->getSnapshot(snapshotNo))->getCOMvel();
+    }
+    
+    /**
+      * Returns system Center of Mass angular momentum for stuntDouble frame from snapshot
+     *
+     */ 
+    Vector3d getCOMw(int snapshotNo){
+      return (snapshotMan_->getSnapshot(snapshotNo))->getCOMw();
+    }
+
+    /**
+     * Returns the previous quaternion of this stuntDouble
+     * @return the quaternion of this stuntDouble
      */    
     Quat4d getPrevQ() {
       return ((snapshotMan_->getPrevSnapshot())->*storage_).aMat[localIndex_].toQuaternion();
     }
        
     /**
-     * Returns the current quaternion of this stuntdouble
-     * @return the quaternion of this stuntdouble
+     * Returns the current quaternion of this stuntDouble
+     * @return the quaternion of this stuntDouble
      */    
     Quat4d getQ() {
       return ((snapshotMan_->getCurrentSnapshot())->*storage_).aMat[localIndex_].toQuaternion();
     }
 
     /**
-     * Returns the quaternion of this stuntdouble in specified snapshot 
-     * @return the quaternion of this stuntdouble
+     * Returns the quaternion of this stuntDouble in specified snapshot 
+     * @return the quaternion of this stuntDouble
      * @param snapshotNo
      */    
     Quat4d getQ(int snapshotNo) {
@@ -406,7 +487,7 @@ namespace oopse{
     }
 
     /**
-     * Sets  the previous quaternion of this stuntdouble
+     * Sets  the previous quaternion of this stuntDouble
      * @param q  new quaternion 
      * @note actual storage data is rotation matrix
      */         
@@ -415,7 +496,7 @@ namespace oopse{
     }
        
     /**
-     * Sets  the current quaternion of this stuntdouble
+     * Sets  the current quaternion of this stuntDouble
      * @param q  new quaternion 
      * @note actual storage data is rotation matrix
      */         
@@ -424,7 +505,7 @@ namespace oopse{
     }
 
     /**
-     * Sets  the quaternion of this stuntdouble in specified snapshot
+     * Sets  the quaternion of this stuntDouble in specified snapshot
      *
      * @param q quaternion to be set 
      * @param snapshotNo 
@@ -435,24 +516,24 @@ namespace oopse{
     }
 
     /**
-     * Returns the previous euler angles of this stuntdouble
-     * @return the euler angles of this stuntdouble
+     * Returns the previous euler angles of this stuntDouble
+     * @return the euler angles of this stuntDouble
      */    
     Vector3d getPrevEuler() {
       return ((snapshotMan_->getPrevSnapshot())->*storage_).aMat[localIndex_].toEulerAngles();
     }
        
     /**
-     * Returns the current euler angles of this stuntdouble
-     * @return the euler angles of this stuntdouble
+     * Returns the current euler angles of this stuntDouble
+     * @return the euler angles of this stuntDouble
      */    
     Vector3d getEuler() {
       return ((snapshotMan_->getCurrentSnapshot())->*storage_).aMat[localIndex_].toEulerAngles();
     }
 
     /**
-     * Returns the euler angles of this stuntdouble in specified snapshot.
-     * @return the euler angles of this stuntdouble
+     * Returns the euler angles of this stuntDouble in specified snapshot.
+     * @return the euler angles of this stuntDouble
      * @param snapshotNo
      */    
     Vector3d getEuler(int snapshotNo) {
@@ -460,7 +541,7 @@ namespace oopse{
     }
 
     /**
-     * Sets  the previous euler angles of this stuntdouble.
+     * Sets  the previous euler angles of this stuntDouble.
      * @param euler  new euler angles 
      * @see #getEuler
      * @note actual storage data is rotation matrix         
@@ -470,7 +551,7 @@ namespace oopse{
     }
        
     /**
-     * Sets  the current euler angles of this stuntdouble
+     * Sets  the current euler angles of this stuntDouble
      * @param euler  new euler angles 
      */         
     void setEuler(const Vector3d& euler) {
@@ -478,7 +559,7 @@ namespace oopse{
     }
 
     /**
-     * Sets  the euler angles  of this stuntdouble in specified snapshot
+     * Sets  the euler angles  of this stuntDouble in specified snapshot
      *
      * @param euler euler angles to be set 
      * @param snapshotNo 
@@ -487,53 +568,80 @@ namespace oopse{
     void setEuler(const Vector3d& euler, int snapshotNo) {
       ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).aMat[localIndex_] = euler;
     }
-       
+
     /**
-     * Returns the previous unit vectors of this stuntdouble
-     * @return the unit vectors of this stuntdouble
+     * Returns the previous dipole vector of this stuntDouble
+     * @return the dipole vector of this stuntDouble
      */    
-    RotMat3x3d getPrevElectroFrame() {
-      return ((snapshotMan_->getPrevSnapshot())->*storage_).electroFrame[localIndex_];
+    Vector3d getPrevDipole() {
+      return ((snapshotMan_->getPrevSnapshot())->*storage_).dipole[localIndex_];
     }
-       
+    
     /**
-     * Returns the current unit vectors of this stuntdouble
-     * @return the unit vectors of this stuntdouble
+     * Returns the current dipole vector of this stuntDouble
+     * @return the dipole vector of this stuntDouble
      */    
-    RotMat3x3d getElectroFrame() {
-      return ((snapshotMan_->getCurrentSnapshot())->*storage_).electroFrame[localIndex_];
+    Vector3d getDipole() {
+      return ((snapshotMan_->getCurrentSnapshot())->*storage_).dipole[localIndex_];
     }
-
+    
     /**
-     * Returns the unit vectors of this stuntdouble in specified snapshot 
+     * Returns the dipole vector of this stuntDouble in specified snapshot 
      *
-     * @return the unit vectors of this stuntdouble
+     * @return the dipole vector of this stuntDouble
      * @param snapshotNo
      */    
-    RotMat3x3d getElectroFrame(int snapshotNo) {
-      return ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).electroFrame[localIndex_];
+    Vector3d getDipole(int snapshotNo) {
+      return ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).dipole[localIndex_];
     }
 
+
     /**
-     * Returns the previous force of this stuntdouble
-     * @return the force of this stuntdouble
+     * Returns the previous quadrupole tensor of this stuntDouble
+     * @return the quadrupole tensor of this stuntDouble
      */    
+    Mat3x3d getPrevQuadrupole() {
+      return ((snapshotMan_->getPrevSnapshot())->*storage_).quadrupole[localIndex_];
+    }
+    
+    /**
+     * Returns the current quadrupole tensor of this stuntDouble
+     * @return the quadrupole tensor of this stuntDouble
+     */    
+    Mat3x3d getQuadrupole() {
+      return ((snapshotMan_->getCurrentSnapshot())->*storage_).quadrupole[localIndex_];
+    }
+    
+    /**
+     * Returns the quadrupole tensor of this stuntDouble in specified snapshot 
+     *
+     * @return the quadrupole tensor of this stuntDouble
+     * @param snapshotNo
+     */    
+    Mat3x3d getQuadrupole(int snapshotNo) {
+      return ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).quadrupole[localIndex_];
+    }
+        
+    /**
+     * Returns the previous force of this stuntDouble
+     * @return the force of this stuntDouble
+     */    
     Vector3d getPrevFrc() {
       return ((snapshotMan_->getPrevSnapshot())->*storage_).force[localIndex_];
     }
        
     /**
-     * Returns the current force of this stuntdouble
-     * @return the force of this stuntdouble
+     * Returns the current force of this stuntDouble
+     * @return the force of this stuntDouble
      */    
     Vector3d getFrc() {
       return ((snapshotMan_->getCurrentSnapshot())->*storage_).force[localIndex_];
     }
 
     /**
-     * Returns the force of this stuntdouble in specified snapshot 
+     * Returns the force of this stuntDouble in specified snapshot 
      *
-     * @return the force of this stuntdouble
+     * @return the force of this stuntDouble
      * @param snapshotNo
      */    
     Vector3d getFrc(int snapshotNo) {
@@ -541,7 +649,7 @@ namespace oopse{
     }
 
     /**
-     * Sets  the previous force of this stuntdouble
+     * Sets  the previous force of this stuntDouble
      *
      * @param frc  new force 
      * @see #getFrc
@@ -551,7 +659,7 @@ namespace oopse{
     }
        
     /**
-     * Sets  the current force of this stuntdouble
+     * Sets  the current force of this stuntDouble
      * @param frc  new force 
      */         
     void setFrc(const Vector3d& frc) {
@@ -559,7 +667,7 @@ namespace oopse{
     }
 
     /**
-     * Sets  the force of this stuntdouble in specified snapshot
+     * Sets  the force of this stuntDouble in specified snapshot
      *
      * @param frc force to be set 
      * @param snapshotNo 
@@ -570,7 +678,7 @@ namespace oopse{
     }
 
     /**
-     * Adds force into the previous force of this stuntdouble
+     * Adds force into the previous force of this stuntDouble
      *
      * @param frc  new force 
      * @see #getFrc
@@ -580,7 +688,7 @@ namespace oopse{
     }
        
     /**
-     * Adds force into the current force of this stuntdouble
+     * Adds force into the current force of this stuntDouble
      * @param frc  new force 
      */         
     void addFrc(const Vector3d& frc) {
@@ -588,7 +696,7 @@ namespace oopse{
     }
 
     /**
-     * Adds force into the force of this stuntdouble in specified snapshot
+     * Adds force into the force of this stuntDouble in specified snapshot
      *
      * @param frc force to be set 
      * @param snapshotNo 
@@ -599,25 +707,25 @@ namespace oopse{
     }
 
     /**
-     * Returns the previous torque of this stuntdouble
-     * @return the torque of this stuntdouble
+     * Returns the previous torque of this stuntDouble
+     * @return the torque of this stuntDouble
      */    
     Vector3d getPrevTrq() {
       return ((snapshotMan_->getPrevSnapshot())->*storage_).torque[localIndex_];
     }
        
     /**
-     * Returns the current torque of this stuntdouble
-     * @return the torque of this stuntdouble
+     * Returns the current torque of this stuntDouble
+     * @return the torque of this stuntDouble
      */    
     Vector3d getTrq() {
       return ((snapshotMan_->getCurrentSnapshot())->*storage_).torque[localIndex_];
     }
 
     /**
-     * Returns the torque of this stuntdouble in specified snapshot 
+     * Returns the torque of this stuntDouble in specified snapshot 
      *
-     * @return the torque of this stuntdouble
+     * @return the torque of this stuntDouble
      * @param snapshotNo
      */    
     Vector3d getTrq(int snapshotNo) {
@@ -625,7 +733,7 @@ namespace oopse{
     }
 
     /**
-     * Sets  the previous torque of this stuntdouble
+     * Sets  the previous torque of this stuntDouble
      *
      * @param trq  new torque 
      * @see #getTrq
@@ -635,7 +743,7 @@ namespace oopse{
     }
        
     /**
-     * Sets  the current torque of this stuntdouble
+     * Sets  the current torque of this stuntDouble
      * @param trq  new torque 
      */         
     void setTrq(const Vector3d& trq) {
@@ -643,7 +751,7 @@ namespace oopse{
     }
 
     /**
-     * Sets  the torque of this stuntdouble in specified snapshot
+     * Sets  the torque of this stuntDouble in specified snapshot
      *
      * @param trq torque to be set 
      * @param snapshotNo 
@@ -654,7 +762,7 @@ namespace oopse{
     }
 
     /**
-     * Adds torque into the previous torque of this stuntdouble
+     * Adds torque into the previous torque of this stuntDouble
      *
      * @param trq  new torque 
      * @see #getTrq
@@ -664,7 +772,7 @@ namespace oopse{
     }
        
     /**
-     * Adds torque into the current torque of this stuntdouble
+     * Adds torque into the current torque of this stuntDouble
      * @param trq  new torque 
      */         
     void addTrq(const Vector3d& trq) {
@@ -672,7 +780,7 @@ namespace oopse{
     }
 
     /**
-     * Adds torque into the torque of this stuntdouble in specified snapshot
+     * Adds torque into the torque of this stuntDouble in specified snapshot
      *
      * @param trq torque to be add 
      * @param snapshotNo 
@@ -683,103 +791,432 @@ namespace oopse{
     }       
 
 
+
     /**
-     * Returns the previous z-angle of this stuntdouble
-     * @return the z-angle of this stuntdouble
+     * Returns the previous particlePot of this stuntDouble
+     * @return the particlePot of this stuntDouble
      */    
-    double getPrevZangle() {
-      return ((snapshotMan_->getPrevSnapshot())->*storage_).zAngle[localIndex_];
+    RealType getPrevParticlePot() {
+      return ((snapshotMan_->getPrevSnapshot())->*storage_).particlePot[localIndex_];
     }
        
     /**
-     * Returns the current z-angle of this stuntdouble
-     * @return the z-angle of this stuntdouble
+     * Returns the current particlePot of this stuntDouble
+     * @return the particlePot of this stuntDouble
      */    
-    double getZangle() {
-      return ((snapshotMan_->getCurrentSnapshot())->*storage_).zAngle[localIndex_];
+    RealType getParticlePot() {
+      return ((snapshotMan_->getCurrentSnapshot())->*storage_).particlePot[localIndex_];
     }
 
     /**
-     * Returns the z-angle of this stuntdouble in specified snapshot
-     * @return the z-angle of this stuntdouble
+     * Returns the particlePot of this stuntDouble in specified snapshot 
+     *
+     * @return the particlePot of this stuntDouble
      * @param snapshotNo
      */    
-    double getZangle(int snapshotNo) {
-      return ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).zAngle[localIndex_];
+    RealType getParticlePot(int snapshotNo) {
+      return ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).particlePot[localIndex_];
     }
 
     /**
-     * Sets  the previous z-angle of this stuntdouble
-     * @param angle  new z-angle 
-     * @see #getZangle
+     * Sets  the previous particlePot of this stuntDouble
+     *
+     * @param particlePot  new particlePot 
+     * @see #getParticlePot
      */         
-    void setPrevZangle(double angle) {
-      ((snapshotMan_->getPrevSnapshot())->*storage_).zAngle[localIndex_] = angle;
+    void setPrevParticlePot(const RealType& particlePot) {
+      ((snapshotMan_->getPrevSnapshot())->*storage_).particlePot[localIndex_] = particlePot;
     }
        
     /**
-     * Sets  the current z-angle of this stuntdouble
-     * @param angle  new z-angle 
+     * Sets  the current particlePot of this stuntDouble
+     * @param particlePot  new particlePot 
      */         
-    void setZangle(double angle) {
-      ((snapshotMan_->getCurrentSnapshot())->*storage_).zAngle[localIndex_] = angle;
+    void setParticlePot(const RealType& particlePot) {
+      ((snapshotMan_->getCurrentSnapshot())->*storage_).particlePot[localIndex_] = particlePot;
     }
 
     /**
-     * Sets  the z-angle of this stuntdouble in specified snapshot
-     * @param angle z-angle to be set 
+     * Sets  the particlePot of this stuntDouble in specified snapshot
+     *
+     * @param particlePot particlePot to be set 
      * @param snapshotNo 
-     * @see #getZangle
+     * @see #getParticlePot
      */         
-    void setZangle(double angle, int snapshotNo) {
-      ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).zAngle[localIndex_] = angle;
+    void setParticlePot(const RealType& particlePot, int snapshotNo) {
+      ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).particlePot[localIndex_] = particlePot;
     }
 
     /**
-     * Adds z-angle into the previous z-angle of this stuntdouble
-     * @param angle  new z-angle 
-     * @see #getZangle
+     * Adds particlePot into the previous particlePot of this stuntDouble
+     *
+     * @param particlePot  new particlePot 
+     * @see #getParticlePot
      */         
-    void addPrevZangle(double angle) {
-      ((snapshotMan_->getPrevSnapshot())->*storage_).zAngle[localIndex_] += angle;
+    void addPrevParticlePot(const RealType& particlePot) {
+      ((snapshotMan_->getPrevSnapshot())->*storage_).particlePot[localIndex_] += particlePot;
     }
        
     /**
-     * Adds z-angle into the current z-angle of this stuntdouble
-     * @param angle  new z-angle 
+     * Adds particlePot into the current particlePot of this stuntDouble
+     * @param particlePot  new particlePot 
      */         
-    void addZangle(double angle) {
-      ((snapshotMan_->getCurrentSnapshot())->*storage_).zAngle[localIndex_] += angle;
+    void addParticlePot(const RealType& particlePot) {
+      ((snapshotMan_->getCurrentSnapshot())->*storage_).particlePot[localIndex_] += particlePot;
     }
 
     /**
-     * Adds z-angle into the z-angle of this stuntdouble in specified snapshot
-     * @param angle z-angle to be add 
+     * Adds particlePot into the particlePot of this stuntDouble in specified snapshot
+     *
+     * @param particlePot particlePot to be add 
      * @param snapshotNo 
-     * @see #getZangle
+     * @see #getParticlePot
      */         
-    void addZangle(double angle, int snapshotNo) {
-      ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).zAngle[localIndex_] += angle;
+    void addParticlePot(const RealType& particlePot, int snapshotNo) {
+      ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).particlePot[localIndex_] += particlePot;
     }       
 
-    /** Set the force of this stuntdouble to zero */
+    /**
+     * Returns the previous fluctuating charge of this stuntDouble
+     * @return the fluctuating charge of this stuntDouble
+     */    
+    RealType getPrevFlucQPos() {
+      return ((snapshotMan_->getPrevSnapshot())->*storage_).flucQPos[localIndex_];
+    }
+       
+    /**
+     * Returns the current fluctuating charge of this stuntDouble
+     * @return the fluctuating charge of this stuntDouble
+     */    
+    RealType getFlucQPos() {
+      return ((snapshotMan_->getCurrentSnapshot())->*storage_).flucQPos[localIndex_];
+    }
+
+    /**
+     * Returns the fluctuating charge of this stuntDouble in specified snapshot
+     * @return the fluctuating charge of this stuntDouble
+     * @param snapshotNo
+     */    
+    RealType getFlucQPos(int snapshotNo) {
+      return ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).flucQPos[localIndex_];
+    }
+
+    /**
+     * Sets  the previous fluctuating charge of this stuntDouble
+     * @param charge  new fluctuating charge 
+     * @see #getFlucQPos
+     */         
+    void setPrevFlucQPos(RealType charge) {
+      ((snapshotMan_->getPrevSnapshot())->*storage_).flucQPos[localIndex_] = charge;
+    }
+       
+    /**
+     * Sets  the current fluctuating charge of this stuntDouble
+     * @param charge  new fluctuating charge 
+     */         
+    void setFlucQPos(RealType charge) {
+      ((snapshotMan_->getCurrentSnapshot())->*storage_).flucQPos[localIndex_] = charge;
+    }
+
+    /**
+     * Sets  the fluctuating charge of this stuntDouble in specified snapshot
+     * @param charge fluctuating charge to be set 
+     * @param snapshotNo 
+     * @see #getFlucQPos
+     */         
+    void setFlucQPos(RealType charge, int snapshotNo) {
+      ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).flucQPos[localIndex_] = charge;
+    }
+
+    /**
+     * Adds fluctuating charge into the previous fluctuating charge of this stuntDouble
+     * @param charge  new fluctuating charge 
+     * @see #getFlucQPos
+     */         
+    void addPrevFlucQPos(RealType charge) {
+      ((snapshotMan_->getPrevSnapshot())->*storage_).flucQPos[localIndex_] += charge;
+    }
+       
+    /**
+     * Adds fluctuating charge into the current fluctuating charge of this stuntDouble
+     * @param charge  new fluctuating charge 
+     */         
+    void addFlucQPos(RealType charge) {
+      ((snapshotMan_->getCurrentSnapshot())->*storage_).flucQPos[localIndex_] += charge;
+    }
+
+    /**
+     * Adds fluctuating charge into the fluctuating charge of this stuntDouble in specified snapshot
+     * @param charge fluctuating charge to be add 
+     * @param snapshotNo 
+     * @see #getFlucQPos
+     */         
+    void addFlucQPos(RealType charge, int snapshotNo) {
+      ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).flucQPos[localIndex_] += charge;
+    }       
+
+
+    /**
+     * Returns the previous charge velocity of this stuntDouble
+     * @return the charge velocity of this stuntDouble
+     */    
+    RealType getPrevFlucQVel() {
+      return ((snapshotMan_->getPrevSnapshot())->*storage_).flucQVel[localIndex_];
+    }
+       
+    /**
+     * Returns the current charge velocity of this stuntDouble
+     * @return the charge velocity of this stuntDouble
+     */    
+    RealType getFlucQVel() {
+      return ((snapshotMan_->getCurrentSnapshot())->*storage_).flucQVel[localIndex_];
+    }
+
+    /**
+     * Returns the charge velocity of this stuntDouble in specified snapshot
+     * @return the charge velocity of this stuntDouble
+     * @param snapshotNo
+     */    
+    RealType getFlucQVel(int snapshotNo) {
+      return ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).flucQVel[localIndex_];
+    }
+
+    /**
+     * Sets  the previous charge velocity of this stuntDouble
+     * @param cvel  new charge velocity 
+     * @see #getFlucQVel
+     */         
+    void setPrevFlucQVel(RealType cvel) {
+      ((snapshotMan_->getPrevSnapshot())->*storage_).flucQVel[localIndex_] = cvel;
+    }
+       
+    /**
+     * Sets  the current charge velocity of this stuntDouble
+     * @param cvel  new charge velocity 
+     */         
+    void setFlucQVel(RealType cvel) {
+      ((snapshotMan_->getCurrentSnapshot())->*storage_).flucQVel[localIndex_] = cvel;
+    }
+
+    /**
+     * Sets  the charge velocity of this stuntDouble in specified snapshot
+     * @param cvel charge velocity to be set 
+     * @param snapshotNo 
+     * @see #getFlucQVel
+     */         
+    void setFlucQVel(RealType cvel, int snapshotNo) {
+      ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).flucQVel[localIndex_] = cvel;
+    }
+
+    /**
+     * Adds charge velocity into the previous charge velocity of this stuntDouble
+     * @param cvel  new charge velocity 
+     * @see #getFlucQVel
+     */         
+    void addPrevFlucQVel(RealType cvel) {
+      ((snapshotMan_->getPrevSnapshot())->*storage_).flucQVel[localIndex_] += cvel;
+    }
+       
+    /**
+     * Adds charge velocity into the current charge velocity of this stuntDouble
+     * @param cvel  new charge velocity 
+     */         
+    void addFlucQVel(RealType cvel) {
+      ((snapshotMan_->getCurrentSnapshot())->*storage_).flucQVel[localIndex_] += cvel;
+    }
+
+    /**
+     * Adds charge velocity into the charge velocity of this stuntDouble in specified snapshot
+     * @param cvel charge velocity to be add 
+     * @param snapshotNo 
+     * @see #getFlucQVel
+     */         
+    void addFlucQVel(RealType cvel, int snapshotNo) {
+      ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).flucQVel[localIndex_] += cvel;
+    }       
+
+
+    /**
+     * Returns the previous charge force of this stuntDouble
+     * @return the charge force of this stuntDouble
+     */    
+    RealType getPrevFlucQFrc() {
+      return ((snapshotMan_->getPrevSnapshot())->*storage_).flucQFrc[localIndex_];
+    }
+       
+    /**
+     * Returns the current charge force of this stuntDouble
+     * @return the charge force of this stuntDouble
+     */    
+    RealType getFlucQFrc() {
+      return ((snapshotMan_->getCurrentSnapshot())->*storage_).flucQFrc[localIndex_];
+    }
+
+    /**
+     * Returns the charge force of this stuntDouble in specified snapshot
+     * @return the charge force of this stuntDouble
+     * @param snapshotNo
+     */    
+    RealType getFlucQFrc(int snapshotNo) {
+      return ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).flucQFrc[localIndex_];
+    }
+
+    /**
+     * Sets  the previous charge force of this stuntDouble
+     * @param cfrc  new charge force 
+     * @see #getFlucQFrc
+     */         
+    void setPrevFlucQFrc(RealType cfrc) {
+      ((snapshotMan_->getPrevSnapshot())->*storage_).flucQFrc[localIndex_] = cfrc;
+    }
+       
+    /**
+     * Sets  the current charge force of this stuntDouble
+     * @param cfrc  new charge force 
+     */         
+    void setFlucQFrc(RealType cfrc) {
+      ((snapshotMan_->getCurrentSnapshot())->*storage_).flucQFrc[localIndex_] = cfrc;
+    }
+
+    /**
+     * Sets  the charge force of this stuntDouble in specified snapshot
+     * @param cfrc charge force to be set 
+     * @param snapshotNo 
+     * @see #getFlucQFrc
+     */         
+    void setFlucQFrc(RealType cfrc, int snapshotNo) {
+      ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).flucQFrc[localIndex_] = cfrc;
+    }
+
+    /**
+     * Adds charge force into the previous charge force of this stuntDouble
+     * @param cfrc   charge force to be added 
+     * @see #getFlucQFrc
+     */         
+    void addPrevFlucQFrc(RealType cfrc) {
+      ((snapshotMan_->getPrevSnapshot())->*storage_).flucQFrc[localIndex_] += cfrc;
+    }
+       
+    /**
+     * Adds charge force into the current charge force of this stuntDouble
+     * @param cfrc   charge force to be added 
+     */         
+    void addFlucQFrc(RealType cfrc) {
+      ((snapshotMan_->getCurrentSnapshot())->*storage_).flucQFrc[localIndex_] += cfrc;
+    }
+
+    /**
+     * Adds charge force into the charge force of this stuntDouble in specified snapshot
+     * @param cfrc charge force to be added
+     * @param snapshotNo 
+     * @see #getFlucQFrc
+     */         
+    void addFlucQFrc(RealType cfrc, int snapshotNo) {
+      ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).flucQFrc[localIndex_] += cfrc;
+    }       
+
+
+    /**
+     * Returns the previous electric field of this stuntDouble
+     * @return the electric field of this stuntDouble
+     */    
+    Vector3d getPrevElectricField() {
+      return ((snapshotMan_->getPrevSnapshot())->*storage_).electricField[localIndex_];
+    }
+       
+    /**
+     * Returns the current electric field of this stuntDouble
+     * @return the electric field of this stuntDouble
+     */    
+    Vector3d getElectricField() {
+      return ((snapshotMan_->getCurrentSnapshot())->*storage_).electricField[localIndex_];
+    }
+
+    /**
+     * Returns the electric field of this stuntDouble in specified snapshot 
+     * @return the electric field of this stuntDouble
+     * @param snapshotNo
+     */    
+    Vector3d getElectricField(int snapshotNo) {
+      return ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).electricField[localIndex_];
+    }
+
+    /**
+     * Sets the previous electric field of this stuntDouble
+     * @param eField  new electric field 
+     * @see #getElectricField
+     */         
+    void setPrevElectricField(const Vector3d& eField) {
+      ((snapshotMan_->getPrevSnapshot())->*storage_).electricField[localIndex_] = eField;
+    }
+       
+    /**
+     * Sets the current electric field of this stuntDouble
+     * @param eField  new electric field 
+     */         
+    void setElectricField(const Vector3d& eField) {
+      ((snapshotMan_->getCurrentSnapshot())->*storage_).electricField[localIndex_] = eField;
+    }
+
+    /**
+     * Sets the electric field of this stuntDouble in specified snapshot
+     * @param eField electric field to be set 
+     * @param snapshotNo 
+     * @see #getElectricField
+     */         
+    void setElectricField(const Vector3d& eField, int snapshotNo) {
+      ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).electricField[localIndex_] = eField;
+    }
+
+    /**
+     * Adds electric field into the previous electric field of this
+     * stuntDouble
+     *
+     * @param eField new electric field 
+     * @see #getElectricField
+     */         
+    void addPrevEelectricField(const Vector3d& eField) {
+      ((snapshotMan_->getPrevSnapshot())->*storage_).electricField[localIndex_] += eField;
+    }
+       
+    /**
+     * Adds electric field into the current electric field of this stuntDouble
+     * @param eField  new electric field 
+     */         
+    void addElectricField(const Vector3d& eField) {
+      ((snapshotMan_->getCurrentSnapshot())->*storage_).electricField[localIndex_] += eField;
+    }
+
+    /**
+     * Adds electric field into the electric field of this stuntDouble in specified snapshot
+     *
+     * @param eField electric field to be added
+     * @param snapshotNo 
+     * @see #getElectricField
+     */         
+    void addElectricField(const Vector3d& eField, int snapshotNo) {
+      ((snapshotMan_->getSnapshot(snapshotNo))->*storage_).electricField[localIndex_] += eField;
+    }       
+
+
+    /** Set the force of this stuntDouble to zero */
     void zeroForcesAndTorques(); 
     /**
-     * Returns the inertia tensor of this stuntdouble
-     * @return the inertia tensor of this stuntdouble
+     * Returns the inertia tensor of this stuntDouble
+     * @return the inertia tensor of this stuntDouble
      */ 
     virtual Mat3x3d getI() = 0;
 
     /**
-     * Returns the gradient of this stuntdouble
-     * @return the gradient of this stuntdouble
+     * Returns the gradient of this stuntDouble
+     * @return the gradient of this stuntDouble
      */ 
-    virtual std::vector<double> getGrad() = 0;
+    virtual std::vector<RealType> getGrad() = 0;
 
     /**
-     * Tests the  if this stuntdouble is a  linear rigidbody
+     * Tests the  if this stuntDouble is a  linear rigidbody
      *
-     * @return true if this stuntdouble is a  linear rigidbody, otherwise return false
+     * @return true if this stuntDouble is a  linear rigidbody, otherwise return false
      * @note atom and directional atom will always return false
      * 
      * @see #linearAxis
@@ -799,8 +1236,8 @@ namespace oopse{
       return linearAxis_;
     }
 
-    /** Returns the mass of this stuntdouble */
-    double getMass() {
+    /** Returns the mass of this stuntDouble */
+    RealType getMass() {
       return mass_;
     }
 
@@ -808,14 +1245,14 @@ namespace oopse{
      * Sets the mass of this stuntdoulbe
      * @param mass the mass to be set
      */         
-    void setMass(double mass) {
+    void setMass(RealType mass) {
       mass_ = mass;
     }
 
-    /** Returns the name of this stuntdouble */
+    /** Returns the name of this stuntDouble */
     virtual std::string getType() = 0;
         
-    /** Sets the name of this stuntdouble*/
+    /** Sets the name of this stuntDouble*/
     virtual void setType(const std::string& name) {}
 
     /**
@@ -843,13 +1280,16 @@ namespace oopse{
     Vector3d body2Lab(const Vector3d& v, int frame){
       return getA(frame).transpose() * v;
     }
+
     /**
      * <p>
-     * The purpose of the Visitor Pattern is to encapsulate an operation that you want to perform on
-     * the elements of a data structure. In this way, you can change the operation being performed 
-     * on a structure without the need of changing the classes of the elements that you are operating
-     * on. Using a Visitor pattern allows you to decouple the classes for the data structure and the 
-     * algorithms used upon them
+     * The purpose of the Visitor Pattern is to encapsulate an
+     * operation that you want to perform on the elements of a data
+     * structure. In this way, you can change the operation being
+     * performed on a structure without the need of changing the
+     * classes of the elements that you are operating on. Using a
+     * Visitor pattern allows you to decouple the classes for the data
+     * structure and the algorithms used upon them
      * </p>
      * @param v visitor
      */      
@@ -909,15 +1349,16 @@ namespace oopse{
 
         
     int globalIndex_;
+    int globalIntegrableObjectIndex_;
     int localIndex_;
 
 
-    double mass_;
+    RealType mass_;
         
   private:
         
     PropertyMap properties_;
   };
 
-}//end namespace oopse
+}//end namespace OpenMD
 #endif //PRIMITIVES_STUNTDOUBLE_HPP