--- trunk/OOPSE/libmdtools/mpiSimulation.cpp	2003/03/27 15:07:29	419
+++ trunk/OOPSE/libmdtools/mpiSimulation.cpp	2003/04/01 16:50:14	441
@@ -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++ ) {
-            atomIndex++;
-            AtomToProcMap[atomIndex] = which_proc;
-          }
-          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
@@ -172,31 +152,46 @@ int* mpiSimulation::divideLabor( void ){
           MolToProcMap[i] = which_proc;
           AtomsPerProc[which_proc] += add_atoms;
           for (j = 0 ; j < add_atoms; j++ ) {
-            atomIndex++;
-            AtomToProcMap[atomIndex] = which_proc;
+	    AtomToProcMap[atomIndex] = which_proc;
+	    atomIndex++;
+          }
+          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++ ) {
-            atomIndex++;
-            AtomToProcMap[atomIndex] = which_proc;
-          }
+	    AtomToProcMap[atomIndex] = which_proc;
+	    atomIndex++;
+           }
           done = 1;
           continue;
         } else {
@@ -234,6 +229,8 @@ int* mpiSimulation::divideLabor( void ){
     
     MPI::COMM_WORLD.Bcast(AtomsPerProc, mpiPlug->numberProcessors,
                           MPI_INT, 0);
+
+
   }
 
 
@@ -285,11 +282,11 @@ int* mpiSimulation::divideLabor( void ){
   local_index = 0;
   for (i = 0; i < mpiPlug->nAtomsGlobal; i++) {
     if (AtomToProcMap[i] == mpiPlug->myNode) {
-      local_index++;
       globalIndex[local_index] = i;
+      local_index++;
     }
   }
- 
+  
   return globalIndex;
 }
 
@@ -299,8 +296,11 @@ void mpiSimulation::mpiRefresh( void ){
   int isError, i;
   int *globalIndex = new int[mpiPlug->myNlocal];
 
-  for(i=0; i<mpiPlug->myNlocal; i++) globalIndex[i] = entryPlug->atoms[i]->getGlobalIndex();
+  // Fortran indexing needs to be increased by 1 in order to get the 2 languages to
+  // not barf
 
+  for(i=0; i<mpiPlug->myNlocal; i++) globalIndex[i] = entryPlug->atoms[i]->getGlobalIndex()+1;
+
   
   isError = 0;
   setFsimParallel( mpiPlug, &(entryPlug->n_atoms), globalIndex, &isError );