--- trunk/OOPSE/libmdtools/mpiSimulation.cpp	2003/03/27 19:21:42	422
+++ trunk/OOPSE/libmdtools/mpiSimulation.cpp	2003/03/28 19:30:59	434
@@ -1,5 +1,5 @@
 #ifdef IS_MPI
-
+#include <iostream>
 #include <cstdlib>
 #include <cstring>
 #include <cmath>
@@ -85,7 +85,7 @@ int* mpiSimulation::divideLabor( void ){
 
   myRandom = new randomSPRNG( baseSeed );
 
-  a = (double)mpiPlug->nMolGlobal / (double)mpiPlug->nAtomsGlobal;
+  a = 3.0 * (double)mpiPlug->nMolGlobal / (double)mpiPlug->nAtomsGlobal;
 
   // Initialize things that we'll send out later:
   for (i = 0; i < mpiPlug->numberProcessors; i++ ) {
@@ -133,28 +133,8 @@ int* mpiSimulation::divideLabor( void ){
         // How many atoms does this processor have?
         
         old_atoms = AtomsPerProc[which_proc];
-
-        // If the processor already had too many atoms, just skip this
-        // processor and try again.
-
-        if (old_atoms >= nTarget) continue;
-
         add_atoms = compStamps[MolComponentType[i]]->getNAtoms();
         new_atoms = old_atoms + add_atoms;
-    
-        // If we can add this molecule to this processor without sending
-        // it above nTarget, then go ahead and do it:
-    
-        if (new_atoms <= nTarget) {
-          MolToProcMap[i] = which_proc;
-          AtomsPerProc[which_proc] += add_atoms;
-          for (j = 0 ; j < add_atoms; j++ ) {
-	    AtomToProcMap[atomIndex] = which_proc;
-	    atomIndex++;
-          }
-          done = 1;
-          continue;
-        }
 
         // If we've been through this loop too many times, we need
         // to just give up and assign the molecule to this processor
@@ -178,19 +158,34 @@ int* mpiSimulation::divideLabor( void ){
           done = 1;
           continue;
         }
+    
+        // If we can add this molecule to this processor without sending
+        // it above nTarget, then go ahead and do it:
+    
+        if (new_atoms <= nTarget) {
+          MolToProcMap[i] = which_proc;
+          AtomsPerProc[which_proc] += add_atoms;
+          for (j = 0 ; j < add_atoms; j++ ) {
+	    AtomToProcMap[atomIndex] = which_proc;
+	    atomIndex++;
+          }
+          done = 1;
+          continue;
+        }
 
-        // The only situation left is where old_atoms < nTarget, but
-        // new_atoms > nTarget.   We want to accept this with some
-        // probability that dies off the farther we are from nTarget
 
+        // The only situation left is when new_atoms > nTarget.  We
+        // want to accept this with some probability that dies off the
+        // farther we are from nTarget
+
         // roughly:  x = new_atoms - nTarget
         //           Pacc(x) = exp(- a * x)
-        // where a = 1 / (average atoms per molecule)
+        // where a = penalty / (average atoms per molecule)
 
         x = (double) (new_atoms - nTarget);
         y = myRandom->getRandom();
-        
-        if (exp(- a * x) > y) {
+      
+        if (y < exp(- a * x)) {
           MolToProcMap[i] = which_proc;
           AtomsPerProc[which_proc] += add_atoms;
           for (j = 0 ; j < add_atoms; j++ ) {
@@ -234,6 +229,8 @@ int* mpiSimulation::divideLabor( void ){
     
     MPI::COMM_WORLD.Bcast(AtomsPerProc, mpiPlug->numberProcessors,
                           MPI_INT, 0);
+
+
   }
 
 
@@ -289,7 +286,7 @@ int* mpiSimulation::divideLabor( void ){
       local_index++;
     }
   }
- 
+  
   return globalIndex;
 }