From e3d2ecce5445499f1d3981ca33467383cbca1abf Mon Sep 17 00:00:00 2001
From: Erik Schnetter <schnetter@gmail.com>
Date: Mon, 10 Sep 2012 22:23:44 -0400
Subject: Carpet: Add "balanced" load balancing mechanism

Add new load balancing mechanism "balanced".
---
 Carpet/Carpet/interface.ccl    |  5 ++-
 Carpet/Carpet/param.ccl        |  6 +++
 Carpet/Carpet/src/Recompose.cc | 88 ++++++++++++++++++++++++++++++++++++++++++
 Carpet/Carpet/src/functions.hh |  4 ++
 Carpet/CarpetLib/interface.ccl |  1 +
 5 files changed, 102 insertions(+), 2 deletions(-)

diff --git a/Carpet/Carpet/interface.ccl b/Carpet/Carpet/interface.ccl
index 93953295a..c2fc042fe 100644
--- a/Carpet/Carpet/interface.ccl
+++ b/Carpet/Carpet/interface.ccl
@@ -15,7 +15,10 @@ uses include header: mpi_string.hh
 
 uses include header: defs.hh
 uses include header: dist.hh
+uses include header: typecase.hh
+uses include header: typeprops.hh
 
+uses include header: balance.hh
 uses include header: bbox.hh
 uses include header: bboxset.hh
 uses include header: fulltree.hh
@@ -30,8 +33,6 @@ uses include header: gf.hh
 uses include header: ggf.hh
 uses include header: gh.hh
 uses include header: th.hh
-uses include header: typecase.hh
-uses include header: typeprops.hh
 
 uses include header: operators.hh
 
diff --git a/Carpet/Carpet/param.ccl b/Carpet/Carpet/param.ccl
index 64cb3fb7c..801f151c2 100644
--- a/Carpet/Carpet/param.ccl
+++ b/Carpet/Carpet/param.ccl
@@ -443,6 +443,7 @@ KEYWORD processor_topology "How to determine the processor topology" STEERABLE=r
   "along-dir" :: "Split the region along one direction only"
   "automatic" :: "Choose the topology automatically"
   "recursive" :: "Choose the topology automatically, using a different algorithm that may lead to better load balancing"
+  "balanced"  :: "Choose the topology automatically, ensuring a maximum load imbalance"
 } "automatic"
 
 CCTK_INT processor_topology_3d_x "Number of processors in x-direction" STEERABLE=recover
@@ -514,6 +515,11 @@ CCTK_REAL ghost_zone_cost "Relative cost of ghost zones for 'recursive' load bal
   0:* :: ""
 } 0.025
 
+CCTK_REAL maximum_imbalance "Maximum load imbalance" STEERABLE=always
+{
+  (0.0:*) :: ""
+} 0.1
+
 BOOLEAN same_number_of_components_on_each_process "Ensure that each process has the same number of components, adding empty dummy components if necessary" STEERABLE=always
 {
 } "yes"
diff --git a/Carpet/Carpet/src/Recompose.cc b/Carpet/Carpet/src/Recompose.cc
index 3dc3517fe..030ae6f76 100644
--- a/Carpet/Carpet/src/Recompose.cc
+++ b/Carpet/Carpet/src/Recompose.cc
@@ -27,6 +27,7 @@
 
 #include <loopcontrol.h>
 
+#include <balance.hh>
 #include <bbox.hh>
 #include <bboxset.hh>
 #include <defs.hh>
