From 9665d6b23d83a7385e71684a65e9f839d9446608 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Franck=20P=C3=A9rignon?= <franck.perignon@imag.fr>
Date: Wed, 5 Feb 2014 09:58:42 +0000
Subject: [PATCH] Review comm management in topo.

---
 HySoP/hysop/fields/discrete.py                |  4 +-
 HySoP/hysop/mpi/topology.py                   | 17 ++++----
 .../hysop/operator/monitors/compute_forces.py | 43 +++++++++----------
 HySoP/hysop/operator/reprojection.py          |  2 +-
 4 files changed, 32 insertions(+), 34 deletions(-)

diff --git a/HySoP/hysop/fields/discrete.py b/HySoP/hysop/fields/discrete.py
index 04b9778c2..db13fff4c 100644
--- a/HySoP/hysop/fields/discrete.py
+++ b/HySoP/hysop/fields/discrete.py
@@ -11,7 +11,6 @@ import numpy as np
 from parmepy.constants import ORDER, PARMES_REAL
 import parmepy.tools.numpywrappers as npw
 import numpy.linalg as la
-from parmepy.mpi import main_comm
 
 
 class DiscreteField(object):
@@ -188,7 +187,7 @@ class DiscreteField(object):
         indNoGhost = self.topology.mesh.iCompute
         for d in range(self.nbComponents):
             result[d] += la.norm(self.data[d][indNoGhost]) ** 2
-        main_comm.Allreduce(result, gResult)
+        self.topology.comm.Allreduce(result, gResult)
         return gResult ** (1. / 2)
 
     def get_data_method(self):
@@ -235,4 +234,3 @@ class DiscreteField(object):
         """ set all components to zero"""
         for dim in xrange(self.nbComponents):
             self.data[dim][...] = 0.0
-
diff --git a/HySoP/hysop/mpi/topology.py b/HySoP/hysop/mpi/topology.py
index 752b5d2da..3386a2705 100644
--- a/HySoP/hysop/mpi/topology.py
+++ b/HySoP/hysop/mpi/topology.py
@@ -72,8 +72,8 @@ class Cartesian(object):
         self.dim = dim
         # Mind the sequential case : dim from withConstr may be equal
         # to 0 if shape[:] = 1
-        ## Communicator used to build the topology
-        self.comm = comm
+        ## (Source) Communicator used to build the topology
+        self._comm_origin = comm
         ## True if everything has been properly set for this topology
         self.isUpToDate = False
         ## Grid of mpi process layout
@@ -122,15 +122,16 @@ class Cartesian(object):
         if(periods is not None):
             assert (periods.size == self.dim)
             self.periods[:] = periods[:]
-        self.topo = self.comm.Create_cart(self.shape[self.cutdir],
-                                          periods=self.periods, reorder=True)
+        self.comm = self._comm_origin.Create_cart(self.shape[self.cutdir],
+                                                  periods=self.periods,
+                                                  reorder=True)
 
         ## Size of the topology (i.e. total number of mpi processes)
-        self.size = self.topo.Get_size()
+        self.size = self.comm.Get_size()
         ## Rank of the current process in the topology
-        self.rank = self.topo.Get_rank()
+        self.rank = self.comm.Get_rank()
         ## Coordinates of the current process
