Skip to content
Snippets Groups Projects

Compare revisions

Changes are shown as if the source revision was being merged into the target revision. Learn more about comparing revisions.

Source

Select target project
No results found
Select Git revision
  • 66-absolute-access-is-probably-not-copied-correctly-after-_eval_subs
  • const_fix
  • fhennig/compiler-warnings
  • fhennig/v2.0-deprecations
  • fma
  • gpu_bufferfield_fix
  • gpu_liveness_opts
  • holzer-master-patch-46757
  • hyteg
  • improved_comm
  • master
  • target_dh_refactoring
  • v2.0-dev
  • vectorization_sqrt_fix
  • zikeliml/124-rework-tutorials
  • zikeliml/Task-96-dotExporterForAST
  • last/Kerncraft
  • last/LLVM
  • last/OpenCL
  • release/0.2.1
  • release/0.2.10
  • release/0.2.11
  • release/0.2.12
  • release/0.2.13
  • release/0.2.14
  • release/0.2.15
  • release/0.2.2
  • release/0.2.3
  • release/0.2.4
  • release/0.2.6
  • release/0.2.7
  • release/0.2.8
  • release/0.2.9
  • release/0.3.0
  • release/0.3.1
  • release/0.3.2
  • release/0.3.3
  • release/0.3.4
  • release/0.4.0
  • release/0.4.1
  • release/0.4.2
  • release/0.4.3
  • release/0.4.4
  • release/1.0
  • release/1.0.1
  • release/1.1
  • release/1.1.1
  • release/1.2
  • release/1.3
  • release/1.3.1
  • release/1.3.2
  • release/1.3.3
  • release/1.3.4
  • release/1.3.5
  • release/1.3.6
  • release/1.3.7
  • release/2.0.dev0
57 results

Target

Select target project
  • anirudh.jonnalagadda/pystencils
  • hyteg/pystencils
  • jbadwaik/pystencils
  • jngrad/pystencils
  • itischler/pystencils
  • ob28imeq/pystencils
  • hoenig/pystencils
  • Bindgen/pystencils
  • hammer/pystencils
  • da15siwa/pystencils
  • holzer/pystencils
  • alexander.reinauer/pystencils
  • ec93ujoh/pystencils
  • Harke/pystencils
  • seitz/pystencils
  • pycodegen/pystencils
16 results
Select Git revision
  • armneon
  • compare_fix
  • const_fix
  • gpu_liveness_opts
  • hyteg
  • improved_comm
  • jan_fix
  • jan_test
  • master
  • mr_parallel_datahandling_fix
  • philox-simd
  • target_dh_refactoring
  • test_martin
  • test_martin2
  • vectorization_sqrt_fix
  • release/0.2.1
  • release/0.2.10
  • release/0.2.11
  • release/0.2.12
  • release/0.2.13
  • release/0.2.14
  • release/0.2.15
  • release/0.2.2
  • release/0.2.3
  • release/0.2.4
  • release/0.2.6
  • release/0.2.7
  • release/0.2.8
  • release/0.2.9
