testing

hysop/backend/device/codegen/base/test.py

@@ -60,3 +60,36 @@ def _test_mesh_info(name, typegen,dim,ghosts,resolution,**kargs):
             **kargs)
 
     return (np_mis, cg_mis)
+
+def make_slice_views(compute_grid_size, 
+        lghosts=None, rghosts=None, step=None):
+    compute_grid_size = np.asarray(compute_grid_size)
+    dim = compute_grid_size.ndim
+
+    if (lghosts is None):
+        lghosts = (0,)*dim
+    elif np.iscalar(lghosts):
+        lghosts = (lghosts,)*dim
+    lghosts = np.asarray(lghosts)
+    
+    if (rghosts is None):
+        rghosts = (0,)*dim
+    elif np.iscalar(rghosts):
+        rghosts = (rghosts,)*dim
+    rghosts = np.asarray(rghosts)
+    
+    if (step is None):
+        step = (1,)*dim
+    elif np.iscalar(step):
+        step = (step,)*dim
+    step = np.asarray(step)
+
+    view  = [slice(lg, gs+lg) for lg, gs in zip(lghosts, compute_grid_size)]
+
+    grid_size = compute_grid_size + lghosts + rghosts
+    grid_shape = grid_size[::-1]
+
+    return view, grid_size, grid_shape, lghosts, rghosts, step
+
+    
+

hysop/backend/device/codegen/kernels/tests/test_directional_remesh.py

@@ -5,7 +5,7 @@ from hysop.deps import np, copy, math
 from hysop.backend.device.opencl import cl
 from hysop.tools.types import check_instance
 from hysop.constants import BoundaryCondition, Precision
-from hysop.backend.device.codegen.base.test import _test_mesh_info , _test_typegen
+from hysop.backend.device.codegen.base.test import _test_mesh_info , _test_typegen, make_slice_views
 from hysop.numerics.odesolvers.runge_kutta import ExplicitRungeKutta
 from hysop.backend.device.codegen.kernels.directional_remesh import DirectionalRemeshKernel
 from hysop.numerics.remesh.remesh import RemeshKernel
@@ -28,58 +28,72 @@ class TestDirectionalRemesh(object):
         dtype = cl_env.precision
         eps = np.finfo(dtype).eps
         
+        # compute grid
         compute_grid_size = np.asarray([64,10,10])
         (A,compute_mesh_info) = _test_mesh_info('base_mesh_info',typegen,3,0,compute_grid_size)
         
         dx     = A['dx'][0][0]
         inv_dx = A['inv_dx'][0][0]
         
+        # velocity bounds
         umax = +10.0
         umin = -10.0
         uinf = max(abs(umax),abs(umin))
         
+        # compute timestep
+        dt = cfl * dx/uinf
+        assert(umin<umax)
+        assert cfl>0
+        assert dt>0
+        
+        # scalar bounds
         S0_min=-1.0
         S0_max=+1.0
         
         S1_min=-10.0
         S1_max=+10.0
 
-        dt = cfl * dx/uinf
-        max_advec = int(np.floor(uinf*dt*inv_dx))
-        assert(umin<umax)
-        assert cfl>0
-        assert dt>0
-        assert max_advec>=0
-        
+        # compute max offsets due to advection and remeshing
         assert self.MAX_MOMENTS % 2 == 0
-        MAX_P=(1+MAX_MOMENTS/2)
-        MAX_ADVEC = int(np.ceil(dt*uinf*inv_dx))
+        MAX_P=(1+self.MAX_MOMENTS/2) 
+        MAX_ADVEC = int(np.floor(dt*uinf*inv_dx))
+        assert MAX_P > 0
+        assert MAX_ADVEC >= 0
 
         min_ghosts_in  = MAX_ADVEC
         min_ghosts_out = MAX_ADVEC + MAX_P
         
-        S0_in_ghosts   = min_ghosts_in + 0
-        S1_in_ghosts   = min_ghosts_in + 1
+        # inject extra ghosts on scalars for testing 
+        S0_in_ghosts   = min_ghosts_in  + 0
+        S1_in_ghosts   = min_ghosts_in  + 1
         S0_out_ghosts  = min_ghosts_out + 0
         S1_out_ghosts  = min_ghosts_out + 2
-
-        S0_grid_ghosts = np.asarray([S0_ghosts_in]*3)
-        S0_grid_size   = compute_grid_size + 2*S0_grid_ghosts
+    
         
-        S1_grid_ghosts = np.asarray([S1_ghosts]*3)
-        S1_grid_size   = compute_grid_size + 2*S1_grid_ghosts
-
-        (B, S0_mesh_info)  = _test_mesh_info('S0_mesh_info', 
-                typegen,3,S0_grid_ghosts,S0_grid_size)
-        (C, S1_mesh_info)  = _test_mesh_info('S1_mesh_info', 
-                typegen,3,S1_grid_ghosts,S1_grid_size)
+        # build 3D view, grid_size, grid_shape, ghosts, mesh_info, ...
+        # size is XYZ, shape is ZYX
+        S0_in_view, S0_in_grid_size, S0_in_grid_shape, S0_in_grid_ghosts, = \
+                make_slice_views(compute_grid_size, S0_in_ghosts, S0_in_ghosts)
+        S0_out_view, S0_out_grid_size, S0_out_grid_shape, S0_out_grid_ghosts, = \
+                make_slice_views(compute_grid_size, S0_out_ghosts, S0_out_ghosts)
+        (B, S0_mesh_info)  = _test_mesh_info('S0_mesh_info_in', 
+                typegen,3,S0_in_grid_ghosts,S0_in_grid_size)
+        (C, S0_mesh_info)  = _test_mesh_info('S0_mesh_info_out', 
+                typegen,3,S0_in_grid_ghosts,S0_in_grid_size)
+        
+        S1_in_view, S1_in_grid_size, S1_in_grid_shape, S1_in_grid_ghosts, = \
+                make_slice_views(compute_grid_size, S1_in_ghosts, S1_in_ghosts)
+        S1_out_view, S1_out_grid_size, S1_out_grid_shape, S1_out_grid_ghosts, = \
+                make_slice_views(compute_grid_size, S1_out_ghosts, S1_out_ghosts)
+        (D, S1_mesh_info)  = _test_mesh_info('S1_mesh_info_in', 
+                typegen,3,S1_in_grid_ghosts,S1_in_grid_size)
+        (E, S1_mesh_info)  = _test_mesh_info('S1_mesh_info_out', 
+                typegen,3,S1_in_grid_ghosts,S1_in_grid_size)
         
-        S0_grid_shape      = S0_grid_size[::-1]
-        S1_grid_shape      = S1_grid_size[::-1]
-        compute_grid_shape = compute_grid_size[::-1]
-
         assert A['dx'][0] == B['dx'][0]
         assert A['dx'][0] == C['dx'][0]
+        assert A['dx'][0] == D['dx'][0]
+        assert A['dx'][0] == E['dx'][0]
 
         compute_grid_bytes = compute_grid_size.size * typegen.FLT_BYTES[typegen.fbtype]