--- trunk/OOPSE/libmdtools/ZConstraint.cpp	2003/08/22 20:04:39	711
+++ trunk/OOPSE/libmdtools/ZConstraint.cpp	2003/08/26 20:37:30	726
@@ -32,11 +32,11 @@ template<typename T> ZConstraint<T>::ZConstraint(SimIn
   if(!data){
     sprintf( painCave.errMsg,
                "ZConstraint Warning: User does not set force substraction policy, "
-               "average force substraction policy is used\n");
+               "PolicyByMass is used\n");
     painCave.isFatal = 0;
     simError();      
 
-    forcePolicy = (ForceSubstractionPolicy*) new PolicyByNumber(this);
+    forcePolicy = (ForceSubstractionPolicy*) new PolicyByMass(this);
   }
   else{
     policy = dynamic_cast<StringData*>(data);
@@ -44,11 +44,11 @@ template<typename T> ZConstraint<T>::ZConstraint(SimIn
     if(!policy){
       sprintf( painCave.errMsg,
                  "ZConstraint Error: Convertion from GenericData to StringData failure, "
-                 "average force substraction policy is used\n");
+                 "PolicyByMass is used\n");
       painCave.isFatal = 0;
       simError();      
 
-      forcePolicy = (ForceSubstractionPolicy*) new PolicyByNumber(this);
+      forcePolicy = (ForceSubstractionPolicy*) new PolicyByMass(this);
     }
     else{
       if(policy->getData() == "BYNUMBER")
@@ -276,8 +276,8 @@ template<typename T> ZConstraint<T>::ZConstraint(SimIn
       massOfZConsMols.push_back(molecules[i].getTotalMass());  
 
       zPos.push_back((*parameters)[searchResult].zPos);
-      cout << "index: "<< (*parameters)[searchResult].zconsIndex 
-             <<"\tzPos = " << (*parameters)[searchResult].zPos << endl; 
+//       cout << "index: "<< (*parameters)[searchResult].zconsIndex 
+//              <<"\tzPos = " << (*parameters)[searchResult].zPos << endl; 
       kz.push_back((*parameters)[searchResult]. kRatio * zForceConst);
       
       molecules[i].getCOM(COM);
@@ -536,11 +536,11 @@ template<typename T> void ZConstraint<T>::calcForce(in
 #ifdef IS_MPI
     if(worldRank == 0){
 #endif
-      std::cerr << "\n"
-		<< "*******************************************\n"
-		<< " about to call zeroOutVel()\n"
-		<< "*******************************************\n"
-		<< "\n";
+//       std::cerr << "\n"
+// 		<< "*******************************************\n"
+// 		<< " about to call zeroOutVel()\n"
+// 		<< "*******************************************\n"
+// 		<< "\n";
 #ifdef IS_MPI
     }
 #endif
@@ -549,11 +549,11 @@ template<typename T> void ZConstraint<T>::calcForce(in
 #ifdef IS_MPI
     if(worldRank == 0){
 #endif
-      std::cerr << "\n"
-		<< "*******************************************\n"
-		<< " finished zeroOutVel()\n"
-		<< "*******************************************\n"
-		<< "\n";
+//       std::cerr << "\n"
+// 		<< "*******************************************\n"
+// 		<< " finished zeroOutVel()\n"
+// 		<< "*******************************************\n"
+// 		<< "\n";
 #ifdef IS_MPI
     }
 #endif
@@ -566,10 +566,10 @@ template<typename T> void ZConstraint<T>::calcForce(in
 #ifdef IS_MPI
   if(worldRank == 0){
 #endif
-    cout << "---------------------------------------------------------------------" <<endl;
-    cout << "current time: " << info->getTime() << endl;
-    cout << "center of mass at z: " << zsys << endl;     
-    cout << "before calcForce, the COMVel of system is " << zSysCOMVel <<endl;
+//     cout << "---------------------------------------------------------------------" <<endl;
+//     cout << "current time: " << info->getTime() << endl;
+//     cout << "center of mass at z: " << zsys << endl;     
+//     cout << "before calcForce, the COMVel of system is " << zSysCOMVel <<endl;
 
 #ifdef IS_MPI
   }
@@ -608,7 +608,7 @@ template<typename T> void ZConstraint<T>::calcForce(in
 #ifdef IS_MPI
   if(worldRank == 0){
 #endif
-    cout << "after calcForce, the COMVel of system is " << zSysCOMVel <<endl;
+ //    cout << "after calcForce, the COMVel of system is " << zSysCOMVel <<endl;
 #ifdef IS_MPI
   }
 #endif
@@ -708,7 +708,7 @@ template<typename T> void ZConstraint<T>::zeroOutVel()
 #ifdef IS_MPI
   if(worldRank == 0){
 #endif
-    cout << "before resetting the COMVel of sytem is " << zSysCOMVel << endl;  
+//     cout << "before resetting the COMVel of sytem is " << zSysCOMVel << endl;  
 #ifdef IS_MPI
   }
 #endif
@@ -782,7 +782,7 @@ template<typename T> void ZConstraint<T>::zeroOutVel()
 #ifdef IS_MPI
   if(worldRank == 0){
 #endif
-    cout << "after resetting the COMVel of moving molecules is " << zSysCOMVel << endl;  
+//     cout << "after resetting the COMVel of moving molecules is " << zSysCOMVel << endl;  
 #ifdef IS_MPI
   }
 #endif
@@ -825,10 +825,11 @@ template<typename T> void ZConstraint<T>::doZconstrain
       } 
       totalFZ_local += fz[i];
 
-//       cout << "Fixed Molecule\tindex: " << indexOfZConsMols[i] 
-//              <<"\tcurrent zpos: " << COM[whichDirection] 
-//              << "\tcurrent fz: " <<fz[i] << endl;
+      //cout << "Fixed Molecule\tindex: " << indexOfZConsMols[i] 
+      //      <<"\tcurrent zpos: " << COM[whichDirection] 
+      //       << "\tcurrent fz: " <<fz[i] << endl;
 
+
     }
     
   } 
@@ -854,8 +855,8 @@ template<typename T> void ZConstraint<T>::doZconstrain
       zconsAtoms = zconsMols[i]->getMyAtoms();  
     
       for(int j =0; j < nAtomOfCurZConsMol; j++) {
-      force[whichDirection] = -fz[i]/ nAtomOfCurZConsMol;
-        //force[whichDirection] = - forcePolicy->getZFOfFixedZMols(zconsMols[i], zconsAtoms[j], fz[i]);
+        //force[whichDirection] = -fz[i]/ nAtomOfCurZConsMol;
+        force[whichDirection] = - forcePolicy->getZFOfFixedZMols(zconsMols[i], zconsAtoms[j], fz[i]);
         zconsAtoms[j]->addFrc(force);
       }
     
@@ -892,8 +893,8 @@ template<typename T> void ZConstraint<T>::doZconstrain
      Atom** unconsAtoms = unconsMols[i]->getMyAtoms();
      
      for(int j = 0; j < unconsMols[i]->getNAtoms(); j++){          
-       force[whichDirection] = totalFZ / (totNumOfUnconsAtoms + nMovingZAtoms);
-       //force[whichDirection] = forcePolicy->getZFOfMovingMols(unconsAtoms[j],totalFZ);
+       //force[whichDirection] = totalFZ / (totNumOfUnconsAtoms + nMovingZAtoms);
+       force[whichDirection] = forcePolicy->getZFOfMovingMols(unconsAtoms[j],totalFZ);
        unconsAtoms[j]->addFrc(force);
      }
      
@@ -906,8 +907,8 @@ template<typename T> void ZConstraint<T>::doZconstrain
       Atom** movingZAtoms = zconsMols[i]->getMyAtoms();     
 
       for(int j = 0; j < zconsMols[i]->getNAtoms(); j++){
-        force[whichDirection] = totalFZ / (totNumOfUnconsAtoms + nMovingZAtoms);
-        //force[whichDirection] = forcePolicy->getZFOfMovingMols(movingZAtoms[j],totalFZ);
+        //force[whichDirection] = totalFZ / (totNumOfUnconsAtoms + nMovingZAtoms);
+        force[whichDirection] = forcePolicy->getZFOfMovingMols(movingZAtoms[j],totalFZ);
         movingZAtoms[j]->addFrc(force);
       }
     }
@@ -942,13 +943,13 @@ template<typename T> void ZConstraint<T>::doHarmonic()
 
     if (states[i] == zcsMoving){
       zconsMols[i]->getCOM(COM);
-      cout << "Moving Molecule\tindex: " << indexOfZConsMols[i] 
-	   << "\tcurrent zpos: " << COM[whichDirection] << endl;
+//       cout << "Moving Molecule\tindex: " << indexOfZConsMols[i] 
+// 	   << "\tcurrent zpos: " << COM[whichDirection] << endl;
+
+      diff = COM[whichDirection] -zPos[i];
     
-    diff = COM[whichDirection] -zPos[i];
-    
       harmonicU = 0.5 * kz[i] * diff * diff;  
-    info->lrPot += harmonicU;
+      info->lrPot += harmonicU;
 
       harmonicF =  - kz[i] * diff;
       totalFZ_local += harmonicF;
@@ -958,8 +959,8 @@ template<typename T> void ZConstraint<T>::doHarmonic()
       Atom** movingZAtoms = zconsMols[i]->getMyAtoms();     
 
        for(int j = 0; j < zconsMols[i]->getNAtoms(); j++){          
-        force[whichDirection] = harmonicF / zconsMols[i]->getNAtoms();
-         //force[whichDirection] = forcePolicy->getHFOfFixedZMols(zconsMols[i], movingZAtoms[j], harmonicF);
+        //force[whichDirection] = harmonicF / zconsMols[i]->getNAtoms();
+        force[whichDirection] = forcePolicy->getHFOfFixedZMols(zconsMols[i], movingZAtoms[j], harmonicF);
          movingZAtoms[j]->addFrc(force);
        }
     }
@@ -982,8 +983,8 @@ template<typename T> void ZConstraint<T>::doHarmonic()
      Atom** unconsAtoms = unconsMols[i]->getMyAtoms();
      
      for(int j = 0; j < unconsMols[i]->getNAtoms(); j++){          
-       force[whichDirection] = - totalFZ /totNumOfUnconsAtoms;
-       //force[whichDirection] = - forcePolicy->getHFOfUnconsMols(unconsAtoms[j], totalFZ);
+       //force[whichDirection] = - totalFZ /totNumOfUnconsAtoms;
+       force[whichDirection] = - forcePolicy->getHFOfUnconsMols(unconsAtoms[j], totalFZ);
        unconsAtoms[j]->addFrc(force);     
      }
   }   
@@ -1023,8 +1024,9 @@ template<typename T> bool ZConstraint<T>::checkZConsSt
 #else
   MPI_Allreduce(&changed_local, &changed, 1, MPI_INT,MPI_SUM, MPI_COMM_WORLD);
 #endif
-  
+
   return (changed > 0);
+
 }
 
 template<typename T> bool ZConstraint<T>::haveFixedZMols(){
@@ -1046,7 +1048,7 @@ template<typename T> bool ZConstraint<T>::haveFixedZMo
   MPI_Allreduce(&havingFixed_local, &havingFixed, 1, MPI_INT,MPI_SUM, MPI_COMM_WORLD);
 #endif
 
-  return havingFixed > 0 ? true : false;
+  return (havingFixed > 0);
 }
 
 
@@ -1072,7 +1074,7 @@ template<typename T> bool ZConstraint<T>::haveMovingZM
   MPI_Allreduce(&havingMoving_local, &havingMoving, 1, MPI_INT,MPI_SUM, MPI_COMM_WORLD);
 #endif
 
-  return havingMoving > 0 ? true : false;
+  return (havingMoving > 0);
   
 }
 
@@ -1207,11 +1209,11 @@ template<typename T> void ZConstraint<T>::PolicyByNumb
 #ifdef IS_MPI
   if(worldRank == 0){
 #endif
-    std::cerr << "\n"
-	      << "*******************************************\n"
-	      << " fiished Policy by numbr()\n"
-	      << "*******************************************\n"
-	      << "\n";
+ //    std::cerr << "\n"
+// 	      << "*******************************************\n"
+// 	      << " fiished Policy by numbr()\n"
+// 	      << "*******************************************\n"
+// 	      << "\n";
 #ifdef IS_MPI
   }
 #endif