Carpet: Allow multi-model simulations

Split MPI_COMM_WORLD according to a new parameter Carpet::model into
different sets, and evolve these sets independently.

Add aliased functions to query the MPI communicators for this
simulation and for all models.

darcs-hash:20050825084739-891bb-7f102b46ec495f884be5197eaf077eef93f89dac.gz

9 files changed, 206 insertions, 3 deletions

diff --git a/Carpet/Carpet/interface.ccl b/Carpet/Carpet/interface.ccl
index 76d8f9b19..fa16cefcc 100644
--- a/Carpet/Carpet/interface.ccl
+++ b/Carpet/Carpet/interface.ccl
@@ -23,6 +23,17 @@ uses include header: operators.hh
 
 INCLUDE HEADER: typecase in carpet_typecase.hh
 
+
+
+# Get access to communicators
+CCTK_INT FUNCTION GetMPICommUniverse (CCTK_POINTER_TO_CONST IN cctkGH)
+PROVIDES FUNCTION GetMPICommUniverse WITH Carpet_GetMPICommUniverse LANGUAGE C
+
+CCTK_INT FUNCTION GetMPICommWorld (CCTK_POINTER_TO_CONST IN cctkGH)
+PROVIDES FUNCTION GetMPICommWorld WITH Carpet_GetMPICommWorld LANGUAGE C
+
+
+
 # Provide the function to en-/disable prolongation
 CCTK_INT FUNCTION EnableProlongating(CCTK_INT IN flag)
 PROVIDES FUNCTION EnableProlongating WITH CarpetEnableProlongating LANGUAGE C
diff --git a/Carpet/Carpet/param.ccl b/Carpet/Carpet/param.ccl
index 3bd189e38..35bb73f1e 100644
--- a/Carpet/Carpet/param.ccl
+++ b/Carpet/Carpet/param.ccl
@@ -147,6 +147,13 @@ CCTK_INT num_maps "Number of maps"
 
 
 
