--- trunk/OOPSE/libmdtools/SimSetup.cpp	2004/05/11 20:33:41	1157
+++ trunk/OOPSE/libmdtools/SimSetup.cpp	2004/05/12 20:54:10	1174
@@ -185,7 +185,9 @@ void SimSetup::makeMolecules(void){
   RigidBodyStamp* currentRigidBody;
   CutoffGroupStamp* currentCutoffGroup;
   CutoffGroup* myCutoffGroup;
-  
+  int nCutoffGroups;// number of cutoff group of a molecule defined in mdl file
+  set<int> cutoffAtomSet; //atoms belong to  cutoffgroup defined at mdl file
+
   bond_pair* theBonds;
   bend_set* theBends;
   torsion_set* theTorsions;
@@ -218,22 +220,23 @@ void SimSetup::makeMolecules(void){
       molInfo.nBends = comp_stamps[stampID]->getNBends();
       molInfo.nTorsions = comp_stamps[stampID]->getNTorsions();
       molInfo.nRigidBodies = comp_stamps[stampID]->getNRigidBodies();
-      molInfo.nCutoffGroups = comp_stamps[stampID]->getNCutoffGroups();
+
+      nCutoffGroups = comp_stamps[stampID]->getNCutoffGroups();
       
       molInfo.myAtoms = &(info[k].atoms[atomOffset]);
 
       if (molInfo.nBonds > 0) 
-        molInfo.myBonds = new (Bond *) [molInfo.nBonds];
+        molInfo.myBonds = new Bond*[molInfo.nBonds];
       else 
         molInfo.myBonds = NULL;
 
       if (molInfo.nBends > 0) 
-        molInfo.myBends = new (Bend *) [molInfo.nBends];
+        molInfo.myBends = new Bend*[molInfo.nBends];
       else 
         molInfo.myBends = NULL;
 
       if (molInfo.nTorsions > 0) 
-        molInfo.myTorsions = new (Torsion *) [molInfo.nTorsions];
+        molInfo.myTorsions = new Torsion *[molInfo.nTorsions];
       else 
         molInfo.myTorsions = NULL;
 
@@ -484,9 +487,12 @@ void SimSetup::makeMolecules(void){
       }
       
 
-      //creat cutoff group for molecule
+      //create cutoff group for molecule
+
+      cutoffAtomSet.clear();
       molInfo.myCutoffGroups.clear();
-      for (j = 0; j < molInfo.nCutoffGroups; j++){
+      
+      for (j = 0; j < nCutoffGroups; j++){
 
         currentCutoffGroup = comp_stamps[stampID]->getCutoffGroup(j);
         nMembers = currentCutoffGroup->getNMembers(); 
@@ -500,15 +506,30 @@ void SimSetup::makeMolecules(void){
 
           // tempI is atom numbering on local processor
           tempI = molI + atomOffset;
-
+          
           myCutoffGroup->addAtom(info[k].atoms[tempI]);          
+
+          cutoffAtomSet.insert(tempI);
         }
 
         molInfo.myCutoffGroups.push_back(myCutoffGroup);
       }//end for (j = 0; j < molInfo.nCutoffGroups; j++)
-      
 
+      //creat a cutoff group for every atom  in current molecule which does not belong to cutoffgroup defined at mdl file
 
+      for(j = 0; j < molInfo.nAtoms; j++){
+
+        if(cutoffAtomSet.find(molInfo.myAtoms[j]->getIndex()) == cutoffAtomSet.end()){
+          myCutoffGroup = new CutoffGroup();
+          myCutoffGroup->addAtom(molInfo.myAtoms[j]);
+          molInfo.myCutoffGroups.push_back(myCutoffGroup);
+        }
+          
+      }
+
+              
+
+
       // After this is all set up, scan through the atoms to 
       // see if they can be added to the integrableObjects:
 
@@ -564,18 +585,6 @@ void SimSetup::makeMolecules(void){
   MPIcheckPoint();
 #endif // is_mpi
 
-  // clean up the forcefield
-
-  if (!globals->haveRcut()){
-
-    the_ff->calcRcut();
-
-  } else {
-    
-    the_ff->setRcut( globals->getRcut() );
-  }
-
-  the_ff->cleanMe();
 }
 
 void SimSetup::initFromBass(void){
@@ -1016,10 +1025,33 @@ void SimSetup::finalInfoCheck(void){
 #endif //is_mpi
 
     double theRcut, theRsw;
+
+    if (globals->haveRcut()) {
+      theRcut = globals->getRcut();
 
+      if (globals->haveRsw()) 
+        theRsw = globals->getRsw();
+      else 
+        theRsw = theRcut;
+      
+      info[i].setDefaultRcut(theRcut, theRsw);
+
+    } else {
+      
+      the_ff->calcRcut();
+      theRcut = info[i].getRcut();
+
+      if (globals->haveRsw()) 
+        theRsw = globals->getRsw();
+      else 
+        theRsw = theRcut;
+      
+      info[i].setDefaultRcut(theRcut, theRsw);
+    }
+
     if (globals->getUseRF()){
       info[i].useReactionField = 1;
-
+      
       if (!globals->haveRcut()){
         sprintf(painCave.errMsg,
                 "SimSetup Warning: No value was set for the cutoffRadius.\n"
@@ -1096,6 +1128,9 @@ void SimSetup::finalInfoCheck(void){
   strcpy(checkPointMsg, "post processing checks out");
   MPIcheckPoint();
 #endif // is_mpi
+
+  // clean up the forcefield
+  the_ff->cleanMe();
 }
   
 void SimSetup::initSystemCoords(void){