test_blocking.py 2.27 KB
Newer Older
1
import numpy as np
Martin Bauer's avatar
Martin Bauer committed
2
3
import sympy as sp

4
5
6
7
8
9
10
11
12
13
14
15
16
17
import pystencils as ps


def jacobi(dst, src):
    assert dst.spatial_dimensions == src.spatial_dimensions
    assert src.index_dimensions == 0 and dst.index_dimensions == 0
    neighbors = []
    for d in range(src.spatial_dimensions):
        neighbors += [src.neighbor(d, offset) for offset in (1, -1)]
    return ps.Assignment(dst.center, sp.Add(*neighbors) / len(neighbors))


def check_equivalence(assignments, src_arr):
    for openmp in (False, True):
18
        for vectorization in [False, {'assume_inner_stride_one': True}]:
19
20
21
            with_blocking = ps.create_kernel(assignments, cpu_blocking=(8, 16, 4), cpu_openmp=openmp,
                                             cpu_vectorize_info=vectorization).compile()
            without_blocking = ps.create_kernel(assignments).compile()
22
            print(f"  openmp {openmp}, vectorization {vectorization}")
23
24
25
26
27
28
29
30
31
            dst_arr = np.zeros_like(src_arr)
            ref_arr = np.zeros_like(src_arr)
            np.copyto(src_arr, np.random.rand(*src_arr.shape))
            with_blocking(src=src_arr, dst=dst_arr)
            without_blocking(src=src_arr, dst=ref_arr)
            np.testing.assert_almost_equal(ref_arr, dst_arr)


def test_jacobi3d_var_size():
32
    src, dst = ps.fields("src, dst: double[3D]", layout='c')
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61

    print("Var Size: Smaller than block sizes")
    arr = np.empty([4, 5, 6])
    check_equivalence(jacobi(dst, src), arr)

    print("Var Size: Large non divisible sizes")
    arr = np.empty([100, 80, 9])
    check_equivalence(jacobi(dst, src), arr)

    print("Var Size: Multiples of block sizes")
    arr = np.empty([8*4, 16*2, 4*3])
    check_equivalence(jacobi(dst, src), arr)


def test_jacobi3d_fixed_size():
    print("Fixed Size: Large non divisible sizes")
    arr = np.empty([10, 10, 9])
    src, dst = ps.fields("src, dst: double[3D]", src=arr, dst=arr)
    check_equivalence(jacobi(dst, src), arr)

    print("Fixed Size: Smaller than block sizes")
    arr = np.empty([4, 5, 6])
    src, dst = ps.fields("src, dst: double[3D]", src=arr, dst=arr)
    check_equivalence(jacobi(dst, src), arr)

    print("Fixed Size: Multiples of block sizes")
    arr = np.empty([8*4, 16*2, 4*3])
    src, dst = ps.fields("src, dst: double[3D]", src=arr, dst=arr)
    check_equivalence(jacobi(dst, src), arr)