+STRING model "Model name for multi-model simulations -- the model name is used to distribute the processors onto the models"
+{
+  "" :: ""
+} "world"
+
+
+
 CCTK_INT prolongation_order_space "Order of prolongation operator in space" STEERABLE=recover
 {
   1 :: "first order (linear)"
diff --git a/Carpet/Carpet/src/MultiModel.cc b/Carpet/Carpet/src/MultiModel.cc
new file mode 100644
index 000000000..35bd17fb3
--- /dev/null
+++ b/Carpet/Carpet/src/MultiModel.cc
@@ -0,0 +1,133 @@
+#include <cassert>
+#include <cstring>
+#include <iostream>
+#include <map>
+#include <string>
+#include <vector>
+
+#include <mpi.h>
+
+#include "cctk.h"
+
+
+
+namespace Carpet
+{
+  
+  using namespace std;
+  
+  
+  
+  void
+  SplitWorld (MPI_Comm const world, string const model, MPI_Comm & comm,
+              bool const verbose)
+  {
+    // Get the total number of processors
+    int num_procs;
+    MPI_Comm_size (world, & num_procs);
+    int my_proc;
+    MPI_Comm_rank (world, & my_proc);
+    
+    // Gather the lengths of the model strings
+    int const length = model.length();
+    vector <int> lengths (num_procs);
+    
+    MPI_Allgather (const_cast <void *> (static_cast <void const *> (& length)),
+                   1, MPI_INT,
+                   static_cast <void *> (& lengths.front()), 1, MPI_INT,
+                   world);
+    
+    // Allocate space for all model strings
+    vector <int> offsets (num_procs + 1);
+    offsets.at(0) = 0;
+    for (int n = 0; n < num_procs; ++ n)
+    {
+      offsets.at(n + 1) = offsets.at(n) + lengths.at(n);
+    }
+    int const total_length = offsets.at(num_procs);
+    
+    // Gather all model strings
+    vector <char> models_buffer (total_length);
+    
+    MPI_Allgatherv (const_cast <void *>
+                    (static_cast <void const *> (model.c_str())),
+                    length, MPI_CHAR,
+                    static_cast <void *> (& models_buffer.front()),
+                    const_cast <int *>
+                    (static_cast <int const *> (& lengths.front())),
+                    const_cast <int *>
+                    (static_cast <int const *> (& offsets.front())),
+                    MPI_CHAR,
+                    world);
+    
+    // Convert model name buffer with C strings to C++ strings
+    vector <string> models (num_procs);
+    for (int n = 0; n < num_procs; ++ n)
+    {
+      models.at(n)
+        = string (& models_buffer.at (offsets.at(n)), lengths.at(n));
+    }
+    
+    // Map model strings to small integers
+    int num_models = 0;
+    vector <int> model_ids (num_procs);
+    std::map <string, int> model_map;
+    for (int n = 0; n < num_procs; ++ n)
+    {
+      if (model_map.find (models.at(n)) != model_map.end())
+      {
+        model_ids.at(n) = model_map[models.at(n)];
+      }
+      else
+      {
+        model_map[models.at(n)] = num_models;
+        model_ids.at(n) = num_models;
+        ++ num_models;
+      }
+    }
+    
+    // Determine processors per model
+    vector <int> num_model_procs (num_models, 0);
+    for (int n = 0; n < num_procs; ++ n)
+    {
+      ++ num_model_procs.at (model_ids.at(n));
+    }
+    
+    vector <vector <int> > model_procs (num_models);
+    for (int m = 0; m < num_models; ++ m)
+    {
+      model_procs.at(m).reserve (num_model_procs.at(m));
+    }
+    for (int n = 0; n < num_procs; ++ n)
+    {
+      model_procs.at (model_ids.at(n)).push_back (n);
+    }
+    for (int m = 0; m < num_models; ++ m)
+    {
+      assert (static_cast<int> (model_procs.at(m).size())
+              == num_model_procs.at(m));
+    }
+    
+    // Create a new communicator for each model
+    MPI_Comm_split (world, model_ids.at(my_proc), my_proc, & comm);
+    
+    if (verbose)
+    {
+      CCTK_INFO ("Multi-Model listing:");
+      for (int m = 0; m < num_models; ++ m)
+      {
+        cout << "   model " << m << ": \"" << models.at(m) << "\"" << endl;
+      }
+      CCTK_INFO ("Multi-Model processor distribution:");
+      for (int n = 0; n < num_procs; ++ n)
+      {
+        int const m = model_ids.at(n);
+        cout << "   processor " << n << ": "
+             << "model " << m << " \"" << models.at(m) << "\"" << endl;
+      }
+      CCTK_VInfo (CCTK_THORNSTRING,
+                  "Multi-Model: This is processor processor %d", my_proc);
+    }
+  }
+  
+} // namespace Carpet
diff --git a/Carpet/Carpet/src/SetupGH.cc b/Carpet/Carpet/src/SetupGH.cc
index 48bc71c58..39ca84187 100644
--- a/Carpet/Carpet/src/SetupGH.cc
+++ b/Carpet/Carpet/src/SetupGH.cc
@@ -356,6 +356,7 @@ namespace Carpet {
   static void enable_storage_for_all_groups (cGH* cgh);
   static void leave_all_modes (cGH* cgh);
 
+
   
   void* SetupGH (tFleshConfig* fc, int convLevel, cGH* cgh)
   {
@@ -366,13 +367,15 @@ namespace Carpet {
     // Not sure what to do with that
     assert (convLevel==0);
     
-    dist::pseudoinit();
-
+    comm_universe = MPI_COMM_WORLD;
+    SplitWorld (comm_universe, model, comm_world, verbose);
+    dist::pseudoinit (comm_world);
+    
     initialise_current_position ();
     
     Waypoint ("Setting up the grid hierarchy");
     
-    cgh->identity = strdup ("");
+    cgh->identity = strdup (model);
     
     // Processor information
     Output ("Carpet is running on %d processors", CCTK_nProcs(cgh));
diff --git a/Carpet/Carpet/src/functions.hh b/Carpet/Carpet/src/functions.hh
index 14ba8c50d..9bf8c8f04 100644
--- a/Carpet/Carpet/src/functions.hh
+++ b/Carpet/Carpet/src/functions.hh
@@ -3,6 +3,8 @@
 
 #include <vector>
 
+#include <mpi.h>
+
 #include "cctk.h"
 #include "cctk_Schedule.h"
 
@@ -39,6 +41,21 @@ namespace Carpet {
    
    
    
+  // Multi-Model
+  void
+  SplitWorld (MPI_Comm const world, string const model, MPI_Comm & comm,
+              bool verbose);
+  
+  extern "C" {
+    CCTK_INT
+    Carpet_GetMPICommUniverse (CCTK_POINTER_TO_CONST cctkGH);
+    CCTK_INT
+    Carpet_GetMPICommWorld (CCTK_POINTER_TO_CONST cctkGH);
+  }
+  
+  
+  
+  
   // Helpers for recomposing the grid hierarchy
   void CheckRegions (const gh::mexts & bbsss,
                      const gh::rbnds & obss,
diff --git a/Carpet/Carpet/src/helpers.cc b/Carpet/Carpet/src/helpers.cc
index af2fddd9e..a5140f7f9 100644
--- a/Carpet/Carpet/src/helpers.cc
+++ b/Carpet/Carpet/src/helpers.cc
@@ -48,6 +48,23 @@ namespace Carpet {
 
 
 
+  // Multi-Model
+  CCTK_INT
+  Carpet_GetMPICommUniverse (CCTK_POINTER_TO_CONST const cctkGH)
+  {
+    assert (comm_universe >= 0);
+    return comm_universe;
+  }
+
+  CCTK_INT
+  Carpet_GetMPICommWorld (CCTK_POINTER_TO_CONST const cctkGH)
+  {
+    assert (comm_world >= 0);
+    return comm_world;
+  }
+  
+
+
   // Communication
 
   int Barrier (const cGH* cgh)
diff --git a/Carpet/Carpet/src/make.code.defn b/Carpet/Carpet/src/make.code.defn
index a42787eac..08e694214 100644
--- a/Carpet/Carpet/src/make.code.defn
+++ b/Carpet/Carpet/src/make.code.defn
@@ -9,6 +9,7 @@ SRCS =	CallFunction.cc				\
 	Cycle.cc				\
 	Evolve.cc				\
 	Initialise.cc				\
+	MultiModel.cc				\
 	Poison.cc				\
 	Recompose.cc				\
 	Restrict.cc				\
diff --git a/Carpet/Carpet/src/variables.cc b/Carpet/Carpet/src/variables.cc
index 8b26a38e7..71520920c 100644
--- a/Carpet/Carpet/src/variables.cc
+++ b/Carpet/Carpet/src/variables.cc
@@ -102,4 +102,10 @@ namespace Carpet {
   // Data for everything
   vector<vector<arrdesc> > arrdata; // [group][map]
   
+  
+  
+  // MPI Communicators
+  MPI_Comm comm_universe = -1;
+  MPI_Comm comm_world = -1;
+  
 } // namespace Carpet
diff --git a/Carpet/Carpet/src/variables.hh b/Carpet/Carpet/src/variables.hh
index 1e5330d80..2d8a6d294 100644
--- a/Carpet/Carpet/src/variables.hh
+++ b/Carpet/Carpet/src/variables.hh
@@ -14,6 +14,8 @@
 
 #include <vector>
 
+#include <mpi.h>
+
 #include "cctk.h"
 
 #include "data.hh"
@@ -140,6 +142,12 @@ namespace Carpet {
   };
   extern vector<vector<arrdesc> > arrdata; // [group][map]
   
+  
+  
+  // MPI Communicators
+  extern MPI_Comm comm_universe;
+  extern MPI_Comm comm_world;
+  
 } // namespace Carpet
 
 #endif // !defined(VARIABLES_HH)