29 results
Show changes
Showing
with 1477 additions and 49 deletions
...@@ -18,14 +18,25 @@ def check_equivalence(assignments, src_arr): ...@@ -18,14 +18,25 @@ def check_equivalence(assignments, src_arr):
for vectorization in [False, {'assume_inner_stride_one': True}]: for vectorization in [False, {'assume_inner_stride_one': True}]:
with_blocking = ps.create_kernel(assignments, cpu_blocking=(8, 16, 4), cpu_openmp=openmp, with_blocking = ps.create_kernel(assignments, cpu_blocking=(8, 16, 4), cpu_openmp=openmp,
cpu_vectorize_info=vectorization).compile() cpu_vectorize_info=vectorization).compile()
with_blocking_only_over_y = ps.create_kernel(assignments, cpu_blocking=(0, 16, 0), cpu_openmp=openmp,
cpu_vectorize_info=vectorization).compile()
without_blocking = ps.create_kernel(assignments).compile() without_blocking = ps.create_kernel(assignments).compile()
print(" openmp {}, vectorization {}".format(openmp, vectorization))
only_omp = ps.create_kernel(assignments, cpu_openmp=2).compile()
print(f" openmp {openmp}, vectorization {vectorization}")
dst_arr = np.zeros_like(src_arr) dst_arr = np.zeros_like(src_arr)
dst2_arr = np.zeros_like(src_arr)
dst3_arr = np.zeros_like(src_arr)
ref_arr = np.zeros_like(src_arr) ref_arr = np.zeros_like(src_arr)
np.copyto(src_arr, np.random.rand(*src_arr.shape)) np.copyto(src_arr, np.random.rand(*src_arr.shape))
with_blocking(src=src_arr, dst=dst_arr) with_blocking(src=src_arr, dst=dst_arr)
with_blocking_only_over_y(src=src_arr, dst=dst2_arr)
without_blocking(src=src_arr, dst=ref_arr) without_blocking(src=src_arr, dst=ref_arr)
only_omp(src=src_arr, dst=dst3_arr)
np.testing.assert_almost_equal(ref_arr, dst_arr) np.testing.assert_almost_equal(ref_arr, dst_arr)
np.testing.assert_almost_equal(ref_arr, dst2_arr)
np.testing.assert_almost_equal(ref_arr, dst3_arr)
def test_jacobi3d_var_size(): def test_jacobi3d_var_size():
...@@ -59,3 +70,11 @@ def test_jacobi3d_fixed_size(): ...@@ -59,3 +70,11 @@ def test_jacobi3d_fixed_size():
arr = np.empty([8*4, 16*2, 4*3]) arr = np.empty([8*4, 16*2, 4*3])
src, dst = ps.fields("src, dst: double[3D]", src=arr, dst=arr) src, dst = ps.fields("src, dst: double[3D]", src=arr, dst=arr)
check_equivalence(jacobi(dst, src), arr) check_equivalence(jacobi(dst, src), arr)
def test_jacobi3d_fixed_field_size():
src, dst = ps.fields("src, dst: double[3, 5, 6]", layout='c')
print("Fixed Field Size: Smaller than block sizes")
arr = np.empty([3, 5, 6])
check_equivalence(jacobi(dst, src), arr)
...@@ -4,14 +4,18 @@ import pystencils as ps ...@@ -4,14 +4,18 @@ import pystencils as ps
def test_blocking_staggered(): def test_blocking_staggered():
f, stag = ps.fields("f, stag(3): double[3D]") f = ps.fields("f: double[3D]")
stag = ps.fields("stag(3): double[3D]", field_type=ps.FieldType.STAGGERED)
terms = [ terms = [
f[0, 0, 0] - f[-1, 0, 0], f[0, 0, 0] - f[-1, 0, 0],
f[0, 0, 0] - f[0, -1, 0], f[0, 0, 0] - f[0, -1, 0],
f[0, 0, 0] - f[0, 0, -1], f[0, 0, 0] - f[0, 0, -1],
] ]
kernel = ps.create_staggered_kernel(stag, terms, cpu_blocking=(3, 16, 8)).compile() assignments = [ps.Assignment(stag.staggered_access(d), terms[i]) for i, d in enumerate(stag.staggered_stencil)]
reference_kernel = ps.create_staggered_kernel(stag, terms).compile() reference_kernel = ps.create_staggered_kernel(assignments)
print(ps.show_code(reference_kernel))
reference_kernel = reference_kernel.compile()
kernel = ps.create_staggered_kernel(assignments, cpu_blocking=(3, 16, 8)).compile()
print(ps.show_code(kernel.ast)) print(ps.show_code(kernel.ast))
f_arr = np.random.rand(80, 33, 19) f_arr = np.random.rand(80, 33, 19)
......
...@@ -2,11 +2,15 @@ import os ...@@ -2,11 +2,15 @@ import os
from tempfile import TemporaryDirectory from tempfile import TemporaryDirectory
import numpy as np import numpy as np
import pytest
import pystencils
from pystencils import Assignment, create_kernel from pystencils import Assignment, create_kernel
from pystencils.boundaries import BoundaryHandling, Neumann, add_neumann_boundary from pystencils.boundaries import BoundaryHandling, Dirichlet, Neumann, add_neumann_boundary
from pystencils.datahandling import SerialDataHandling from pystencils.datahandling import SerialDataHandling
from pystencils.enums import Target
from pystencils.slicing import slice_from_direction from pystencils.slicing import slice_from_direction
from pystencils.timeloop import TimeLoop
def test_kernel_vs_copy_boundary(): def test_kernel_vs_copy_boundary():
...@@ -83,5 +87,160 @@ def test_kernel_vs_copy_boundary(): ...@@ -83,5 +87,160 @@ def test_kernel_vs_copy_boundary():
np.testing.assert_almost_equal(python_copy_result, handling_result) np.testing.assert_almost_equal(python_copy_result, handling_result)
with TemporaryDirectory() as tmp_dir: with TemporaryDirectory() as tmp_dir:
pytest.importorskip('pyevtk')
boundary_handling.geometry_to_vtk(file_name=os.path.join(tmp_dir, 'test_output1'), ghost_layers=False) boundary_handling.geometry_to_vtk(file_name=os.path.join(tmp_dir, 'test_output1'), ghost_layers=False)
boundary_handling.geometry_to_vtk(file_name=os.path.join(tmp_dir, 'test_output2'), ghost_layers=True) boundary_handling.geometry_to_vtk(file_name=os.path.join(tmp_dir, 'test_output2'), ghost_layers=True)
boundaries = list(boundary_handling._boundary_object_to_boundary_info.keys()) + ['domain']
boundary_handling.geometry_to_vtk(file_name=os.path.join(tmp_dir, 'test_output3'),
boundaries=boundaries[0], ghost_layers=False)
def test_boundary_gpu():
pytest.importorskip('cupy')
dh = SerialDataHandling(domain_size=(7, 7), default_target=Target.GPU)
src = dh.add_array('src')
dh.fill("src", 0.0, ghost_layers=True)
dh.fill("src", 1.0, ghost_layers=False)
src_cpu = dh.add_array('src_cpu', gpu=False)
dh.fill("src_cpu", 0.0, ghost_layers=True)
dh.fill("src_cpu", 1.0, ghost_layers=False)
boundary_stencil = [(1, 0), (-1, 0), (0, 1), (0, -1)]
boundary_handling_cpu = BoundaryHandling(dh, src_cpu.name, boundary_stencil,
name="boundary_handling_cpu", target=Target.CPU)
boundary_handling = BoundaryHandling(dh, src.name, boundary_stencil,
name="boundary_handling_gpu", target=Target.GPU)
neumann = Neumann()
for d in ('N', 'S', 'W', 'E'):
boundary_handling.set_boundary(neumann, slice_from_direction(d, dim=2))
boundary_handling_cpu.set_boundary(neumann, slice_from_direction(d, dim=2))
boundary_handling.prepare()
boundary_handling_cpu.prepare()
boundary_handling_cpu()
dh.all_to_gpu()
boundary_handling()
dh.all_to_cpu()
np.testing.assert_almost_equal(dh.cpu_arrays["src_cpu"], dh.cpu_arrays["src"])
def test_boundary_utility():
dh = SerialDataHandling(domain_size=(7, 7))
src = dh.add_array('src')
dh.fill("src", 0.0, ghost_layers=True)
boundary_stencil = [(1, 0), (-1, 0), (0, 1), (0, -1)]
boundary_handling = BoundaryHandling(dh, src.name, boundary_stencil,
name="boundary_handling", target=Target.CPU)
neumann = Neumann()
dirichlet = Dirichlet(2)
for d in ('N', 'S', 'W', 'E'):
boundary_handling.set_boundary(neumann, slice_from_direction(d, dim=2))
boundary_handling.set_boundary(neumann, (slice(2, 4, None), slice(2, 4, None)))
boundary_handling.prepare()
assert boundary_handling.get_flag(boundary_handling.boundary_objects[0]) == 2
assert boundary_handling.shape == dh.shape
assert boundary_handling.flag_array_name == 'boundary_handlingFlags'
mask_neumann = boundary_handling.get_mask((slice(0, 7), slice(0, 7)), boundary_handling.boundary_objects[0])
np.testing.assert_almost_equal(mask_neumann[1:3, 1:3], 2)
mask_domain = boundary_handling.get_mask((slice(0, 7), slice(0, 7)), "domain")
assert np.sum(mask_domain) == 7 ** 2 - 4
def set_sphere(x, y):
mid = (4, 4)
radius = 2
return (x - mid[0]) ** 2 + (y - mid[1]) ** 2 < radius ** 2
boundary_handling.set_boundary(dirichlet, mask_callback=set_sphere, force_flag_value=4)
mask_dirichlet = boundary_handling.get_mask((slice(0, 7), slice(0, 7)), boundary_handling.boundary_objects[1])
assert np.sum(mask_dirichlet) == 48
assert boundary_handling.set_boundary("domain") == 1
assert boundary_handling.set_boundary(dirichlet, mask_callback=set_sphere, force_flag_value=8, replace=False) == 4
assert boundary_handling.set_boundary(dirichlet, force_flag_value=16, replace=False) == 4
assert boundary_handling.set_boundary_where_flag_is_set(boundary_handling.boundary_objects[0], 16) == 16
def test_add_fix_steps():
dh = SerialDataHandling(domain_size=(7, 7))
src = dh.add_array('src')
dh.fill("src", 0.0, ghost_layers=True)
dh.fill("src", 1.0, ghost_layers=False)
boundary_stencil = [(1, 0), (-1, 0), (0, 1), (0, -1)]
boundary_handling = BoundaryHandling(dh, src.name, boundary_stencil,
name="boundary_handling", target=pystencils.Target.CPU)
neumann = Neumann()
for d in ('N', 'S', 'W', 'E'):
boundary_handling.set_boundary(neumann, slice_from_direction(d, dim=2))
timeloop = TimeLoop(steps=1)
boundary_handling.add_fixed_steps(timeloop)
timeloop.run()
assert np.sum(dh.cpu_arrays['src']) == 7 * 7 + 7 * 4
def test_boundary_data_setter():
dh = SerialDataHandling(domain_size=(7, 7))
src = dh.add_array('src')
dh.fill("src", 0.0, ghost_layers=True)
dh.fill("src", 1.0, ghost_layers=False)
boundary_stencil = [(1, 0), (-1, 0), (0, 1), (0, -1)]
boundary_handling = BoundaryHandling(dh, src.name, boundary_stencil,
name="boundary_handling", target=Target.CPU)
neumann = Neumann()
for d in 'N':
boundary_handling.set_boundary(neumann, slice_from_direction(d, dim=2))
boundary_handling.prepare()
for b in dh.iterate(ghost_layers=True):
index_array_bd = b[boundary_handling._index_array_name]
data_setter = index_array_bd.boundary_object_to_data_setter[boundary_handling.boundary_objects[0]]
y_pos = data_setter.boundary_cell_positions(1)
assert all(y_pos == 5.5)
assert np.all(data_setter.link_offsets() == [0, -1])
assert np.all(data_setter.link_positions(1) == 6.)
@pytest.mark.parametrize('with_indices', ('with_indices', False))
def test_dirichlet(with_indices):
value = (1, 20, 3) if with_indices else 1
dh = SerialDataHandling(domain_size=(7, 7))
src = dh.add_array('src', values_per_cell=3 if with_indices else 1)
dh.cpu_arrays.src[...] = np.random.rand(*src.shape)
boundary_stencil = [(1, 0), (-1, 0), (0, 1), (0, -1)]
boundary_handling = BoundaryHandling(dh, src.name, boundary_stencil)
dirichlet = Dirichlet(value)
assert dirichlet.name == 'Dirichlet'
dirichlet.name = "wall"
assert dirichlet.name == 'wall'
for d in ('N', 'S', 'W', 'E'):
boundary_handling.set_boundary(dirichlet, slice_from_direction(d, dim=2))
boundary_handling()
assert all([np.allclose(a, np.array(value)) for a in dh.cpu_arrays.src[1:-2, 0]])
assert all([np.allclose(a, np.array(value)) for a in dh.cpu_arrays.src[1:-2, -1]])
assert all([np.allclose(a, np.array(value)) for a in dh.cpu_arrays.src[0, 1:-2]])
assert all([np.allclose(a, np.array(value)) for a in dh.cpu_arrays.src[-1, 1:-2]])
import numpy as np
from itertools import product
import pystencils.boundaries.createindexlist as cil
import pytest
@pytest.mark.parametrize('single_link', [False, True])
@pytest.mark.skipif(not cil.cython_funcs_available, reason='Cython functions are not available')
def test_equivalence_cython_python_version(single_link):
# D2Q9
stencil_2d = tuple((x, y) for x, y in product([-1, 0, 1], [-1, 0, 1]))
# D3Q19
stencil_3d = tuple(
(x, y, z) for x, y, z in product([-1, 0, 1], [-1, 0, 1], [-1, 0, 1]) if abs(x) + abs(y) + abs(z) < 3)
for dtype in [int, np.int16, np.uint32]:
fluid_mask = dtype(1)
mask = dtype(2)
flag_field_2d = np.ones([15, 16], dtype=dtype) * fluid_mask
flag_field_3d = np.ones([15, 16, 17], dtype=dtype) * fluid_mask
flag_field_2d[0, :] = mask
flag_field_2d[-1, :] = mask
flag_field_2d[7, 7] = mask
flag_field_3d[0, :, :] = mask
flag_field_3d[-1, :, :] = mask
flag_field_3d[7, 7, 7] = mask
result_python_2d = cil._create_index_list_python(flag_field_2d, mask, fluid_mask,
stencil_2d, single_link, True, 1)
result_python_3d = cil._create_index_list_python(flag_field_3d, mask, fluid_mask,
stencil_3d, single_link, True, 1)
result_cython_2d = cil.create_boundary_index_list(flag_field_2d, stencil_2d, mask,
fluid_mask, 1, True, single_link)
result_cython_3d = cil.create_boundary_index_list(flag_field_3d, stencil_3d, mask,
fluid_mask, 1, True, 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('single_link', [False, True])
@pytest.mark.skipif(not cil.cython_funcs_available, reason='Cython functions are not available')
def test_equivalence_cell_idx_list_cython_python_version(single_link):
# D2Q9
stencil_2d = tuple((x, y) for x, y in product([-1, 0, 1], [-1, 0, 1]))
# D3Q19
stencil_3d = tuple(
(x, y, z) for x, y, z in product([-1, 0, 1], [-1, 0, 1], [-1, 0, 1]) if abs(x) + abs(y) + abs(z) < 3)
for dtype in [int, np.int16, np.uint32]:
fluid_mask = dtype(1)
mask = dtype(2)
flag_field_2d = np.ones([15, 16], dtype=dtype) * fluid_mask
flag_field_3d = np.ones([15, 16, 17], dtype=dtype) * fluid_mask
flag_field_2d[0, :] = mask
flag_field_2d[-1, :] = mask
flag_field_2d[7, 7] = mask
flag_field_3d[0, :, :] = mask
flag_field_3d[-1, :, :] = mask
flag_field_3d[7, 7, 7] = mask
result_python_2d = cil._create_index_list_python(flag_field_2d, mask, fluid_mask,
stencil_2d, single_link, False)
result_python_3d = cil._create_index_list_python(flag_field_3d, mask, fluid_mask,
stencil_3d, single_link, False)
result_cython_2d = cil.create_boundary_index_list(flag_field_2d, stencil_2d, mask, fluid_mask, None,
False, single_link)
result_cython_3d = cil.create_boundary_index_list(flag_field_3d, stencil_3d, mask, fluid_mask, None,
False, 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 = tuple((x, y) for x, y in product([-1, 0, 1], [-1, 0, 1]))
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
...@@ -2,6 +2,7 @@ ...@@ -2,6 +2,7 @@
import numpy as np import numpy as np
import pystencils as ps
from pystencils import Assignment, Field, FieldType, create_kernel from pystencils import Assignment, Field, FieldType, create_kernel
from pystencils.field import create_numpy_array_with_layout, layout_string_to_tuple from pystencils.field import create_numpy_array_with_layout, layout_string_to_tuple
from pystencils.slicing import ( from pystencils.slicing import (
...@@ -19,9 +20,9 @@ def _generate_fields(dt=np.uint64, num_directions=1, layout='numpy'): ...@@ -19,9 +20,9 @@ def _generate_fields(dt=np.uint64, num_directions=1, layout='numpy'):
fields = [] fields = []
for size in field_sizes: for size in field_sizes:
field_layout = layout_string_to_tuple(layout, len(size)) field_layout = layout_string_to_tuple(layout, len(size))
src_arr = create_numpy_array_with_layout(size, field_layout) src_arr = create_numpy_array_with_layout(size, field_layout, dtype=dt)
array_data = np.reshape(np.arange(1, int(np.prod(size)+1)), size) array_data = np.reshape(np.arange(1, int(np.prod(size) + 1)), size)
# Use flat iterator to input data into the array # Use flat iterator to input data into the array
src_arr.flat = add_ghost_layers(array_data, index_dimensions=1 if num_directions > 1 else 0).astype(dt).flat src_arr.flat = add_ghost_layers(array_data, index_dimensions=1 if num_directions > 1 else 0).astype(dt).flat
dst_arr = np.zeros(src_arr.shape, dtype=dt) dst_arr = np.zeros(src_arr.shape, dtype=dt)
...@@ -40,13 +41,18 @@ def test_full_scalar_field(): ...@@ -40,13 +41,18 @@ def test_full_scalar_field():
field_type=FieldType.BUFFER, dtype=src_arr.dtype) field_type=FieldType.BUFFER, dtype=src_arr.dtype)
pack_eqs = [Assignment(buffer.center(), src_field.center())] pack_eqs = [Assignment(buffer.center(), src_field.center())]
pack_code = create_kernel(pack_eqs, data_type={'src_field': src_arr.dtype, 'buffer': buffer.dtype}) config = ps.CreateKernelConfig(data_type={'src_field': src_arr.dtype, 'buffer': buffer.dtype})
pack_code = create_kernel(pack_eqs, config=config)
code = ps.get_code_str(pack_code)
ps.show_code(pack_code)
pack_kernel = pack_code.compile() pack_kernel = pack_code.compile()
pack_kernel(buffer=buffer_arr, src_field=src_arr) pack_kernel(buffer=buffer_arr, src_field=src_arr)
unpack_eqs = [Assignment(dst_field.center(), buffer.center())] unpack_eqs = [Assignment(dst_field.center(), buffer.center())]
unpack_code = create_kernel(unpack_eqs, data_type={'dst_field': dst_arr.dtype, 'buffer': buffer.dtype})
config = ps.CreateKernelConfig(data_type={'dst_field': dst_arr.dtype, 'buffer': buffer.dtype})
unpack_code = create_kernel(unpack_eqs, config=config)
unpack_kernel = unpack_code.compile() unpack_kernel = unpack_code.compile()
unpack_kernel(dst_field=dst_arr, buffer=buffer_arr) unpack_kernel(dst_field=dst_arr, buffer=buffer_arr)
...@@ -70,14 +76,18 @@ def test_field_slice(): ...@@ -70,14 +76,18 @@ def test_field_slice():
field_type=FieldType.BUFFER, dtype=src_arr.dtype) field_type=FieldType.BUFFER, dtype=src_arr.dtype)
pack_eqs = [Assignment(buffer.center(), src_field.center())] pack_eqs = [Assignment(buffer.center(), src_field.center())]
pack_code = create_kernel(pack_eqs, data_type={'src_field': src_arr.dtype, 'buffer': buffer.dtype})
config = ps.CreateKernelConfig(data_type={'src_field': src_arr.dtype, 'buffer': buffer.dtype})
pack_code = create_kernel(pack_eqs, config=config)
pack_kernel = pack_code.compile() pack_kernel = pack_code.compile()
pack_kernel(buffer=bufferArr, src_field=src_arr[pack_slice]) pack_kernel(buffer=bufferArr, src_field=src_arr[pack_slice])
# Unpack into ghost layer of dst_field in N direction # Unpack into ghost layer of dst_field in N direction
unpack_eqs = [Assignment(dst_field.center(), buffer.center())] unpack_eqs = [Assignment(dst_field.center(), buffer.center())]
unpack_code = create_kernel(unpack_eqs, data_type={'dst_field': dst_arr.dtype, 'buffer': buffer.dtype})
config = ps.CreateKernelConfig(data_type={'dst_field': dst_arr.dtype, 'buffer': buffer.dtype})
unpack_code = create_kernel(unpack_eqs, config=config)
unpack_kernel = unpack_code.compile() unpack_kernel = unpack_code.compile()
unpack_kernel(buffer=bufferArr, dst_field=dst_arr[unpack_slice]) unpack_kernel(buffer=bufferArr, dst_field=dst_arr[unpack_slice])
...@@ -102,7 +112,8 @@ def test_all_cell_values(): ...@@ -102,7 +112,8 @@ def test_all_cell_values():
eq = Assignment(buffer(idx), src_field(idx)) eq = Assignment(buffer(idx), src_field(idx))
pack_eqs.append(eq) pack_eqs.append(eq)
pack_code = create_kernel(pack_eqs, data_type={'src_field': src_arr.dtype, 'buffer': buffer.dtype}) config = ps.CreateKernelConfig(data_type={'src_field': src_arr.dtype, 'buffer': buffer.dtype})
pack_code = create_kernel(pack_eqs, config=config)
pack_kernel = pack_code.compile() pack_kernel = pack_code.compile()
pack_kernel(buffer=bufferArr, src_field=src_arr) pack_kernel(buffer=bufferArr, src_field=src_arr)
...@@ -112,7 +123,8 @@ def test_all_cell_values(): ...@@ -112,7 +123,8 @@ def test_all_cell_values():
eq = Assignment(dst_field(idx), buffer(idx)) eq = Assignment(dst_field(idx), buffer(idx))
unpack_eqs.append(eq) unpack_eqs.append(eq)
unpack_code = create_kernel(unpack_eqs, data_type={'dst_field': dst_arr.dtype, 'buffer': buffer.dtype}) config = ps.CreateKernelConfig(data_type={'dst_field': dst_arr.dtype, 'buffer': buffer.dtype})
unpack_code = create_kernel(unpack_eqs, config=config)
unpack_kernel = unpack_code.compile() unpack_kernel = unpack_code.compile()
unpack_kernel(buffer=bufferArr, dst_field=dst_arr) unpack_kernel(buffer=bufferArr, dst_field=dst_arr)
...@@ -138,7 +150,8 @@ def test_subset_cell_values(): ...@@ -138,7 +150,8 @@ def test_subset_cell_values():
eq = Assignment(buffer(buffer_idx), src_field(cell_idx)) eq = Assignment(buffer(buffer_idx), src_field(cell_idx))
pack_eqs.append(eq) pack_eqs.append(eq)
pack_code = create_kernel(pack_eqs, data_type={'src_field': src_arr.dtype, 'buffer': buffer.dtype}) config = ps.CreateKernelConfig(data_type={'src_field': src_arr.dtype, 'buffer': buffer.dtype})
pack_code = create_kernel(pack_eqs, config=config)
pack_kernel = pack_code.compile() pack_kernel = pack_code.compile()
pack_kernel(buffer=bufferArr, src_field=src_arr) pack_kernel(buffer=bufferArr, src_field=src_arr)
...@@ -148,7 +161,8 @@ def test_subset_cell_values(): ...@@ -148,7 +161,8 @@ def test_subset_cell_values():
eq = Assignment(dst_field(cell_idx), buffer(buffer_idx)) eq = Assignment(dst_field(cell_idx), buffer(buffer_idx))
unpack_eqs.append(eq) unpack_eqs.append(eq)
unpack_code = create_kernel(unpack_eqs, data_type={'dst_field': dst_arr.dtype, 'buffer': buffer.dtype}) config = ps.CreateKernelConfig(data_type={'dst_field': dst_arr.dtype, 'buffer': buffer.dtype})
unpack_code = create_kernel(unpack_eqs, config=config)
unpack_kernel = unpack_code.compile() unpack_kernel = unpack_code.compile()
unpack_kernel(buffer=bufferArr, dst_field=dst_arr) unpack_kernel(buffer=bufferArr, dst_field=dst_arr)
...@@ -173,7 +187,8 @@ def test_field_layouts(): ...@@ -173,7 +187,8 @@ def test_field_layouts():
eq = Assignment(buffer(idx), src_field(idx)) eq = Assignment(buffer(idx), src_field(idx))
pack_eqs.append(eq) pack_eqs.append(eq)
pack_code = create_kernel(pack_eqs, data_type={'src_field': src_arr.dtype, 'buffer': buffer.dtype}) config = ps.CreateKernelConfig(data_type={'src_field': src_arr.dtype, 'buffer': buffer.dtype})
pack_code = create_kernel(pack_eqs, config=config)
pack_kernel = pack_code.compile() pack_kernel = pack_code.compile()
pack_kernel(buffer=bufferArr, src_field=src_arr) pack_kernel(buffer=bufferArr, src_field=src_arr)
...@@ -183,6 +198,62 @@ def test_field_layouts(): ...@@ -183,6 +198,62 @@ def test_field_layouts():
eq = Assignment(dst_field(idx), buffer(idx)) eq = Assignment(dst_field(idx), buffer(idx))
unpack_eqs.append(eq) unpack_eqs.append(eq)
unpack_code = create_kernel(unpack_eqs, data_type={'dst_field': dst_arr.dtype, 'buffer': buffer.dtype}) config = ps.CreateKernelConfig(data_type={'dst_field': dst_arr.dtype, 'buffer': buffer.dtype})
unpack_code = create_kernel(unpack_eqs, config=config)
unpack_kernel = unpack_code.compile() unpack_kernel = unpack_code.compile()
unpack_kernel(buffer=bufferArr, dst_field=dst_arr) unpack_kernel(buffer=bufferArr, dst_field=dst_arr)
def test_iteration_slices():
num_cell_values = 19
dt = np.uint64
fields = _generate_fields(dt=dt, num_directions=num_cell_values)
for (src_arr, dst_arr, bufferArr) in fields:
spatial_dimensions = len(src_arr.shape) - 1
# src_field = Field.create_from_numpy_array("src_field", src_arr, index_dimensions=1)
# dst_field = Field.create_from_numpy_array("dst_field", dst_arr, index_dimensions=1)
src_field = Field.create_generic("src_field", spatial_dimensions, index_shape=(num_cell_values,), dtype=dt)
dst_field = Field.create_generic("dst_field", spatial_dimensions, index_shape=(num_cell_values,), dtype=dt)
buffer = Field.create_generic("buffer", spatial_dimensions=1, index_dimensions=1,
field_type=FieldType.BUFFER, dtype=src_arr.dtype)
pack_eqs = []
# Since we are packing all cell values for all cells, then
# the buffer index is equivalent to the field index
for idx in range(num_cell_values):
eq = Assignment(buffer(idx), src_field(idx))
pack_eqs.append(eq)
dim = src_field.spatial_dimensions
# Pack only the leftmost slice, only every second cell
pack_slice = (slice(None, None, 2),) * (dim - 1) + (0,)
# Fill the entire array with data
src_arr[(slice(None, None, 1),) * dim] = np.arange(num_cell_values)
dst_arr.fill(0)
config = ps.CreateKernelConfig(iteration_slice=pack_slice,
data_type={'src_field': src_arr.dtype, 'buffer': buffer.dtype})
pack_code = create_kernel(pack_eqs, config=config)
pack_kernel = pack_code.compile()
pack_kernel(buffer=bufferArr, src_field=src_arr)
unpack_eqs = []
for idx in range(num_cell_values):
eq = Assignment(dst_field(idx), buffer(idx))
unpack_eqs.append(eq)
config = ps.CreateKernelConfig(iteration_slice=pack_slice,
data_type={'dst_field': dst_arr.dtype, 'buffer': buffer.dtype})
unpack_code = create_kernel(unpack_eqs, config=config)
unpack_kernel = unpack_code.compile()
unpack_kernel(buffer=bufferArr, dst_field=dst_arr)
# Check if only every second entry of the leftmost slice has been copied
np.testing.assert_equal(dst_arr[pack_slice], src_arr[pack_slice])
np.testing.assert_equal(dst_arr[(slice(1, None, 2),) * (dim - 1) + (0,)], 0)
np.testing.assert_equal(dst_arr[(slice(None, None, 1),) * (dim - 1) + (slice(1, None),)], 0)
"""Tests for the (un)packing (from)to buffers on a CUDA GPU.""" """Tests for the (un)packing (from)to buffers on a CUDA GPU."""
from dataclasses import replace
import numpy as np import numpy as np
import pytest import pytest
from pystencils import Assignment, Field, FieldType import pystencils
from pystencils import Assignment, Field, FieldType, Target, CreateKernelConfig, create_kernel, fields
from pystencils.bit_masks import flag_cond
from pystencils.field import create_numpy_array_with_layout, layout_string_to_tuple from pystencils.field import create_numpy_array_with_layout, layout_string_to_tuple
from pystencils.gpucuda import create_cuda_kernel, make_python_function
from pystencils.slicing import ( from pystencils.slicing import (
add_ghost_layers, get_ghost_region_slice, get_slice_before_ghost_layer) add_ghost_layers, get_ghost_region_slice, get_slice_before_ghost_layer)
from pystencils.stencil import direction_string_to_offset from pystencils.stencil import direction_string_to_offset
try: try:
# noinspection PyUnresolvedReferences # noinspection PyUnresolvedReferences
import pycuda.autoinit import cupy as cp
import pycuda.gpuarray as gpuarray
except ImportError: except ImportError:
pass pass
...@@ -22,6 +23,7 @@ FIELD_SIZES = [(4, 3), (9, 3, 7)] ...@@ -22,6 +23,7 @@ FIELD_SIZES = [(4, 3), (9, 3, 7)]
def _generate_fields(dt=np.uint8, stencil_directions=1, layout='numpy'): def _generate_fields(dt=np.uint8, stencil_directions=1, layout='numpy'):
pytest.importorskip('cupy')
field_sizes = FIELD_SIZES field_sizes = FIELD_SIZES
if stencil_directions > 1: if stencil_directions > 1:
field_sizes = [s + (stencil_directions,) for s in field_sizes] field_sizes = [s + (stencil_directions,) for s in field_sizes]
...@@ -36,15 +38,15 @@ def _generate_fields(dt=np.uint8, stencil_directions=1, layout='numpy'): ...@@ -36,15 +38,15 @@ def _generate_fields(dt=np.uint8, stencil_directions=1, layout='numpy'):
src_arr.flat = add_ghost_layers(array_data, src_arr.flat = add_ghost_layers(array_data,
index_dimensions=1 if stencil_directions > 1 else 0).astype(dt).flat index_dimensions=1 if stencil_directions > 1 else 0).astype(dt).flat
gpu_src_arr = gpuarray.to_gpu(src_arr) gpu_src_arr = cp.asarray(src_arr)
gpu_dst_arr = gpuarray.zeros_like(gpu_src_arr) gpu_dst_arr = cp.zeros_like(gpu_src_arr)
gpu_buffer_arr = gpuarray.zeros(np.prod(src_arr.shape), dtype=dt) size = int(np.prod(src_arr.shape))
gpu_buffer_arr = cp.zeros(size, dtype=dt)
fields.append((src_arr, gpu_src_arr, gpu_dst_arr, gpu_buffer_arr)) fields.append((src_arr, gpu_src_arr, gpu_dst_arr, gpu_buffer_arr))
return fields return fields
@pytest.mark.gpu
def test_full_scalar_field(): def test_full_scalar_field():
"""Tests fully (un)packing a scalar field (from)to a GPU buffer.""" """Tests fully (un)packing a scalar field (from)to a GPU buffer."""
fields = _generate_fields() fields = _generate_fields()
...@@ -56,16 +58,20 @@ def test_full_scalar_field(): ...@@ -56,16 +58,20 @@ def test_full_scalar_field():
pack_eqs = [Assignment(buffer.center(), src_field.center())] pack_eqs = [Assignment(buffer.center(), src_field.center())]
pack_types = {'src_field': gpu_src_arr.dtype, 'buffer': gpu_buffer_arr.dtype} pack_types = {'src_field': gpu_src_arr.dtype, 'buffer': gpu_buffer_arr.dtype}
pack_code = create_cuda_kernel(pack_eqs, type_info=pack_types)
pack_kernel = make_python_function(pack_code) config = CreateKernelConfig(target=pystencils.Target.GPU, data_type=pack_types)
pack_ast = create_kernel(pack_eqs, config=config)
pack_kernel = pack_ast.compile()
pack_kernel(buffer=gpu_buffer_arr, src_field=gpu_src_arr) pack_kernel(buffer=gpu_buffer_arr, src_field=gpu_src_arr)
unpack_eqs = [Assignment(dst_field.center(), buffer.center())] unpack_eqs = [Assignment(dst_field.center(), buffer.center())]
unpack_types = {'dst_field': gpu_dst_arr.dtype, 'buffer': gpu_buffer_arr.dtype} unpack_types = {'dst_field': gpu_dst_arr.dtype, 'buffer': gpu_buffer_arr.dtype}
unpack_code = create_cuda_kernel(unpack_eqs, type_info=unpack_types)
unpack_kernel = make_python_function(unpack_code) config = CreateKernelConfig(target=pystencils.Target.GPU, data_type=unpack_types)
unpack_ast = create_kernel(unpack_eqs, config=config)
unpack_kernel = unpack_ast.compile()
unpack_kernel(dst_field=gpu_dst_arr, buffer=gpu_buffer_arr) unpack_kernel(dst_field=gpu_dst_arr, buffer=gpu_buffer_arr)
dst_arr = gpu_dst_arr.get() dst_arr = gpu_dst_arr.get()
...@@ -73,7 +79,6 @@ def test_full_scalar_field(): ...@@ -73,7 +79,6 @@ def test_full_scalar_field():
np.testing.assert_equal(src_arr, dst_arr) np.testing.assert_equal(src_arr, dst_arr)
@pytest.mark.gpu
def test_field_slice(): def test_field_slice():
"""Tests (un)packing slices of a scalar field (from)to a buffer.""" """Tests (un)packing slices of a scalar field (from)to a buffer."""
fields = _generate_fields() fields = _generate_fields()
...@@ -91,17 +96,21 @@ def test_field_slice(): ...@@ -91,17 +96,21 @@ def test_field_slice():
pack_eqs = [Assignment(buffer.center(), src_field.center())] pack_eqs = [Assignment(buffer.center(), src_field.center())]
pack_types = {'src_field': gpu_src_arr.dtype, 'buffer': gpu_buffer_arr.dtype} pack_types = {'src_field': gpu_src_arr.dtype, 'buffer': gpu_buffer_arr.dtype}
pack_code = create_cuda_kernel(pack_eqs, type_info=pack_types)
pack_kernel = make_python_function(pack_code) config = CreateKernelConfig(target=pystencils.Target.GPU, data_type=pack_types)
pack_ast = create_kernel(pack_eqs, config=config)
pack_kernel = pack_ast.compile()
pack_kernel(buffer=gpu_buffer_arr, src_field=gpu_src_arr[pack_slice]) pack_kernel(buffer=gpu_buffer_arr, src_field=gpu_src_arr[pack_slice])
# Unpack into ghost layer of dst_field in N direction # Unpack into ghost layer of dst_field in N direction
unpack_eqs = [Assignment(dst_field.center(), buffer.center())] unpack_eqs = [Assignment(dst_field.center(), buffer.center())]
unpack_types = {'dst_field': gpu_dst_arr.dtype, 'buffer': gpu_buffer_arr.dtype} unpack_types = {'dst_field': gpu_dst_arr.dtype, 'buffer': gpu_buffer_arr.dtype}
unpack_code = create_cuda_kernel(unpack_eqs, type_info=unpack_types)
unpack_kernel = make_python_function(unpack_code) config = CreateKernelConfig(target=pystencils.Target.GPU, data_type=unpack_types)
unpack_ast = create_kernel(unpack_eqs, config=config)
unpack_kernel = unpack_ast.compile()
unpack_kernel(buffer=gpu_buffer_arr, dst_field=gpu_dst_arr[unpack_slice]) unpack_kernel(buffer=gpu_buffer_arr, dst_field=gpu_dst_arr[unpack_slice])
dst_arr = gpu_dst_arr.get() dst_arr = gpu_dst_arr.get()
...@@ -109,7 +118,6 @@ def test_field_slice(): ...@@ -109,7 +118,6 @@ def test_field_slice():
np.testing.assert_equal(src_arr[pack_slice], dst_arr[unpack_slice]) np.testing.assert_equal(src_arr[pack_slice], dst_arr[unpack_slice])
@pytest.mark.gpu
def test_all_cell_values(): def test_all_cell_values():
"""Tests (un)packing all cell values of the a field (from)to a buffer.""" """Tests (un)packing all cell values of the a field (from)to a buffer."""
num_cell_values = 7 num_cell_values = 7
...@@ -128,8 +136,11 @@ def test_all_cell_values(): ...@@ -128,8 +136,11 @@ def test_all_cell_values():
pack_eqs.append(eq) pack_eqs.append(eq)
pack_types = {'src_field': gpu_src_arr.dtype, 'buffer': gpu_buffer_arr.dtype} pack_types = {'src_field': gpu_src_arr.dtype, 'buffer': gpu_buffer_arr.dtype}
pack_code = create_cuda_kernel(pack_eqs, type_info=pack_types)
pack_kernel = make_python_function(pack_code) config = CreateKernelConfig(target=pystencils.Target.GPU, data_type=pack_types)
pack_code = create_kernel(pack_eqs, config=config)
pack_kernel = pack_code.compile()
pack_kernel(buffer=gpu_buffer_arr, src_field=gpu_src_arr) pack_kernel(buffer=gpu_buffer_arr, src_field=gpu_src_arr)
unpack_eqs = [] unpack_eqs = []
...@@ -139,8 +150,10 @@ def test_all_cell_values(): ...@@ -139,8 +150,10 @@ def test_all_cell_values():
unpack_eqs.append(eq) unpack_eqs.append(eq)
unpack_types = {'dst_field': gpu_dst_arr.dtype, 'buffer': gpu_buffer_arr.dtype} unpack_types = {'dst_field': gpu_dst_arr.dtype, 'buffer': gpu_buffer_arr.dtype}
unpack_code = create_cuda_kernel(unpack_eqs, type_info=unpack_types)
unpack_kernel = make_python_function(unpack_code) config = CreateKernelConfig(target=pystencils.Target.GPU, data_type=unpack_types)
unpack_ast = create_kernel(unpack_eqs, config=config)
unpack_kernel = unpack_ast.compile()
unpack_kernel(buffer=gpu_buffer_arr, dst_field=gpu_dst_arr) unpack_kernel(buffer=gpu_buffer_arr, dst_field=gpu_dst_arr)
dst_arr = gpu_dst_arr.get() dst_arr = gpu_dst_arr.get()
...@@ -148,9 +161,8 @@ def test_all_cell_values(): ...@@ -148,9 +161,8 @@ def test_all_cell_values():
np.testing.assert_equal(src_arr, dst_arr) np.testing.assert_equal(src_arr, dst_arr)
@pytest.mark.gpu
def test_subset_cell_values(): def test_subset_cell_values():
"""Tests (un)packing a subset of cell values of the a field (from)to a buffer.""" """Tests (un)packing a subset of cell values of a field (from)to a buffer."""
num_cell_values = 7 num_cell_values = 7
# Cell indices of the field to be (un)packed (from)to the buffer # Cell indices of the field to be (un)packed (from)to the buffer
cell_indices = [1, 3, 5, 6] cell_indices = [1, 3, 5, 6]
...@@ -169,8 +181,9 @@ def test_subset_cell_values(): ...@@ -169,8 +181,9 @@ def test_subset_cell_values():
pack_eqs.append(eq) pack_eqs.append(eq)
pack_types = {'src_field': gpu_src_arr.dtype, 'buffer': gpu_buffer_arr.dtype} pack_types = {'src_field': gpu_src_arr.dtype, 'buffer': gpu_buffer_arr.dtype}
pack_code = create_cuda_kernel(pack_eqs, type_info=pack_types) config = CreateKernelConfig(target=pystencils.Target.GPU, data_type=pack_types)
pack_kernel = make_python_function(pack_code) pack_ast = create_kernel(pack_eqs, config=config)
pack_kernel = pack_ast.compile()
pack_kernel(buffer=gpu_buffer_arr, src_field=gpu_src_arr) pack_kernel(buffer=gpu_buffer_arr, src_field=gpu_src_arr)
unpack_eqs = [] unpack_eqs = []
...@@ -180,8 +193,10 @@ def test_subset_cell_values(): ...@@ -180,8 +193,10 @@ def test_subset_cell_values():
unpack_eqs.append(eq) unpack_eqs.append(eq)
unpack_types = {'dst_field': gpu_dst_arr.dtype, 'buffer': gpu_buffer_arr.dtype} unpack_types = {'dst_field': gpu_dst_arr.dtype, 'buffer': gpu_buffer_arr.dtype}
unpack_code = create_cuda_kernel(unpack_eqs, type_info=unpack_types) config = CreateKernelConfig(target=pystencils.Target.GPU, data_type=unpack_types)
unpack_kernel = make_python_function(unpack_code) unpack_ast = create_kernel(unpack_eqs, config=config)
unpack_kernel = unpack_ast.compile()
unpack_kernel(buffer=gpu_buffer_arr, dst_field=gpu_dst_arr) unpack_kernel(buffer=gpu_buffer_arr, dst_field=gpu_dst_arr)
dst_arr = gpu_dst_arr.get() dst_arr = gpu_dst_arr.get()
...@@ -190,7 +205,6 @@ def test_subset_cell_values(): ...@@ -190,7 +205,6 @@ def test_subset_cell_values():
np.testing.assert_equal(dst_arr, mask_arr.filled(int(0))) np.testing.assert_equal(dst_arr, mask_arr.filled(int(0)))
@pytest.mark.gpu
def test_field_layouts(): def test_field_layouts():
num_cell_values = 7 num_cell_values = 7
for layout_str in ['numpy', 'fzyx', 'zyxf', 'reverse_numpy']: for layout_str in ['numpy', 'fzyx', 'zyxf', 'reverse_numpy']:
...@@ -209,8 +223,10 @@ def test_field_layouts(): ...@@ -209,8 +223,10 @@ def test_field_layouts():
pack_eqs.append(eq) pack_eqs.append(eq)
pack_types = {'src_field': gpu_src_arr.dtype, 'buffer': gpu_buffer_arr.dtype} pack_types = {'src_field': gpu_src_arr.dtype, 'buffer': gpu_buffer_arr.dtype}
pack_code = create_cuda_kernel(pack_eqs, type_info=pack_types) config = CreateKernelConfig(target=pystencils.Target.GPU, data_type=pack_types)
pack_kernel = make_python_function(pack_code) pack_ast = create_kernel(pack_eqs, config=config)
pack_kernel = pack_ast.compile()
pack_kernel(buffer=gpu_buffer_arr, src_field=gpu_src_arr) pack_kernel(buffer=gpu_buffer_arr, src_field=gpu_src_arr)
unpack_eqs = [] unpack_eqs = []
...@@ -220,6 +236,99 @@ def test_field_layouts(): ...@@ -220,6 +236,99 @@ def test_field_layouts():
unpack_eqs.append(eq) unpack_eqs.append(eq)
unpack_types = {'dst_field': gpu_dst_arr.dtype, 'buffer': gpu_buffer_arr.dtype} unpack_types = {'dst_field': gpu_dst_arr.dtype, 'buffer': gpu_buffer_arr.dtype}
unpack_code = create_cuda_kernel(unpack_eqs, type_info=unpack_types) config = CreateKernelConfig(target=pystencils.Target.GPU, data_type=unpack_types)
unpack_kernel = make_python_function(unpack_code) unpack_ast = create_kernel(unpack_eqs, config=config)
unpack_kernel = unpack_ast.compile()
unpack_kernel(buffer=gpu_buffer_arr, dst_field=gpu_dst_arr) unpack_kernel(buffer=gpu_buffer_arr, dst_field=gpu_dst_arr)
def test_buffer_indexing():
src_field, dst_field = fields(f'pdfs_src(19), pdfs_dst(19) :double[3D]')
mask_field = fields(f'mask : uint32 [3D]')
buffer = Field.create_generic('buffer', spatial_dimensions=1, field_type=FieldType.BUFFER,
dtype="float64",
index_shape=(19,))
src_field_size = src_field.spatial_shape
mask_field_size = mask_field.spatial_shape
up = Assignment(buffer(0), flag_cond(1, mask_field.center, src_field[0, 1, 0](1)))
iteration_slice = tuple(slice(None, None, 2) for _ in range(3))
config = CreateKernelConfig(target=Target.GPU)
config = replace(config, iteration_slice=iteration_slice, ghost_layers=0)
ast = create_kernel(up, config=config)
parameters = ast.get_parameters()
spatial_shape_symbols = [p.symbol for p in parameters if p.is_field_shape]
# The loop counters as well as the resolved field access should depend on one common spatial shape
if spatial_shape_symbols[0] in mask_field_size:
for s in spatial_shape_symbols:
assert s in mask_field_size
if spatial_shape_symbols[0] in src_field_size:
for s in spatial_shape_symbols:
assert s in src_field_size
assert len(spatial_shape_symbols) <= 3
@pytest.mark.parametrize('gpu_indexing', ("block", "line"))
def test_iteration_slices(gpu_indexing):
num_cell_values = 19
dt = np.uint64
fields = _generate_fields(dt=dt, stencil_directions=num_cell_values)
for (src_arr, gpu_src_arr, gpu_dst_arr, gpu_buffer_arr) in fields:
src_field = Field.create_from_numpy_array("src_field", gpu_src_arr, index_dimensions=1)
dst_field = Field.create_from_numpy_array("dst_field", gpu_src_arr, index_dimensions=1)
buffer = Field.create_generic("buffer", spatial_dimensions=1, index_dimensions=1,
field_type=FieldType.BUFFER, dtype=src_arr.dtype)
pack_eqs = []
# Since we are packing all cell values for all cells, then
# the buffer index is equivalent to the field index
for idx in range(num_cell_values):
eq = Assignment(buffer(idx), src_field(idx))
pack_eqs.append(eq)
dim = src_field.spatial_dimensions
# Pack only the leftmost slice, only every second cell
pack_slice = (slice(None, None, 2),) * (dim - 1) + (0,)
# Fill the entire array with data
src_arr[(slice(None, None, 1),) * dim] = np.arange(num_cell_values)
gpu_src_arr.set(src_arr)
gpu_dst_arr.fill(0)
config = CreateKernelConfig(target=Target.GPU, iteration_slice=pack_slice,
data_type={'src_field': gpu_src_arr.dtype, 'buffer': gpu_buffer_arr.dtype},
gpu_indexing=gpu_indexing)
pack_code = create_kernel(pack_eqs, config=config)
pack_kernel = pack_code.compile()
pack_kernel(buffer=gpu_buffer_arr, src_field=gpu_src_arr)
unpack_eqs = []
for idx in range(num_cell_values):
eq = Assignment(dst_field(idx), buffer(idx))
unpack_eqs.append(eq)
config = CreateKernelConfig(target=Target.GPU, iteration_slice=pack_slice,
data_type={'dst_field': gpu_dst_arr.dtype, 'buffer': gpu_buffer_arr.dtype},
gpu_indexing=gpu_indexing)
unpack_code = create_kernel(unpack_eqs, config=config)
unpack_kernel = unpack_code.compile()
unpack_kernel(buffer=gpu_buffer_arr, dst_field=gpu_dst_arr)
dst_arr = gpu_dst_arr.get()
src_arr = gpu_src_arr.get()
# Check if only every second entry of the leftmost slice has been copied
np.testing.assert_equal(dst_arr[pack_slice], src_arr[pack_slice])
np.testing.assert_equal(dst_arr[(slice(1, None, 2),) * (dim - 1) + (0,)], 0)
np.testing.assert_equal(dst_arr[(slice(None, None, 1),) * (dim - 1) + (slice(1, None),)], 0)
# -*- coding: utf-8 -*-
#
# Copyright © 2019 Stephan Seitz <stephan.seitz@fau.de>
#
# Distributed under terms of the GPLv3 license.
"""
"""
import itertools
import numpy as np
import pytest
import sympy as sp
import pystencils as ps
from pystencils import Field, x_vector
from pystencils.astnodes import ConditionalFieldAccess
from pystencils.simp import sympy_cse
def add_fixed_constant_boundary_handling(assignments, with_cse):
common_shape = next(iter(set().union(itertools.chain.from_iterable(
[a.atoms(Field.Access) for a in assignments]
)))).field.spatial_shape
ndim = len(common_shape)
def is_out_of_bound(access, shape):
return sp.Or(*[sp.Or(a < 0, a >= s) for a, s in zip(access, shape)])
safe_assignments = [ps.Assignment(
assignment.lhs, assignment.rhs.subs({
a: ConditionalFieldAccess(a, is_out_of_bound(sp.Matrix(a.offsets) + x_vector(ndim), common_shape))
for a in assignment.rhs.atoms(Field.Access) if not a.is_absolute_access
})) for assignment in assignments.all_assignments]
# subs = [{a: ConditionalFieldAccess(a, is_out_of_bound(
# sp.Matrix(a.offsets) + x_vector(ndim), common_shape))
# for a in assignment.rhs.atoms(Field.Access) if not a.is_absolute_access
# } for assignment in assignments.all_assignments]
# print(subs)
if with_cse:
safe_assignments = sympy_cse(ps.AssignmentCollection(safe_assignments))
return safe_assignments
else:
return ps.AssignmentCollection(safe_assignments)
@pytest.mark.parametrize('dtype', ('float64', 'float32'))
@pytest.mark.parametrize('with_cse', (False, 'with_cse'))
def test_boundary_check(dtype, with_cse):
f, g = ps.fields(f"f, g : {dtype}[2D]")
stencil = ps.Assignment(g[0, 0], (f[1, 0] + f[-1, 0] + f[0, 1] + f[0, -1]) / 4)
f_arr = np.random.rand(10, 10).astype(dtype=dtype)
g_arr = np.zeros_like(f_arr)
assignments = add_fixed_constant_boundary_handling(ps.AssignmentCollection([stencil]), with_cse)
config = ps.CreateKernelConfig(data_type=dtype, default_number_float=dtype, ghost_layers=0)
kernel_checked = ps.create_kernel(assignments, config=config).compile()
# ps.show_code(kernel_checked)
# No SEGFAULT, please!!
kernel_checked(f=f_arr, g=g_arr)
import numpy as np
import sympy as sp
import pytest
import pystencils as ps
from pystencils.alignedarray import aligned_zeros
from pystencils.astnodes import Block, Conditional, SympyAssignment
from pystencils.backends.simd_instruction_sets import get_supported_instruction_sets, get_vector_instruction_set
from pystencils.enums import Target
from pystencils.cpu.vectorization import vec_all, vec_any
from pystencils.node_collection import NodeCollection
supported_instruction_sets = get_supported_instruction_sets() if get_supported_instruction_sets() else []
@pytest.mark.parametrize('instruction_set', supported_instruction_sets)
@pytest.mark.parametrize('dtype', ('float32', 'float64'))
def test_vec_any(instruction_set, dtype):
if instruction_set in ['sve', 'sve2', 'sme', 'rvv']:
width = 4 # we don't know the actual value
else:
width = get_vector_instruction_set(dtype, instruction_set)['width']
data_arr = np.zeros((4 * width, 4 * width), dtype=dtype)
data_arr[3:9, 1:3 * width - 1] = 1.0
data = ps.fields(f"data: {dtype}[2D]", data=data_arr)
c = [
SympyAssignment(sp.Symbol("t1"), vec_any(data.center() > 0.0)),
Conditional(vec_any(data.center() > 0.0), Block([SympyAssignment(data.center(), 2.0)]))
]
assignmets = NodeCollection(c)
ast = ps.create_kernel(assignments=assignmets, target=ps.Target.CPU,
cpu_vectorize_info={'instruction_set': instruction_set})
kernel = ast.compile()
kernel(data=data_arr)
if instruction_set in ['sve', 'sve2', 'sme', 'rvv']:
# we only know that the first value has changed
np.testing.assert_equal(data_arr[3:9, :3 * width - 1], 2.0)
else:
np.testing.assert_equal(data_arr[3:9, :3 * width], 2.0)
@pytest.mark.parametrize('instruction_set', supported_instruction_sets)
@pytest.mark.parametrize('dtype', ('float32', 'float64'))
def test_vec_all(instruction_set, dtype):
if instruction_set in ['sve', 'sve2', 'sme', 'rvv']:
width = 1000 # we don't know the actual value, need something guaranteed larger than vector
else:
width = get_vector_instruction_set(dtype, instruction_set)['width']
data_arr = np.zeros((4 * width, 4 * width), dtype=dtype)
data_arr[3:9, 1:3 * width - 1] = 1.0
data = ps.fields(f"data: {dtype}[2D]", data=data_arr)
c = [Conditional(vec_all(data.center() > 0.0), Block([SympyAssignment(data.center(), 2.0)]))]
assignmets = NodeCollection(c)
ast = ps.create_kernel(assignmets, target=Target.CPU,
cpu_vectorize_info={'instruction_set': instruction_set})
kernel = ast.compile()
kernel(data=data_arr)
if instruction_set in ['sve', 'sve2', 'sme', 'rvv']:
# we only know that some values in the middle have been replaced
assert np.all(data_arr[3:9, :2] <= 1.0)
assert np.any(data_arr[3:9, 2:] == 2.0)
else:
np.testing.assert_equal(data_arr[3:9, :1], 0.0)
np.testing.assert_equal(data_arr[3:9, 1:width], 1.0)
np.testing.assert_equal(data_arr[3:9, width:2 * width], 2.0)
np.testing.assert_equal(data_arr[3:9, 2 * width:3 * width - 1], 1.0)
np.testing.assert_equal(data_arr[3:9, 3 * width - 1:], 0.0)
@pytest.mark.skipif(not supported_instruction_sets, reason='cannot detect CPU instruction set')
def test_boolean_before_loop():
t1, t2 = sp.symbols('t1, t2')
f_arr = np.ones((10, 10))
g_arr = np.zeros_like(f_arr)
f, g = ps.fields("f, g : double[2D]", f=f_arr, g=g_arr)
a = [
ps.Assignment(t1, t2 > 0),
ps.Assignment(g[0, 0],
sp.Piecewise((f[0, 0], t1), (42, True)))
]
ast = ps.create_kernel(a, cpu_vectorize_info={'instruction_set': supported_instruction_sets[-1]})
kernel = ast.compile()
kernel(f=f_arr, g=g_arr, t2=1.0)
# print(g)
np.testing.assert_array_equal(g_arr, 1.0)
kernel(f=f_arr, g=g_arr, t2=-1.0)
np.testing.assert_array_equal(g_arr, 42.0)
@pytest.mark.parametrize('instruction_set', supported_instruction_sets)
@pytest.mark.parametrize('dtype', ('float32', 'float64'))
@pytest.mark.parametrize('nontemporal', [False, True])
@pytest.mark.parametrize('aligned', [False, True])
def test_vec_maskstore(instruction_set, dtype, nontemporal, aligned):
data_arr = (aligned_zeros if aligned else np.zeros)((16, 16), dtype=dtype)
data_arr[3:-3, 3:-3] = 1.0
data = ps.fields(f"data: {dtype}[2D]", data=data_arr)
c = [Conditional(data.center() < 1.0, Block([SympyAssignment(data.center(), 2.0)]))]
assignmets = NodeCollection(c)
config = ps.CreateKernelConfig(cpu_vectorize_info={'instruction_set': instruction_set,
'nontemporal': nontemporal,
'assume_aligned': aligned},
default_number_float=dtype)
ast = ps.create_kernel(assignmets, config=config)
if 'maskStore' in ast.instruction_set:
instruction = 'maskStream' if nontemporal and 'maskStream' in ast.instruction_set else (
'maskStoreA' if aligned and 'maskStoreA' in ast.instruction_set else 'maskStore')
assert ast.instruction_set[instruction].split('{')[0] in ps.get_code_str(ast)
print(ps.get_code_str(ast))
kernel = ast.compile()
kernel(data=data_arr)
np.testing.assert_equal(data_arr[:3, :], 2.0)
np.testing.assert_equal(data_arr[-3:, :], 2.0)
np.testing.assert_equal(data_arr[:, :3], 2.0)
np.testing.assert_equal(data_arr[:, -3:], 2.0)
np.testing.assert_equal(data_arr[3:-3, 3:-3], 1.0)
@pytest.mark.parametrize('instruction_set', supported_instruction_sets)
@pytest.mark.parametrize('dtype', ('float32', 'float64'))
@pytest.mark.parametrize('nontemporal', [False, True])
def test_vec_maskscatter(instruction_set, dtype, nontemporal):
data_arr = np.zeros((16, 16), dtype=dtype)
data_arr[3:-3, 3:-3] = 1.0
data = ps.fields(f"data: {dtype}[2D]")
c = [Conditional(data.center() < 1.0, Block([SympyAssignment(data.center(), 2.0)]))]
assignmets = NodeCollection(c)
config = ps.CreateKernelConfig(cpu_vectorize_info={'instruction_set': instruction_set,
'nontemporal': nontemporal},
default_number_float=dtype)
if 'maskStoreS' not in get_vector_instruction_set(dtype, instruction_set) \
and not instruction_set.startswith('sve'):
with pytest.warns(UserWarning) as warn:
ast = ps.create_kernel(assignmets, config=config)
assert 'Could not vectorize loop' in warn[0].message.args[0]
else:
with pytest.warns(None) as warn:
ast = ps.create_kernel(assignmets, config=config)
assert len(warn) == 0
instruction = 'maskStreamS' if nontemporal and 'maskStreamS' in ast.instruction_set else 'maskStoreS'
assert ast.instruction_set[instruction].split('{')[0] in ps.get_code_str(ast)
print(ps.get_code_str(ast))
kernel = ast.compile()
kernel(data=data_arr)
np.testing.assert_equal(data_arr[:3, :], 2.0)
np.testing.assert_equal(data_arr[-3:, :], 2.0)
np.testing.assert_equal(data_arr[:, :3], 2.0)
np.testing.assert_equal(data_arr[:, -3:], 2.0)
np.testing.assert_equal(data_arr[3:-3, 3:-3], 1.0)
from collections import defaultdict
import numpy as np
import pytest
from pystencils import CreateKernelConfig, Target, Backend
from pystencils.typing import BasicType
def test_config():
# targets
config = CreateKernelConfig(target=Target.CPU)
assert config.target == Target.CPU
assert config.backend == Backend.C
config = CreateKernelConfig(target=Target.GPU)
assert config.target == Target.GPU
assert config.backend == Backend.CUDA
# typing
config = CreateKernelConfig(data_type=np.float64)
assert isinstance(config.data_type, defaultdict)
assert config.data_type.default_factory() == BasicType('float64')
assert config.default_number_float == BasicType('float64')
assert config.default_number_int == BasicType('int64')
config = CreateKernelConfig(data_type=np.float32)
assert isinstance(config.data_type, defaultdict)
assert config.data_type.default_factory() == BasicType('float32')
assert config.default_number_float == BasicType('float32')
assert config.default_number_int == BasicType('int64')
config = CreateKernelConfig(data_type=np.float32, default_number_float=np.float64)
assert isinstance(config.data_type, defaultdict)
assert config.data_type.default_factory() == BasicType('float32')
assert config.default_number_float == BasicType('float64')
assert config.default_number_int == BasicType('int64')
config = CreateKernelConfig(data_type=np.float32, default_number_float=np.float64, default_number_int=np.int16)
assert isinstance(config.data_type, defaultdict)
assert config.data_type.default_factory() == BasicType('float32')
assert config.default_number_float == BasicType('float64')
assert config.default_number_int == BasicType('int16')
config = CreateKernelConfig(data_type='float64')
assert isinstance(config.data_type, defaultdict)
assert config.data_type.default_factory() == BasicType('float64')
assert config.default_number_float == BasicType('float64')
assert config.default_number_int == BasicType('int64')
config = CreateKernelConfig(data_type={'a': np.float64, 'b': np.float32})
assert isinstance(config.data_type, defaultdict)
assert config.data_type.default_factory() == BasicType('float64')
assert config.default_number_float == BasicType('float64')
assert config.default_number_int == BasicType('int64')
config = CreateKernelConfig(data_type={'a': np.float32, 'b': np.int32})
assert isinstance(config.data_type, defaultdict)
assert config.data_type.default_factory() == BasicType('float32')
assert config.default_number_float == BasicType('float32')
assert config.default_number_int == BasicType('int64')
def test_config_target_as_string():
with pytest.raises(ValueError):
CreateKernelConfig(target='cpu')
def test_config_backend_as_string():
with pytest.raises(ValueError):
CreateKernelConfig(backend='C')
def test_config_python_types():
with pytest.raises(ValueError):
CreateKernelConfig(data_type=float)
def test_config_python_types2():
with pytest.raises(ValueError):
CreateKernelConfig(data_type={'a': float})
def test_config_python_types3():
with pytest.raises(ValueError):
CreateKernelConfig(default_number_float=float)
def test_config_python_types4():
with pytest.raises(ValueError):
CreateKernelConfig(default_number_int=int)
def test_config_python_types5():
with pytest.raises(ValueError):
CreateKernelConfig(data_type="float")
def test_config_python_types6():
with pytest.raises(ValueError):
CreateKernelConfig(default_number_float="float")
def test_config_python_types7():
dtype = defaultdict(lambda: 'float', {'a': np.float64, 'b': np.int64})
with pytest.raises(ValueError):
CreateKernelConfig(data_type=dtype)
def test_config_python_types8():
dtype = defaultdict(lambda: float, {'a': np.float64, 'b': np.int64})
with pytest.raises(ValueError):
CreateKernelConfig(data_type=dtype)
def test_config_python_types9():
dtype = defaultdict(lambda: 'float32', {'a': 'float', 'b': np.int64})
with pytest.raises(ValueError):
CreateKernelConfig(data_type=dtype)
def test_config_python_types10():
dtype = defaultdict(lambda: 'float32', {'a': float, 'b': np.int64})
with pytest.raises(ValueError):
CreateKernelConfig(data_type=dtype)
import numpy as np
import sympy as sp
import pystencils as ps
import pystencils.config
def test_create_kernel_config():
c = pystencils.config.CreateKernelConfig()
assert c.backend == ps.Backend.C
assert c.target == ps.Target.CPU
c = pystencils.config.CreateKernelConfig(target=ps.Target.GPU)
assert c.backend == ps.Backend.CUDA
c = pystencils.config.CreateKernelConfig(backend=ps.Backend.CUDA)
assert c.target == ps.Target.CPU
assert c.backend == ps.Backend.CUDA
def test_kernel_decorator_config():
config = pystencils.config.CreateKernelConfig()
a, b, c = ps.fields(a=np.ones(100), b=np.ones(100), c=np.ones(100))
@ps.kernel_config(config)
def test():
a[0] @= b[0] + c[0]
ps.create_kernel(**test)
def test_kernel_decorator2():
h = sp.symbols("h")
dtype = "float64"
src, dst = ps.fields(f"src, src_tmp: {dtype}[3D]")
@ps.kernel
def kernel_func():
dst[0, 0, 0] @= (src[1, 0, 0] + src[-1, 0, 0]
+ src[0, 1, 0] + src[0, -1, 0]
+ src[0, 0, 1] + src[0, 0, -1]) / (6 * h ** 2)
# assignments = ps.assignment_from_stencil(stencil, src, dst, normalization_factor=2)
ast = ps.create_kernel(kernel_func)
code = ps.get_code_str(ast)
from subprocess import CalledProcessError from subprocess import CalledProcessError
import pycuda.driver
import pytest import pytest
import sympy
import pystencils import pystencils
import pystencils.cpu.cpujit import pystencils.cpu.cpujit
import pystencils.gpucuda.cudajit
from pystencils.backends.cbackend import CBackend from pystencils.backends.cbackend import CBackend
from pystencils.backends.cuda_backend import CudaBackend from pystencils.backends.cuda_backend import CudaBackend
from pystencils.enums import Target
class ScreamingBackend(CBackend): class ScreamingBackend(CBackend):
...@@ -25,18 +23,29 @@ class ScreamingGpuBackend(CudaBackend): ...@@ -25,18 +23,29 @@ class ScreamingGpuBackend(CudaBackend):
return normal_code.upper() return normal_code.upper()
def test_custom_backends(): def test_custom_backends_cpu():
z, x, y = pystencils.fields("z, y, x: [2d]") z, y, x = pystencils.fields("z, y, x: [2d]")
normal_assignments = pystencils.AssignmentCollection([pystencils.Assignment( normal_assignments = pystencils.AssignmentCollection([pystencils.Assignment(
z[0, 0], x[0, 0] * sympy.log(x[0, 0] * y[0, 0]))], []) z[0, 0], x[0, 0] * x[0, 0] * y[0, 0])], [])
ast = pystencils.create_kernel(normal_assignments, target='cpu') ast = pystencils.create_kernel(normal_assignments, target=Target.CPU)
print(pystencils.show_code(ast, ScreamingBackend())) pystencils.show_code(ast, ScreamingBackend())
with pytest.raises(CalledProcessError): with pytest.raises(CalledProcessError):
pystencils.cpu.cpujit.make_python_function(ast, custom_backend=ScreamingBackend()) pystencils.cpu.cpujit.make_python_function(ast, custom_backend=ScreamingBackend())
ast = pystencils.create_kernel(normal_assignments, target='gpu')
print(pystencils.show_code(ast, ScreamingGpuBackend())) def test_custom_backends_gpu():
with pytest.raises(pycuda.driver.CompileError): pytest.importorskip('cupy')
pystencils.gpucuda.cudajit.make_python_function(ast, custom_backend=ScreamingGpuBackend()) import cupy
import pystencils.gpu.gpujit
z, x, y = pystencils.fields("z, y, x: [2d]")
normal_assignments = pystencils.AssignmentCollection([pystencils.Assignment(
z[0, 0], x[0, 0] * x[0, 0] * y[0, 0])], [])
ast = pystencils.create_kernel(normal_assignments, target=Target.GPU)
pystencils.show_code(ast, ScreamingGpuBackend())
with pytest.raises((cupy.cuda.compiler.JitifyException, cupy.cuda.compiler.CompileException)):
pystencils.gpu.gpujit.make_python_function(ast, custom_backend=ScreamingGpuBackend())
File added
File added
File added
File added
File added
File added
tests/test_data/lenna.png