@@ -1305,6 +1306,8 @@ namespace Carpet {
       SplitRegionsMaps_Automatic (cctkGH, superregss, regss);
     } else if (CCTK_EQUALS (processor_topology, "recursive")) {
       SplitRegionsMaps_Recursively (cctkGH, superregss, regss);
+    } else if (CCTK_EQUALS (processor_topology, "balanced")) {
+      SplitRegionsMaps_Balanced (cctkGH, superregss, regss);
     } else if (CCTK_EQUALS (processor_topology, "manual")) {
       assert (0);
 //       SplitRegionsMaps_AsSpecified (cctkGH, superregss, regss);
@@ -1845,6 +1848,91 @@ namespace Carpet {
   
   
   
+  void
+  SplitRegionsMaps_Balanced (cGH const *const cctkGH,
+                             vector<vector<region_t> >& superregss,
+                             vector<vector<region_t> >& regss)
+  {
+    DECLARE_CCTK_PARAMETERS;
+    
+     // Find map numbers
+    int const nmaps = superregss.size();
+    vector<int> map_index(maps, -1);
+    for (int m=0; m<nmaps; ++m) {
+      assert(not superregss.AT(m).empty());
+      if (not superregss.AT(m).empty()) {
+        assert(map_index.AT(superregss.AT(m).AT(0).map) == -1);
+        map_index.AT(superregss.AT(m).AT(0).map) = m;
+      }
+    }
+    
+    // Combine regions from all maps
+    vector<region_t> regs;
+    for (int m=0; m<nmaps; ++m) {
+      combine_regions(superregss.AT(m), regs);
+    }
+
+    // Create superregion tree structure
+    {
+      int const ncomps = regs.size();
+      for (int c=0; c<ncomps; ++c) {
+        assert(regs.AT(c).processor == -1);
+        assert(not regs.AT(c).processors);
+        pseudoregion_t const preg(regs.AT(c).extent, -1);
+        regs.AT(c).processors = new ipfulltree(preg);
+      }
+    }
+    
+    // Set up new superregions
+    // Note: The tree structure is now reachable from both the
+    // superregions and the regions
+    for (int m=0; m<nmaps; ++m) {
+      superregss.AT(m).clear();
+    }
+    {
+      int const ncomps = regs.size();
+      for (int c=0; c<ncomps; ++c) {
+        int const m = map_index.AT(regs.AT(c).map);
+        superregss.AT(m).push_back(regs.AT(c));
+      }
+    }
+    
+    // Split onto processes
+    vector<vector<region_t> > split_regss;
+    int const nprocs = CCTK_nProcs(cctkGH);
+    int const nworkers = nprocs;
+    balance(regs, split_regss, nworkers,
+            maximum_imbalance, same_number_of_components_on_each_process);
+    
+    // Assign process numbers, and make the tree structure
+    // inaccessible from the regions
+    for (int p=0; p<nprocs; ++p) {
+      int const ncomps = split_regss.AT(p).size();
+      for (int c=0; c<ncomps; ++c) {
+        // Set process number in tree structure
+        assert(split_regss.AT(p).AT(c).processors);
+        assert(split_regss.AT(p).AT(c).processors->is_leaf());
+        assert(split_regss.AT(p).AT(c).processors->payload().component == -1);
+        split_regss.AT(p).AT(c).processors->payload().component = p;
+        // Set process number in region
+        split_regss.AT(p).AT(c).processors = NULL;
+        assert(split_regss.AT(p).AT(c).processor == -1);
+        split_regss.AT(p).AT(c).processor = p;
+      }
+    }
+    
+    // Distribute by map, implicitly assigning component numbers
+    for (int p=0; p<nprocs; ++p) {
+      int const ncomps = split_regss.AT(p).size();
+      for (int c=0; c<ncomps; ++c) {
+        int const m = map_index.AT(split_regss.AT(p).AT(c).map);
+        regss.AT(m).push_back(split_regss.AT(p).AT(c));
+      }
+    }
+  }
+  
+  
+  
   //////////////////////////////////////////////////////////////////////////////
   
   
diff --git a/Carpet/Carpet/src/functions.hh b/Carpet/Carpet/src/functions.hh
index d513d154b..163381ce4 100644
--- a/Carpet/Carpet/src/functions.hh
+++ b/Carpet/Carpet/src/functions.hh
@@ -200,6 +200,10 @@ namespace Carpet {
   SplitRegionsMaps_Recursively (cGH const * cctkGH,
                                 vector<vector<region_t> > & superregss,
                                 vector<vector<region_t> > & regss);
+  void
+  SplitRegionsMaps_Balanced (cGH const * cctkGH,
+                             vector<vector<region_t> > & superregss,
+                             vector<vector<region_t> > & regss);
   
   void
   MakeMultigridBoxes (cGH const * cctkGH,
diff --git a/Carpet/CarpetLib/interface.ccl b/Carpet/CarpetLib/interface.ccl
index 66ebc9f5d..12cba170c 100644
--- a/Carpet/CarpetLib/interface.ccl
+++ b/Carpet/CarpetLib/interface.ccl
@@ -9,6 +9,7 @@ includes header: dist.hh in dist.hh
 includes header: typecase.hh in typecase.hh
 includes header: typeprops.hh in typeprops.hh
 
+includes header: balance.hh in balance.hh
 includes header: bbox.hh in bbox.hh
 includes header: bboxset.hh in bboxset.hh
 includes header: fulltree.hh in fulltree.hh
-- 
cgit v1.2.3