-        self.reduced_coords = np.asarray(self.topo.Get_coords(self.rank),
+        self.reduced_coords = np.asarray(self.comm.Get_coords(self.rank),
                                          dtype=PARMES_INDEX)
         ## Coordinates of the current process
         ## What is different between proc_coords and reduced_coords?
@@ -144,7 +145,7 @@ class Cartesian(object):
         ## (warning : direction in the grid of process).
         self.neighbours = np.zeros((2, self.dim), dtype=PARMES_INTEGER)
         for direction in range(self.dim):
-            self.neighbours[:, direction] = self.topo.Shift(direction, 1)
+            self.neighbours[:, direction] = self.comm.Shift(direction, 1)
 
         ## ghost layer (default = 0)
         if(ghosts is None):
diff --git a/HySoP/hysop/operator/monitors/compute_forces.py b/HySoP/hysop/operator/monitors/compute_forces.py
index 38099234b..021e9ab10 100644
--- a/HySoP/hysop/operator/monitors/compute_forces.py
+++ b/HySoP/hysop/operator/monitors/compute_forces.py
@@ -24,11 +24,13 @@ class DragAndLift(Monitoring):
     The present class implements formula (52) of Plouhmans2002.
     Integral inside the obstacle is not taken into account.
     """
-    def __init__(self, velocity, vorticity, nu, coefForce, topo, volumeOfControl,
-                 obstacles=None, frequency=1, filename=None, safeOutput=None):
+    def __init__(self, velocity, vorticity, nu, coefForce, topo,
+                 volumeOfControl, obstacles=None, frequency=1,
+                 filename=None, safeOutput=True):
         """
         @param velocity field
-        @param vorticity field@
+        @param vorticity field
+        @param nu viscosity
         @param the topology on which forces will be computed
         @param a volume of control
         (parmepy.domain.obstacle.controlBox.ControlBox object)
@@ -77,9 +79,9 @@ class DragAndLift(Monitoring):
         self.vd = None
         # discrete vorticity field
         self.wd = None
-        # viscosity
+        ## viscosity
         self.nu = nu
-        # Normalizing coefficient for forces 
+        # Normalizing coefficient for forces
         # (based on the physics of the flow)
         self.coefForce = coefForce
         # Output file
@@ -94,17 +96,14 @@ class DragAndLift(Monitoring):
                 d = os.path.dirname(self.filename)
                 if not os.path.exists(d):
                     os.makedirs(d)
-            if safeOutput is None:
-                ## Defines how often
-                ## output file is written :
-                ## True --> open/close file everytime
-                ## data are written.
-                ## False --> open at init and close
-                ## during finalize. Cost less but if simu
-                ## crashes, data are lost.
-                self.safeOutput = True
-            else:
-                self.safeOutput = safeOutput
+            ## Defines how often
+            ## output file is written :
+            ## True --> open/close file everytime
+            ## data are written.
+            ## False --> open at init and close
+            ## during finalize. Cost less but if simu
+            ## crashes, data are lost.
+            self.safeOutput = safeOutput
             if self.safeOutput:
                 self._writeFile = self._fullwrite
             else:
@@ -191,9 +190,9 @@ class DragAndLift(Monitoring):
 
         # Print results, if required
         if self._writeOuput and self.topo.rank == 0 and (ite % self.freq) == 0:
-            dataout = npw.zeros((1,4))
-            dataout[0,0] = simulation.time
-            dataout[0,1:] = self.force
+            dataout = npw.zeros((1, 4))
+            dataout[0, 0] = simulation.time
+            dataout[0, 1:] = self.force
             self._writeFile(dataout)
 
         return self.force
@@ -243,8 +242,8 @@ class DragAndLift(Monitoring):
         self._laplacian.fd_scheme.computeIndices(sl)
         for j in i2:
             self._work[...] = self._laplacian(vdata[j], self._work)
-            res[i1] += self._coeff * self.nu * normal[i1] * np.sum(coords[j]
-                                                         * self._work[sl])
+            res[i1] += self._coeff * self.nu * normal[i1]\
+                * np.sum(coords[j] * self._work[sl])
             res[j] -= self._coeff * self.nu * normal[i1] * coords[i1] * \
                 np.sum(self._work[sl])
         self._fd_scheme.computeIndices(sl)
@@ -572,7 +571,7 @@ class Forces(Monitoring):
 
             # Reduction operation
             if main_rank == 0:
-                comm = main_comm
+                comm = self.topo.comm
                 comm.Reduce([localForce, MPI.DOUBLE], [force, MPI.DOUBLE], 
                             op=MPI.SUM, root=0)
                 force = force * self.coef
diff --git a/HySoP/hysop/operator/reprojection.py b/HySoP/hysop/operator/reprojection.py
index 51f6085a3..c6d0f5355 100644
--- a/HySoP/hysop/operator/reprojection.py
+++ b/HySoP/hysop/operator/reprojection.py
@@ -11,7 +11,7 @@ from parmepy.numerics.finite_differences import FD_C_4, FD_C_2
 from parmepy.numerics.differential_operations import GradV
 import parmepy.tools.numpywrappers as npw
 from parmepy.numerics.updateGhosts import UpdateGhosts
-from parmepy.mpi import main_rank, main_comm, main_size, MPI
+from parmepy.mpi import main_rank, MPI
 from parmepy.tools.timers import timed_function
 
 
-- 
GitLab