Added test case

pystencils/boundaries/createindexlist.py

@@ -49,7 +49,7 @@ def _create_index_list_python(flag_field_arr, boundary_mask,
 
     # boundary cells are extracted via np.where. To ensure continous memory access in the compute kernel these cells
     # have to be sorted.
-    boundary_cells = np.transpose(np.where(flag_field_arr == boundary_mask))
+    boundary_cells = np.transpose(np.nonzero(flag_field_arr == boundary_mask))
     for i in range(len(flag_field_arr.shape)):
         boundary_cells = boundary_cells[boundary_cells[:, i].argsort(kind='mergesort')]
 

pystencils_tests/test_boundary_indexlist_creation.py

@@ -2,6 +2,8 @@ import numpy as np
 from itertools import product
 import pystencils.boundaries.createindexlist as cil
 
+from lbmpy.stencils import get_stencil
+
 import pytest
 
 @pytest.mark.parametrize('single_link', [False, True])
@@ -75,3 +77,42 @@ def test_equivalence_cell_idx_list_cython_python_version(single_link):
 
         np.testing.assert_equal(result_python_2d, result_cython_2d)
         np.testing.assert_equal(result_python_3d, result_cython_3d)
+
+@pytest.mark.parametrize('inner_or_boundary', [False, True])
+def test_normal_calculation(inner_or_boundary):
+    stencil = get_stencil("D2Q9")
+    domain_size = (32, 32)
+    dtype = np.uint32
+    fluid_mask = dtype(1)
+    mask = dtype(2)
+    flag_field = np.ones([domain_size[0], domain_size[1]], dtype=dtype) * fluid_mask
+
+    radius_inner = domain_size[0] // 4
+    radius_outer = domain_size[0] // 2
+    y_mid = domain_size[1] / 2
+    x_mid = domain_size[0] / 2
+
+    for x in range(0, domain_size[0]):
+        for y in range(0, domain_size[1]):
+            if (y - y_mid) ** 2 + (x - x_mid) ** 2 < radius_inner ** 2:
+                flag_field[x, y] = mask
+            if (x - x_mid) ** 2 + (y - y_mid) ** 2 > radius_outer ** 2:
+                flag_field[x, y] = mask
+
+    args_no_gl = (flag_field, mask, fluid_mask, np.array(stencil, dtype=np.int32), True)
+    index_list = cil._create_index_list_python(*args_no_gl, inner_or_boundary=inner_or_boundary, nr_of_ghost_layers=1)
+
+    checkmask = mask if inner_or_boundary else fluid_mask
+
+    for cell in index_list:
+        idx = cell[2]
+        cell = tuple((cell[0], cell[1]))
+        sum_cells = np.zeros(len(cell))
+        for dir_idx, direction in enumerate(stencil):
+            neighbor_cell = tuple([cell_i + dir_i for cell_i, dir_i in zip(cell, direction)])
+            if any(not 0 <= e < upper for e, upper in zip(neighbor_cell, flag_field.shape)):
+                continue
+            if flag_field[neighbor_cell] & checkmask:
+                sum_cells += np.array(direction)
+
+        assert np.argmax(np.inner(sum_cells, stencil)) == idx