463 KiB

tests/test_data/test_vessel2d_mask.png

7.43 KiB

import os import os
from tempfile import TemporaryDirectory from tempfile import TemporaryDirectory
from pathlib import Path
import numpy as np import numpy as np
import pystencils as ps import pystencils as ps
from pystencils import create_data_handling, create_kernel from pystencils import create_data_handling, create_kernel
from pystencils.gpu.gpu_array_handler import GPUArrayHandler
from pystencils.enums import Target
try:
import pytest
except ImportError:
import unittest.mock
pytest = unittest.mock.MagicMock()
try:
import cupy.cuda.runtime
device_numbers = range(cupy.cuda.runtime.getDeviceCount())
except ImportError:
device_numbers = []
SCRIPT_FOLDER = Path(__file__).parent.absolute()
INPUT_FOLDER = SCRIPT_FOLDER / "test_data"
def basic_iteration(dh): def basic_iteration(dh):
...@@ -17,7 +35,7 @@ def basic_iteration(dh): ...@@ -17,7 +35,7 @@ def basic_iteration(dh):
def access_and_gather(dh, domain_size): def access_and_gather(dh, domain_size):
dh.add_array('f1', dtype=np.dtype(np.int32)) dh.add_array('f1', dtype=np.dtype(np.int8))
dh.add_array_like('f2', 'f1') dh.add_array_like('f2', 'f1')
dh.add_array('v1', values_per_cell=3, dtype=np.int64, ghost_layers=2) dh.add_array('v1', values_per_cell=3, dtype=np.int64, ghost_layers=2)
dh.add_array_like('v2', 'v1') dh.add_array_like('v2', 'v1')
...@@ -28,7 +46,7 @@ def access_and_gather(dh, domain_size): ...@@ -28,7 +46,7 @@ def access_and_gather(dh, domain_size):
# Check symbolic field properties # Check symbolic field properties
assert dh.fields.f1.index_dimensions == 0 assert dh.fields.f1.index_dimensions == 0
assert dh.fields.f1.spatial_dimensions == len(domain_size) assert dh.fields.f1.spatial_dimensions == len(domain_size)
assert dh.fields.f1.dtype.numpy_dtype == np.int32 assert dh.fields.f1.dtype.numpy_dtype == np.int8
assert dh.fields.v1.index_dimensions == 1 assert dh.fields.v1.index_dimensions == 1
assert dh.fields.v1.spatial_dimensions == len(domain_size) assert dh.fields.v1.spatial_dimensions == len(domain_size)
...@@ -73,14 +91,10 @@ def access_and_gather(dh, domain_size): ...@@ -73,14 +91,10 @@ def access_and_gather(dh, domain_size):
def synchronization(dh, test_gpu=False): def synchronization(dh, test_gpu=False):
field_name = 'comm_field_test' field_name = 'comm_field_test'
if test_gpu: if test_gpu:
try: pytest.importorskip("cupy")
from pycuda import driver
import pycuda.autoinit
except ImportError:
return
field_name += 'Gpu' field_name += 'Gpu'
dh.add_array(field_name, ghost_layers=1, dtype=np.int32, cpu=True, gpu=test_gpu) dh.add_array(field_name, ghost_layers=1, dtype=np.int8, cpu=True, gpu=test_gpu)
# initialize everything with 1 # initialize everything with 1
for b in dh.iterate(ghost_layers=1): for b in dh.iterate(ghost_layers=1):
...@@ -90,8 +104,10 @@ def synchronization(dh, test_gpu=False): ...@@ -90,8 +104,10 @@ def synchronization(dh, test_gpu=False):
if test_gpu: if test_gpu:
dh.to_gpu(field_name) dh.to_gpu(field_name)
dh.synchronization_function_gpu(field_name)()
else:
dh.synchronization_function_cpu(field_name)()
dh.synchronization_function(field_name, target='gpu' if test_gpu else 'cpu')()
if test_gpu: if test_gpu:
dh.to_cpu(field_name) dh.to_cpu(field_name)
...@@ -100,17 +116,16 @@ def synchronization(dh, test_gpu=False): ...@@ -100,17 +116,16 @@ def synchronization(dh, test_gpu=False):
np.testing.assert_equal(42, b[field_name]) np.testing.assert_equal(42, b[field_name])
def kernel_execution_jacobi(dh, test_gpu=False): def kernel_execution_jacobi(dh, target):
if test_gpu:
try:
from pycuda import driver
import pycuda.autoinit
except ImportError:
print("Skipping kernel_execution_jacobi GPU version, because pycuda not available")
return
test_gpu = target == Target.GPU
dh.add_array('f', gpu=test_gpu) dh.add_array('f', gpu=test_gpu)
dh.add_array('tmp', gpu=test_gpu) dh.add_array('tmp', gpu=test_gpu)
if test_gpu:
assert dh.is_on_gpu('f')
assert dh.is_on_gpu('tmp')
stencil_2d = [(1, 0), (-1, 0), (0, 1), (0, -1)] stencil_2d = [(1, 0), (-1, 0), (0, 1), (0, -1)]
stencil_3d = [(1, 0, 0), (-1, 0, 0), (0, 1, 0), (0, -1, 0), (0, 0, 1), (0, 0, -1)] stencil_3d = [(1, 0, 0), (-1, 0, 0), (0, 1, 0), (0, -1, 0), (0, 0, 1), (0, 0, -1)]
stencil = stencil_2d if dh.dim == 2 else stencil_3d stencil = stencil_2d if dh.dim == 2 else stencil_3d
...@@ -119,7 +134,7 @@ def kernel_execution_jacobi(dh, test_gpu=False): ...@@ -119,7 +134,7 @@ def kernel_execution_jacobi(dh, test_gpu=False):
def jacobi(): def jacobi():
dh.fields.tmp.center @= sum(dh.fields.f.neighbors(stencil)) / len(stencil) dh.fields.tmp.center @= sum(dh.fields.f.neighbors(stencil)) / len(stencil)
kernel = create_kernel(jacobi, target='gpu' if test_gpu else 'cpu').compile() kernel = create_kernel(jacobi, config=ps.CreateKernelConfig(target=target)).compile()
for b in dh.iterate(ghost_layers=1): for b in dh.iterate(ghost_layers=1):
b['f'].fill(42) b['f'].fill(42)
dh.run_kernel(kernel) dh.run_kernel(kernel)
...@@ -128,6 +143,7 @@ def kernel_execution_jacobi(dh, test_gpu=False): ...@@ -128,6 +143,7 @@ def kernel_execution_jacobi(dh, test_gpu=False):
def vtk_output(dh): def vtk_output(dh):
pytest.importorskip('pyevtk')
dh.add_array('scalar_field') dh.add_array('scalar_field')
dh.add_array('vector_field', values_per_cell=dh.dim) dh.add_array('vector_field', values_per_cell=dh.dim)
dh.add_array('multiple_scalar_field', values_per_cell=9) dh.add_array('multiple_scalar_field', values_per_cell=9)
...@@ -196,18 +212,200 @@ def test_access_and_gather(): ...@@ -196,18 +212,200 @@ def test_access_and_gather():
def test_kernel(): def test_kernel():
for domain_shape in [(4, 5), (3, 4, 5)]: for domain_shape in [(4, 5), (3, 4, 5)]:
dh = create_data_handling(domain_size=domain_shape, periodicity=True) dh = create_data_handling(domain_size=domain_shape, periodicity=True)
kernel_execution_jacobi(dh, test_gpu=True) assert all(dh.periodicity)
kernel_execution_jacobi(dh, Target.CPU)
reduction(dh) reduction(dh)
try: try:
import pycuda import cupy
dh = create_data_handling(domain_size=domain_shape, periodicity=True) dh = create_data_handling(domain_size=domain_shape, periodicity=True)
kernel_execution_jacobi(dh, test_gpu=False) kernel_execution_jacobi(dh, Target.GPU)
except ImportError: except ImportError:
pass pass
@pytest.mark.parametrize('target', (Target.CPU, Target.GPU))
def test_kernel_param(target):
for domain_shape in [(4, 5), (3, 4, 5)]:
if target == Target.GPU:
pytest.importorskip('cupy')
dh = create_data_handling(domain_size=domain_shape, periodicity=True, default_target=target)
kernel_execution_jacobi(dh, target)
reduction(dh)
def test_vtk_output(): def test_vtk_output():
pytest.importorskip('pyevtk')
for domain_shape in [(4, 5), (3, 4, 5)]: for domain_shape in [(4, 5), (3, 4, 5)]:
dh = create_data_handling(domain_size=domain_shape, periodicity=True) dh = create_data_handling(domain_size=domain_shape, periodicity=True)
vtk_output(dh) vtk_output(dh)
def test_add_arrays():
domain_shape = (3, 4, 5)
field_description = 'x, y(9)'
dh = create_data_handling(domain_size=domain_shape, default_ghost_layers=0, default_layout='numpy')
x_, y_ = dh.add_arrays(field_description)
x, y = ps.fields(field_description + ': [3,4,5]')
assert x_ == x
assert y_ == y
assert x == dh.fields['x']
assert y == dh.fields['y']
@pytest.mark.parametrize('shape', [(17, 12), (7, 11, 18)])
@pytest.mark.parametrize('layout', ['zyxf', 'fzyx'])
def test_add_arrays_with_layout(shape, layout):
pytest.importorskip('cupy')
dh = create_data_handling(domain_size=shape, default_layout=layout, default_target=ps.Target.GPU)
f1 = dh.add_array("f1", values_per_cell=19)
dh.fill(f1.name, 1.0)
assert dh.cpu_arrays[f1.name].shape == dh.gpu_arrays[f1.name].shape
assert dh.cpu_arrays[f1.name].strides == dh.gpu_arrays[f1.name].strides
assert dh.cpu_arrays[f1.name].dtype == dh.gpu_arrays[f1.name].dtype
def test_get_kwarg():
domain_shape = (10, 10)
field_description = 'src, dst'
dh = create_data_handling(domain_size=domain_shape, default_ghost_layers=1)
src, dst = dh.add_arrays(field_description)
dh.fill("src", 1.0, ghost_layers=True)
dh.fill("dst", 0.0, ghost_layers=True)
with pytest.raises(ValueError):
dh.add_array('src', values_per_cell=1)
ur = ps.Assignment(src.center, dst.center)
kernel = ps.create_kernel(ur).compile()
kw = dh.get_kernel_kwargs(kernel)
assert np.all(kw[0]['src'] == dh.cpu_arrays['src'])
assert np.all(kw[0]['dst'] == dh.cpu_arrays['dst'])
def test_add_custom_data():
pytest.importorskip('cupy')
import cupy as cp
def cpu_data_create_func():
return np.ones((2, 2), dtype=np.float64)
def gpu_data_create_func():
return cp.zeros((2, 2), dtype=np.float64)
def cpu_to_gpu_transfer_func(gpuarr, cpuarray):
gpuarr.set(cpuarray)
def gpu_to_cpu_transfer_func(gpuarr, cpuarray):
cpuarray[:] = gpuarr.get()
dh = create_data_handling(domain_size=(10, 10))
dh.add_custom_data('custom_data',
cpu_data_create_func,
gpu_data_create_func,
cpu_to_gpu_transfer_func,
gpu_to_cpu_transfer_func)
assert np.all(dh.custom_data_cpu['custom_data'] == 1)
assert np.all(dh.custom_data_gpu['custom_data'].get() == 0)
dh.to_cpu(name='custom_data')
dh.to_gpu(name='custom_data')
assert 'custom_data' in dh.custom_data_names
def test_log():
dh = create_data_handling(domain_size=(10, 10))
dh.log_on_root()
assert dh.is_root
assert dh.world_rank == 0
def test_save_data():
domain_shape = (2, 2)
dh = create_data_handling(domain_size=domain_shape, default_ghost_layers=1)
dh.add_array("src", values_per_cell=9)
dh.fill("src", 1.0, ghost_layers=True)
dh.add_array("dst", values_per_cell=9)
dh.fill("dst", 1.0, ghost_layers=True)
dh.save_all(str(INPUT_FOLDER) + '/datahandling_save_test')
def test_load_data():
domain_shape = (2, 2)
dh = create_data_handling(domain_size=domain_shape, default_ghost_layers=1)
dh.add_array("src", values_per_cell=9)
dh.fill("src", 0.0, ghost_layers=True)
dh.add_array("dst", values_per_cell=9)
dh.fill("dst", 0.0, ghost_layers=True)
dh.load_all(str(INPUT_FOLDER) + '/datahandling_load_test')
assert np.all(dh.cpu_arrays['src']) == 1
assert np.all(dh.cpu_arrays['dst']) == 1
domain_shape = (3, 3)
dh = create_data_handling(domain_size=domain_shape, default_ghost_layers=1)
dh.add_array("src", values_per_cell=9)
dh.fill("src", 0.0, ghost_layers=True)
dh.add_array("dst", values_per_cell=9)
dh.fill("dst", 0.0, ghost_layers=True)
dh.add_array("dst2", values_per_cell=9)
dh.fill("dst2", 0.0, ghost_layers=True)
dh.load_all(str(INPUT_FOLDER) + '/datahandling_load_test')
assert np.all(dh.cpu_arrays['src']) == 0
assert np.all(dh.cpu_arrays['dst']) == 0
assert np.all(dh.cpu_arrays['dst2']) == 0
@pytest.mark.parametrize("device_number", device_numbers)
def test_array_handler(device_number):
size = (2, 2)
pytest.importorskip('cupy')
array_handler = GPUArrayHandler(device_number)
zero_array = array_handler.zeros(size)
cpu_array = np.empty(size)
array_handler.download(zero_array, cpu_array)
assert np.all(cpu_array) == 0
ones_array = array_handler.ones(size)
cpu_array = np.empty(size)
array_handler.download(ones_array, cpu_array)
assert np.all(cpu_array) == 1
empty = array_handler.empty(size)
assert empty.strides == (16, 8)
empty = array_handler.empty(shape=size, order="F")
assert empty.strides == (8, 16)
random_array = array_handler.randn(size)
cpu_array = np.empty((20, 40), dtype=np.float64)
gpu_array = array_handler.to_gpu(cpu_array)
assert cpu_array.base is None
assert gpu_array.base is None
assert gpu_array.strides == cpu_array.strides
cpu_array2 = np.empty((20, 40), dtype=np.float64)
cpu_array2 = cpu_array2.swapaxes(0, 1)
gpu_array2 = array_handler.to_gpu(cpu_array2)
assert cpu_array2.base is not None
assert gpu_array2.base is not None
assert gpu_array2.strides == cpu_array2.strides