From ed280f6f495169c2694231a3bb4e588f4951397a Mon Sep 17 00:00:00 2001
From: Jean-Baptiste Keck <Jean-Baptiste.Keck@imag.fr>
Date: Mon, 8 Apr 2019 23:22:57 +0200
Subject: [PATCH] bug in python advection
---
examples/multiresolution/scalar_advection.py | 2 +-
.../python/operator/directional/advection_dir.py | 4 +++-
hysop/fields/cartesian_discrete_field.py | 14 ++++++++++----
hysop/fields/ghost_exchangers.py | 2 ++
hysop/operator/base/advection_dir.py | 5 ++---
5 files changed, 18 insertions(+), 9 deletions(-)
diff --git a/examples/multiresolution/scalar_advection.py b/examples/multiresolution/scalar_advection.py
index 866a60b14..fa9d3a586 100644
--- a/examples/multiresolution/scalar_advection.py
+++ b/examples/multiresolution/scalar_advection.py
@@ -214,7 +214,7 @@ if __name__=='__main__':
parser = MultiResolutionScalarAdvectionArgParser()
parser.set_defaults(box_origin=(0.0,), box_length=(2*np.pi,),
- tstart=0.0, tend=2*np.pi, npts=(16,16,16,),
+ tstart=0.0, tend=2*np.pi, npts=(16,),
dump_freq=10, cfl=0.5, velocity=(1.0,),
ndim=3, compute_precision='fp64')
diff --git a/hysop/backend/host/python/operator/directional/advection_dir.py b/hysop/backend/host/python/operator/directional/advection_dir.py
index 2819b3a75..1c211b53e 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 = False
+DEBUG = True
class PythonDirectionalAdvection(DirectionalAdvectionBase, HostDirectionalOperator):
counter = 0
@@ -169,6 +169,8 @@ class PythonDirectionalAdvection(DirectionalAdvectionBase, HostDirectionalOperat
sout.accumulate_ghosts(directions=sout.dim-1, ghosts=self.remesh_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 fd64e3a9f..d019383d0 100644
--- a/hysop/fields/cartesian_discrete_field.py
+++ b/hysop/fields/cartesian_discrete_field.py
@@ -1189,12 +1189,12 @@ class CartesianDiscreteScalarFieldView(CartesianDiscreteScalarFieldViewContainer
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 all(g<=mg for (g,mg) in zip(ghosts, self.ghosts))
if any(ghosts[i]>0 for i in directions):
topology_state = self.topology_state
- key = (topology_state, ghosts, ghost_op, ghost_mask,
- exchange_method, components, directions)
+ key = (topology_state, ghosts,
+ ghost_op, ghost_mask,
+ exchange_method, components, directions)
if (key not in self._dfield._ghost_exchangers):
self._dfield._ghost_exchangers[key] = \
self.build_ghost_exchanger(ghosts=ghosts,
@@ -1220,7 +1220,8 @@ class CartesianDiscreteScalarFieldView(CartesianDiscreteScalarFieldViewContainer
return first_not_None(new_evt, evt)
def build_ghost_exchanger(self, name=None, components=None, directions=None,
- data=None, ghosts=None, ghost_op=None, ghost_mask=None, exchange_method=None):
+ data=None, ghosts=None,
+ ghost_op=None, ghost_mask=None, exchange_method=None):
"""
Build a ghost exchanger for cartesian discrete fields, possibly on different data.
"""
@@ -1231,15 +1232,20 @@ class CartesianDiscreteScalarFieldView(CartesianDiscreteScalarFieldViewContainer
check_instance(exchange_method, ExchangeMethod)
name = first_not_None(name, '{}_{}_{}'.format(self.name, ghosts, ghost_op))
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)
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 all(g<=mg for (g,mg) in zip(ghosts, self.ghosts))
if all(ghosts[i]==0 for i in directions):
return None
+ if not data:
+ return None
d0 = data[0]
if isinstance(d0, (np.ndarray, HostArray)):
diff --git a/hysop/fields/ghost_exchangers.py b/hysop/fields/ghost_exchangers.py
index 5b3d1cd82..b480090c8 100644
--- a/hysop/fields/ghost_exchangers.py
+++ b/hysop/fields/ghost_exchangers.py
@@ -287,6 +287,8 @@ 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/operator/base/advection_dir.py b/hysop/operator/base/advection_dir.py
index cda2d33ff..51a134a88 100644
--- a/hysop/operator/base/advection_dir.py
+++ b/hysop/operator/base/advection_dir.py
@@ -209,7 +209,6 @@ class DirectionalAdvectionBase(object):
min_velocity_ghosts[...] = 1
min_velocity_ghosts[ndir] = advection_ghosts
v_requirements.min_ghosts = npw.maximum(v_requirements.min_ghosts, min_velocity_ghosts).astype(HYSOP_INTEGER)
- #v_requirements.max_ghosts = npw.minimum(v_requirements.max_ghosts, min_velocity_ghosts).astype(HYSOP_INTEGER)
scalar_advection_ghosts = self.velocity_advection_ghosts(scalar_cfl)
remesh_ghosts = self.scalars_out_cache_ghosts(scalar_cfl, self.remesh_kernel)
@@ -227,7 +226,7 @@ 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)
+ #_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:
@@ -238,7 +237,7 @@ 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)
+ #_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