--- branches/development/src/brains/SimInfo.cpp	2011/05/26 13:55:04	1569
+++ branches/development/src/brains/SimInfo.cpp	2011/06/16 22:00:08	1583
@@ -125,13 +125,8 @@ namespace OpenMD {
     //equal to the total number of atoms minus number of atoms belong to 
     //cutoff group defined in meta-data file plus the number of cutoff 
     //groups defined in meta-data file
-    std::cerr << "nGA = " << nGlobalAtoms_ << "\n";
-    std::cerr << "nCA = " << nCutoffAtoms << "\n";
-    std::cerr << "nG = " << nGroups << "\n";
 
     nGlobalCutoffGroups_ = nGlobalAtoms_ - nCutoffAtoms + nGroups;
-
-    std::cerr << "nGCG = " << nGlobalCutoffGroups_ << "\n";
     
     //every free atom (atom does not belong to rigid bodies) is an 
     //integrable object therefore the total number of integrable objects 
@@ -274,7 +269,26 @@ namespace OpenMD {
 #endif
     return fdf_;
   }
+  
+  unsigned int SimInfo::getNLocalCutoffGroups(){
+    int nLocalCutoffAtoms = 0;
+    Molecule* mol;
+    MoleculeIterator mi;
+    CutoffGroup* cg;
+    Molecule::CutoffGroupIterator ci;
     
+    for (mol = beginMolecule(mi); mol != NULL; mol  = nextMolecule(mi)) {
+      
+      for (cg = mol->beginCutoffGroup(ci); cg != NULL; 
+           cg = mol->nextCutoffGroup(ci)) {
+        nLocalCutoffAtoms += cg->getNumAtom();
+        
+      }        
+    }
+    
+    return nAtoms_ - nLocalCutoffAtoms + nCutoffGroups_;
+  }
+    
   void SimInfo::calcNdfRaw() {
     int ndfRaw_local;
 
@@ -745,7 +759,7 @@ namespace OpenMD {
       if ( simParams_->getAccumulateBoxDipole() ) {
 	calcBoxDipole_ = true;       
       }
-
+    
     set<AtomType*>::iterator i;
     set<AtomType*> atomTypes;
     atomTypes = getSimulatedAtomTypes();    
@@ -758,18 +772,28 @@ namespace OpenMD {
       usesMetallic |= (*i)->isMetal();
       usesDirectional |= (*i)->isDirectional();
     }
-
+    
 #ifdef IS_MPI    
     int temp;
     temp = usesDirectional;
     MPI_Allreduce(&temp, &usesDirectionalAtoms_, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
-
+    
     temp = usesMetallic;
     MPI_Allreduce(&temp, &usesMetallicAtoms_, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD);    
-
+    
     temp = usesElectrostatic;
     MPI_Allreduce(&temp, &usesElectrostaticAtoms_, 1, MPI_INT, MPI_LOR, MPI_COMM_WORLD); 
+#else
+
+    usesDirectionalAtoms_ = usesDirectional;
+    usesMetallicAtoms_ = usesMetallic;
+    usesElectrostaticAtoms_ = usesElectrostatic;
+
 #endif
+    
+    requiresPrepair_ = usesMetallicAtoms_ ? true : false; 
+    requiresSkipCorrection_ = usesElectrostaticAtoms_ ? true : false;
+    requiresSelfCorrection_ = usesElectrostaticAtoms_ ? true : false;    
   }
 
 
@@ -824,9 +848,16 @@ namespace OpenMD {
     Atom* atom;
     RealType totalMass;
 
-    //to avoid memory reallocation, reserve enough space for massFactors_
+    /**
+     * The mass factor is the relative mass of an atom to the total
+     * mass of the cutoff group it belongs to.  By default, all atoms
+     * are their own cutoff groups, and therefore have mass factors of
+     * 1.  We need some special handling for massless atoms, which
+     * will be treated as carrying the entire mass of the cutoff
+     * group.
+     */
     massFactors_.clear();
-    massFactors_.reserve(getNCutoffGroups());
+    massFactors_.resize(getNAtoms(), 1.0);
     
     for(mol = beginMolecule(mi); mol != NULL; mol = nextMolecule(mi)) {        
       for (cg = mol->beginCutoffGroup(ci); cg != NULL; 
@@ -835,10 +866,10 @@ namespace OpenMD {
 	totalMass = cg->getMass();
 	for(atom = cg->beginAtom(ai); atom != NULL; atom = cg->nextAtom(ai)) {
 	  // Check for massless groups - set mfact to 1 if true
-	  if (totalMass != 0)
-	    massFactors_.push_back(atom->getMass()/totalMass);
+	  if (totalMass != 0) 
+            massFactors_[atom->getLocalIndex()] = atom->getMass()/totalMass;
 	  else
-	    massFactors_.push_back( 1.0 );
+	    massFactors_[atom->getLocalIndex()] = 1.0;
 	}
       }       
     }
@@ -865,14 +896,6 @@ namespace OpenMD {
     int* oneThreeList = oneThreeInteractions_.getPairList();
     int* oneFourList = oneFourInteractions_.getPairList();
 
-    //setFortranSim( &fInfo_, &nGlobalAtoms_, &nAtoms_, &identArray_[0], 
-    //               &nExclude, excludeList, 
-    //               &nOneTwo, oneTwoList,
-    //               &nOneThree, oneThreeList,
-    //               &nOneFour, oneFourList,
-    //               &molMembershipArray[0], &mfact[0], &nCutoffGroups_, 
-    //               &fortranGlobalGroupMembership[0], &isError); 
-    
     topologyDone_ = true;
   }