From 67c68c5a5cf1a0a1b6ae827639c7bd8117b7830c Mon Sep 17 00:00:00 2001
From: Jean-Baptiste Keck <Jean-Baptiste.Keck@imag.fr>
Date: Tue, 9 Apr 2019 17:35:32 +0200
Subject: [PATCH] fixed ghosts bug in accumulation

---
 .../operator/directional/advection_dir.py     |  3 ++-
 .../operator/directional/advection_dir.py     | 14 ++++++-----
 hysop/fields/cartesian_discrete_field.py      | 23 +++++++++++--------
 hysop/fields/ghost_exchangers.py              |  2 --
 hysop/numerics/interpolation/polynomial.py    |  2 +-
 hysop/operator/base/advection_dir.py          |  2 --
 6 files changed, 25 insertions(+), 21 deletions(-)

diff --git a/hysop/backend/device/opencl/operator/directional/advection_dir.py b/hysop/backend/device/opencl/operator/directional/advection_dir.py
index c7e3ff9c4..51cf8b3b8 100644
--- a/hysop/backend/device/opencl/operator/directional/advection_dir.py
+++ b/hysop/backend/device/opencl/operator/directional/advection_dir.py
@@ -158,8 +158,9 @@ class OpenClDirectionalAdvection(DirectionalAdvectionBase, OpenClDirectionalOper
         dsoutputs = self.dadvected_fields_out
         kl = OpenClKernelListLauncher(name='accumulate_and_exchange_ghosts')
         for sout in dsoutputs.values():
+            ghosts = tuple(sout.ghosts[:-1])+(self.remesh_ghosts,)
             kl += sout.accumulate_ghosts(directions=sout.dim-1,
-                                         ghosts=remesh_ghosts,
+                                         ghosts=ghosts,
                                          build_launcher=True)
             kl += sout.exchange_ghosts(build_launcher=True)
         self.accumulate_and_exchange = kl
diff --git a/hysop/backend/host/python/operator/directional/advection_dir.py b/hysop/backend/host/python/operator/directional/advection_dir.py
index 1c211b53e..f746d9bb6 100644
--- a/hysop/backend/host/python/operator/directional/advection_dir.py
+++ b/hysop/backend/host/python/operator/directional/advection_dir.py
@@ -12,7 +12,7 @@ from hysop.operator.base.advection_dir import DirectionalAdvectionBase
 from hysop.methods import Interpolation
 from hysop.numerics.odesolvers.runge_kutta import ExplicitRungeKutta, Euler, RK2, RK3, RK4
 
-DEBUG = True
+DEBUG = False
 
 class PythonDirectionalAdvection(DirectionalAdvectionBase, HostDirectionalOperator):
     counter = 0
@@ -150,14 +150,17 @@ class PythonDirectionalAdvection(DirectionalAdvectionBase, HostDirectionalOperat
             dump(Sin, 'Sin before remesh')
             self._compute_remesh()
             print 'S (before accumulation)'
-            print Sout.collect_data()
+            print Sout[0].sbuffer[Sout[0].local_slices(ghosts=(0,self.remesh_ghosts))]
             dump(Sin, 'Sout (after remesh)')
             for sout in dsoutputs.values():
-                sout.accumulate_ghosts(directions=sout.dim-1, ghosts=self.remesh_ghosts)
+                print 'Accumulate {}'.format(sout.short_description())
+                ghosts = tuple(sout.ghosts[:-1])+(self.remesh_ghosts,)
+                sout.accumulate_ghosts(directions=sout.dim-1, ghosts=ghosts)
             print 'S (after accumulation, before ghost exchange)'
             print Sout.collect_data()
             dump(Sin, 'Sout (after accumulation)')
             for sout in dsoutputs.values():
+                print 'Exchange {}'.format(sout.short_description())
                 sout.exchange_ghosts()
             print 'S (after ghost exchange)'
             print Sout.collect_data()
@@ -166,11 +169,10 @@ class PythonDirectionalAdvection(DirectionalAdvectionBase, HostDirectionalOperat
             self._compute_advection(dt)
             self._compute_remesh()
             for sout in dsoutputs.values():
-                sout.accumulate_ghosts(directions=sout.dim-1, ghosts=self.remesh_ghosts)
+                ghosts = tuple(sout.ghosts[:-1])+(self.remesh_ghosts,)
+                sout.accumulate_ghosts(directions=sout.dim-1, ghosts=ghosts)
             for sout in dsoutputs.values():
                 sout.exchange_ghosts()
-        import sys
-        sys.exit(1)
         
         self.counter += 1
 
diff --git a/hysop/fields/cartesian_discrete_field.py b/hysop/fields/cartesian_discrete_field.py
index d019383d0..a7457025f 100644
--- a/hysop/fields/cartesian_discrete_field.py
+++ b/hysop/fields/cartesian_discrete_field.py
@@ -1175,6 +1175,9 @@ class CartesianDiscreteScalarFieldView(CartesianDiscreteScalarFieldViewContainer
         Defaults to full ghosts exchange, including diagonals (ie. overwrite operation).
         """
         assert ('data' not in kwds)
+        msg='Passing ghosts as an integer is not supported anymore, use a tuple of size dim instead.'
+        if isinstance(ghosts, (int,long)):
+            raise RuntimeError(msg)
 
         directions = to_tuple(first_not_None(directions, range(self.dim)), cast=int)
         components = to_tuple(first_not_None(components, range(self.nb_components)), cast=int)
@@ -1184,11 +1187,10 @@ class CartesianDiscreteScalarFieldView(CartesianDiscreteScalarFieldViewContainer
         ghost_mask      = first_not_None(ghost_mask, GhostMask.FULL)
         exchange_method = first_not_None(exchange_method, ExchangeMethod.ISEND_IRECV)
 
-        assert len(directions) <= self.dim
+        assert len(ghosts) == self.dim, msg
+        assert all(g<=mg for (g,mg) in zip(ghosts, self.ghosts))
         assert len(directions)==len(set(directions))
-        assert len(ghosts) in (1,self.dim), ghosts
-        if len(ghosts)==1:
-            ghosts = tuple(ghosts[0] if (i in directions) else 0 for i in xrange(self.dim))
+        assert 0 < len(directions) <= self.dim
 
         if any(ghosts[i]>0 for i in directions):
             topology_state = self.topology_state
@@ -1225,6 +1227,10 @@ class CartesianDiscreteScalarFieldView(CartesianDiscreteScalarFieldViewContainer
         """
         Build a ghost exchanger for cartesian discrete fields, possibly on different data.
         """
+        msg='Passing ghosts as an integer is not supported anymore, use a tuple of size dim instead.'
+        if isinstance(ghosts, (int,long)):
+            raise RuntimeError(msg)
+
         ghost_op        = first_not_None(ghost_op, GhostOperation.EXCHANGE)
         ghost_mask      = first_not_None(ghost_mask, GhostMask.FULL)
         exchange_method = first_not_None(exchange_method, ExchangeMethod.ISEND_IRECV)
@@ -1234,13 +1240,12 @@ class CartesianDiscreteScalarFieldView(CartesianDiscreteScalarFieldViewContainer
         data       = to_tuple(first_not_None(data, self.data))
         
         directions = to_tuple(first_not_None(directions, range(self.dim)))
-        assert len(set(directions))==len(directions)
+        ghosts = first_not_None(ghosts, self.ghosts)
         
-        ghosts     = first_not_None(ghosts, self.ghosts)
-        if len(ghosts)==1:
-            ghosts = tuple(ghosts[0] if (i in directions) else 0 for i in xrange(self.dim))
-        assert len(ghosts)==self.dim
+        assert len(ghosts) == self.dim, msg
         assert all(g<=mg for (g,mg) in zip(ghosts, self.ghosts))
+        assert len(directions)==len(set(directions))
+        assert 0 < len(directions) <= self.dim
         
         if all(ghosts[i]==0 for i in directions):
             return None
diff --git a/hysop/fields/ghost_exchangers.py b/hysop/fields/ghost_exchangers.py
index b480090c8..5b3d1cd82 100644
--- a/hysop/fields/ghost_exchangers.py
+++ b/hysop/fields/ghost_exchangers.py
@@ -287,8 +287,6 @@ class CartesianDiscreteFieldGhostExchanger(GhostExchanger):
                                                               ghost_mask=ghost_mask)
         self.inner_ghosts = mesh.get_local_inner_ghost_slices(ghosts=ghosts,
                                                               ghost_mask=ghost_mask)
-        for slc in self.outer_ghosts:
-            print slc
         self.boundary_layers = mesh.get_boundary_layer_slices(ghosts=ghosts,
                                                               ghost_mask=ghost_mask)
         self.all_inner_ghost_slices = mesh.get_all_local_inner_ghost_slices(ghosts=ghosts)
diff --git a/hysop/numerics/interpolation/polynomial.py b/hysop/numerics/interpolation/polynomial.py
index 770cdde96..a57a2973a 100644
--- a/hysop/numerics/interpolation/polynomial.py
+++ b/hysop/numerics/interpolation/polynomial.py
@@ -8,7 +8,7 @@ try:
     has_flint = True
 except ImportError:
     import warnings
-    from hysop.tools.warnings import HysopPerformanceWarning
+    from hysop.tools.warning import HysopPerformanceWarning
     msg='Failed to import python-flint module, falling back to slow sympy solver.'
     warnings.warn(msg, HysopPerformanceWarning)
 
diff --git a/hysop/operator/base/advection_dir.py b/hysop/operator/base/advection_dir.py
index 51a134a88..0b21c7a4c 100644
--- a/hysop/operator/base/advection_dir.py
+++ b/hysop/operator/base/advection_dir.py
@@ -226,7 +226,6 @@ class DirectionalAdvectionBase(object):
                 _s_dx = _s_topo.mesh.space_step
             assert (_s_dx == s_dx).all()
             _s_requirements.min_ghosts = npw.maximum(_s_requirements.min_ghosts, min_scalar_ghosts).astype(HYSOP_INTEGER)
-            #_s_requirements.max_ghosts = npw.minimum(_s_requirements.max_ghosts, min_scalar_ghosts).astype(HYSOP_INTEGER)
 
         if is_inplace:
             for sfield in self.advected_fields_in:
@@ -237,7 +236,6 @@ class DirectionalAdvectionBase(object):
                     _s_dx = _s_topo.mesh.space_step
                 assert (_s_dx == s_dx).all()
                 _s_requirements.min_ghosts = npw.maximum(_s_requirements.min_ghosts, min_scalar_ghosts).astype(HYSOP_INTEGER)
-                #_s_requirements.max_ghosts = npw.minimum(_s_requirements.max_ghosts, min_scalar_ghosts).astype(HYSOP_INTEGER)
 
         self.scalar_cfl = scalar_cfl
         self.min_velocity_ghosts = min_velocity_ghosts
-- 
GitLab