From ce8d55378b5d0df250a005c3c62aee550fa832ba Mon Sep 17 00:00:00 2001
From: Rafael Ravedutti <rafaelravedutti@gmail.com>
Date: Wed, 27 Mar 2024 10:58:30 +0100
Subject: [PATCH] First fixes for waLBerla communication

Signed-off-by: Rafael Ravedutti <rafaelravedutti@gmail.com>
---
 runtime/domain/block_forest.cpp      |  3 +++
 runtime/domain/block_forest.hpp      |  4 ----
 src/pairs/sim/comm.py                | 22 ++++++++++------------
 src/pairs/sim/domain_partitioning.py | 25 ++++++++++++++++++++++++-
 4 files changed, 37 insertions(+), 17 deletions(-)

diff --git a/runtime/domain/block_forest.cpp b/runtime/domain/block_forest.cpp
index 81a4166..1ce021e 100644
--- a/runtime/domain/block_forest.cpp
+++ b/runtime/domain/block_forest.cpp
@@ -1,3 +1,4 @@
+#include <map>
 #include <mpi.h>
 #include <vector>
 //---
@@ -21,6 +22,8 @@
 namespace pairs {
 
 void BlockForest::updateNeighborhood() {
+    std::map<int, std::vector<walberla::math::AABB>> neighborhood;
+    std::map<int, std::vector<walberla::BlockID>> blocks_pushed;
     auto me = walberla::mpi::MPIManager::instance()->rank();
     this->nranks = 0;
     this->total_aabbs = 0;
diff --git a/runtime/domain/block_forest.hpp b/runtime/domain/block_forest.hpp
index 673369d..acfa0c5 100644
--- a/runtime/domain/block_forest.hpp
+++ b/runtime/domain/block_forest.hpp
@@ -1,5 +1,3 @@
-#include <map>
-//---
 #include <blockforest/BlockForest.h>
 #include <blockforest/Initialization.h>
 #include <blockforest/loadbalancing/DynamicCurve.h>
@@ -27,8 +25,6 @@ private:
     PairsSimulation *ps;
     std::shared_ptr<walberla::BlockForest> forest;
     std::shared_ptr<walberla::blockforest::InfoCollection> info;
-    std::map<int, std::vector<walberla::math::AABB>> neighborhood;
-    std::map<int, std::vector<walberla::BlockID>> blocks_pushed;
     std::vector<int> ranks;
     std::vector<int> naabbs;
     std::vector<int> aabb_offsets;
diff --git a/src/pairs/sim/comm.py b/src/pairs/sim/comm.py
index 08fb143..fbd8de7 100644
--- a/src/pairs/sim/comm.py
+++ b/src/pairs/sim/comm.py
@@ -93,7 +93,7 @@ class Comm:
             CommunicateData(self, step, prop_list)
             UnpackGhostParticles(self, step, prop_list)
 
-            step_nrecv = sum([self.nrecv[j] for j in self.dom_part.step_indexes(step)])
+            step_nrecv = self.dom_part.reduce_sum_step_indexes(step, self.nrecv)
             Assign(self.sim, self.sim.nghost, self.sim.nghost + step_nrecv)
 
     @pairs_inline
@@ -306,9 +306,9 @@ class PackGhostParticles(Lowerable):
         send_mult = self.comm.send_mult
         self.sim.module_name(f"pack_ghost_particles{self.step}_" + "_".join([str(p.id()) for p in self.prop_list]))
 
-        step_indexes = self.comm.dom_part.step_indexes(self.step)
-        start = self.comm.send_offsets[step_indexes[0]]
-        for i in For(self.sim, start, ScalarOp.inline(start + sum([self.comm.nsend[j] for j in step_indexes]))):
+        start = self.comm.send_offsets[self.comm.dom_part.first_step_index(self.step)]
+        end = ScalarOp.inline(start + self.comm.dom_part.reduce_sum_step_indexes(self.step, self.comm.nsend))
+        for i in For(self.sim, start, end):
             p_offset = 0
             m = send_map[i]
             for p in self.prop_list:
@@ -347,9 +347,9 @@ class UnpackGhostParticles(Lowerable):
         recv_buffer.set_stride(1, self.get_elems_per_particle())
         self.sim.module_name(f"unpack_ghost_particles{self.step}_" + "_".join([str(p.id()) for p in self.prop_list]))
 
-        step_indexes = self.comm.dom_part.step_indexes(self.step)
-        start = self.comm.recv_offsets[step_indexes[0]]
-        for i in For(self.sim, start, ScalarOp.inline(start + sum([self.comm.nrecv[j] for j in step_indexes]))):
+        start = self.comm.recv_offsets[self.comm.dom_part.first_step_index(self.step)]
+        end = ScalarOp.inline(start + self.comm.dom_part.reduce_sum_step_indexes(self.step, self.comm.nrecv))
+        for i in For(self.sim, start, end):
             p_offset = 0
             for p in self.prop_list:
                 if not Types.is_scalar(p.type()):
@@ -422,9 +422,7 @@ class UnpackAllGhostParticles(Lowerable):
         recv_buffer.set_stride(1, self.get_elems_per_particle())
         self.sim.module_name(f"unpack_all_ghost_particles" + "_".join([str(p.id()) for p in self.prop_list]))
 
-        nrecv_size = sum([len(dom_part.step_indexes(s)) for s in range(dom_part.number_of_steps())])
-        nrecv_all = sum([self.comm.nrecv[j] for j in range(nrecv_size)])
-
+        nrecv_all = self.comm.dom_part.reduce_sum_all_steps(self.comm.nrecv)
         for i in For(self.sim, 0, nrecv_all):
             p_offset = 0
             for p in self.prop_list:
@@ -516,7 +514,8 @@ class ChangeSizeAfterExchange(Lowerable):
     def lower(self):
         self.sim.module_name(f"change_size_after_exchange{self.step}")
         self.sim.check_resize(self.sim.particle_capacity, self.sim.nlocal)
-        Assign(self.sim, self.sim.nlocal, self.sim.nlocal + sum([self.comm.nrecv[j] for j in self.comm.dom_part.step_indexes(self.step)]))
+        nrecv = self.comm.dom_part.reduce_sum_step_indexes(self.step, self.comm.nrecv)
+        Assign(self.sim, self.sim.nlocal, self.sim.nlocal + nrecv)
 
 
 class PackContactHistoryData(Lowerable):
@@ -611,7 +610,6 @@ class UnpackContactHistoryData(Lowerable):
         contact_used = self.sim._contact_history.contact_used
         self.sim.module_name(f"unpack_contact_history{self.step}")
 
-        step_indexes = self.comm.dom_part.step_indexes(self.step)
         nelems_per_contact = sum([Types.number_of_elements(self.sim, cp.type()) \
                                   for cp in self.sim.contact_properties]) + 1
 
diff --git a/src/pairs/sim/domain_partitioning.py b/src/pairs/sim/domain_partitioning.py
index 83875af..6cf3dbf 100644
--- a/src/pairs/sim/domain_partitioning.py
+++ b/src/pairs/sim/domain_partitioning.py
@@ -29,6 +29,16 @@ class DimensionRanges:
     def step_indexes(self, step):
         return [step * 2 + 0, step * 2 + 1]
 
+    def first_step_index(self, step):
+        return self.step_indexes(step)[0]
+
+    def reduce_sum_all_steps(self, array):
+        total_size = sum([len(self.step_indexes(s)) for s in range(self.number_of_steps())])
+        return sum([array[i] for i in range(total_size)])
+
+    def reduce_sum_step_indexes(self, step, array):
+       return sum([array[i] for i in self.step_indexes(step)])
+
     def initialize(self):
         grid_array = [(self.sim.grid.min(d), self.sim.grid.max(d)) for d in range(self.sim.ndims())]
         Call_Void(self.sim, "pairs->initDomain", [param for delim in grid_array for param in delim])
@@ -96,7 +106,20 @@ class BlockForest:
         return 1
 
     def step_indexes(self, step):
-        return [step]
+        yield from For(self.sim, 0, self.nranks)
+
+    def first_step_index(self, step):
+        return 0
+
+    def reduce_sum_all_steps(self, array):
+        return self.reduce_sum_step_indexes(0, array)
+
+    def reduce_sum_step_indexes(self, step, array):
+        nsend_sum = self.sim.add_temp_var(0)
+        for i in For(self.sim, 0, self.nranks):
+            Assign(self.sim, nsend_sum, nsend_sum + array[i])
+
+        return nsend_sum
 
     def initialize(self):
         grid_array = [(self.sim.grid.min(d), self.sim.grid.max(d)) for d in range(self.sim.ndims())]
-- 
GitLab