diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml index 59877cae1d7517d041ab799c064d89b3c7d3ccb5..d44dd00060d68ca69050f0a9f804c3b49e5d4c65 100644 --- a/.gitlab-ci.yml +++ b/.gitlab-ci.yml @@ -33,15 +33,15 @@ Python 2.7 POCL: reports: junit: test/pytest.xml -Python 3.6 POCL: +Python 3 POCL: script: - - export PY_EXE=python3.6 + - export PY_EXE=python3 - export PYOPENCL_TEST=portable - export EXTRA_INSTALL="pybind11 numpy mako" - curl -L -O -k https://gitlab.tiker.net/inducer/ci-support/raw/master/build-and-test-py-project.sh - ". ./build-and-test-py-project.sh" tags: - - python3.6 + - python3 - pocl except: - tags @@ -49,31 +49,15 @@ Python 3.6 POCL: reports: junit: test/pytest.xml -Python 3.7 POCL: +Python 3 Titan X: script: - - export PY_EXE=python3.7 - - export PYOPENCL_TEST=portable - - export EXTRA_INSTALL="pybind11 numpy mako" - - curl -L -O -k https://gitlab.tiker.net/inducer/ci-support/raw/master/build-and-test-py-project.sh - - ". ./build-and-test-py-project.sh" - tags: - - python3.7 - - pocl - except: - - tags - artifacts: - reports: - junit: test/pytest.xml - -Python 3.7 Titan X: - script: - - py_version=3.7 + - py_version=3 - export PYOPENCL_TEST=nvi:titan - EXTRA_INSTALL="pybind11 numpy mako" - curl -L -O -k https://gitlab.tiker.net/inducer/ci-support/raw/master/build-and-test-py-project.sh - ". ./build-and-test-py-project.sh" tags: - - python3.7 + - python3 - nvidia-titan-x except: - tags @@ -105,7 +89,7 @@ Documentation: - curl -L -O -k https://gitlab.tiker.net/inducer/ci-support/raw/master/build-docs.sh - ". ./build-docs.sh" tags: - - python3.5 + - python3 only: - master @@ -114,7 +98,7 @@ Flake8: - curl -L -O -k https://gitlab.tiker.net/inducer/ci-support/raw/master/prepare-and-run-flake8.sh - ". ./prepare-and-run-flake8.sh sumpy test" tags: - - python3.5 + - python3 except: - tags diff --git a/sumpy/visualization.py b/sumpy/visualization.py index 3a4dc8df66512fdf1da8e5e3c40fa1794e0bcfad..bf5f016b4c089be7723a2ff92c91530c998603c2 100644 --- a/sumpy/visualization.py +++ b/sumpy/visualization.py @@ -64,35 +64,35 @@ def separate_by_real_and_imag(data, real_only): def make_field_plotter_from_bbox(bbox, h, extend_factor=0): - """ - :arg bbox: a tuple (low, high) of points represented as 1D numpy arrays - indicating the low and high ends of the extent of a bounding box. - :arg h: Either a number or a sequence of numbers indicating the desired - (approximate) grid spacing in all or each of the dimensions. If a - sequence, the length must match the number of dimensions. - :arg extend_factor: A floating point number indicating by what percentage - the plot area should be grown compared to *bbox*. - """ - low, high = bbox - - extent = (high-low) * (1 + extend_factor) - center = 0.5*(high+low) - - dimensions = len(center) - from numbers import Number - if isinstance(h, Number): - h = (h,)*dimensions - else: - if len(h) != dimensions: - raise ValueError("length of 'h' must match number of dimensions") + """ + :arg bbox: a tuple (low, high) of points represented as 1D numpy arrays + indicating the low and high ends of the extent of a bounding box. + :arg h: Either a number or a sequence of numbers indicating the desired + (approximate) grid spacing in all or each of the dimensions. If a + sequence, the length must match the number of dimensions. + :arg extend_factor: A floating point number indicating by what percentage + the plot area should be grown compared to *bbox*. + """ + low, high = bbox + + extent = (high-low) * (1 + extend_factor) + center = 0.5*(high+low) + + dimensions = len(center) + from numbers import Number + if isinstance(h, Number): + h = (h,)*dimensions + else: + if len(h) != dimensions: + raise ValueError("length of 'h' must match number of dimensions") - from math import ceil + from math import ceil - npoints = tuple( - int(ceil(extent[i] / h[i])) - for i in range(dimensions)) + npoints = tuple( + int(ceil(extent[i] / h[i])) + for i in range(dimensions)) - return FieldPlotter(center, extent, npoints) + return FieldPlotter(center, extent, npoints) class FieldPlotter(object): diff --git a/test/test_misc.py b/test/test_misc.py index 7b35b2ddc3743d0eb9b7e80daf38c617ea73940d..33e8535a8cbce064263958962421d9b43e77ee6c 100644 --- a/test/test_misc.py +++ b/test/test_misc.py @@ -32,6 +32,8 @@ import pyopencl as cl # noqa: F401 from pyopencl.tools import ( # noqa pytest_generate_tests_for_pyopencl as pytest_generate_tests) +from pytools import Record + from sumpy.kernel import (LaplaceKernel, HelmholtzKernel, BiharmonicKernel, YukawaKernel) @@ -135,6 +137,122 @@ def test_order_finder(knl): print(orders) +# {{{ expansion toys p2e2e2p test cases + +def approx_convergence_factor(orders, errors): + poly = np.polyfit(orders, np.log(errors), deg=1) + return np.exp(poly[0]) + + +class P2E2E2PTestCase(Record): + + @property + def dim(self): + return len(self.source) + + @staticmethod + def eval(expr, source, center1, center2, target): + from pymbolic import parse, evaluate + context = { + "s": source, + "c1": center1, + "c2": center2, + "t": target, + "norm": la.norm} + + return evaluate(parse(expr), context) + + def __init__(self, + source, center1, center2, target, expansion1, expansion2, conv_factor): + + if isinstance(conv_factor, str): + conv_factor = self.eval(conv_factor, source, center1, center2, target) + + Record.__init__(self, + source=source, + center1=center1, + center2=center2, + target=target, + expansion1=expansion1, + expansion2=expansion2, + conv_factor=conv_factor) + + +P2E2E2P_TEST_CASES = ( + # local to local, 3D + P2E2E2PTestCase( + source=np.array([3., 4., 5.]), + center1=np.array([1., 0., 0.]), + center2=np.array([1., 3., 0.]), + target=np.array([1., 1., 1.]), + expansion1=t.local_expand, + expansion2=t.local_expand, + conv_factor="norm(t-c1)/norm(s-c1)"), + + # multipole to multipole, 3D + P2E2E2PTestCase( + source=np.array([1., 1., 1.]), + center1=np.array([1., 0., 0.]), + center2=np.array([1., 0., 3.]), + target=np.array([3., 4., 5.]), + expansion1=t.multipole_expand, + expansion2=t.multipole_expand, + conv_factor="norm(s-c2)/norm(t-c2)"), + + # multipole to local, 3D + P2E2E2PTestCase( + source=np.array([-2., 2., 1.]), + center1=np.array([-2., 5., 3.]), + center2=np.array([0., 0., 0.]), + target=np.array([0., 0., -1]), + expansion1=t.multipole_expand, + expansion2=t.local_expand, + conv_factor="norm(t-c2)/(norm(c2-c1)-norm(c1-s))"), +) + +# }}} + + +ORDERS_P2E2E2P = (3, 4, 5) +RTOL_P2E2E2P = 1e-2 + + +@pytest.mark.parametrize("case", P2E2E2P_TEST_CASES) +def test_toy_p2e2e2p(ctx_getter, case): + dim = case.dim + + src = case.source.reshape(dim, -1) + tgt = case.target.reshape(dim, -1) + + if not 0 <= case.conv_factor <= 1: + raise ValueError( + "convergence factor not in valid range: %e" % case.conv_factor) + + from sumpy.expansion.local import VolumeTaylorLocalExpansion + from sumpy.expansion.multipole import VolumeTaylorMultipoleExpansion + + cl_ctx = ctx_getter() + ctx = t.ToyContext(cl_ctx, + LaplaceKernel(dim), + VolumeTaylorMultipoleExpansion, + VolumeTaylorLocalExpansion) + + errors = [] + + src_pot = t.PointSources(ctx, src, weights=np.array([1.])) + pot_actual = src_pot.eval(tgt).item() + + for order in ORDERS_P2E2E2P: + expn = case.expansion1(src_pot, case.center1, order=order) + expn2 = case.expansion2(expn, case.center2, order=order) + pot_p2e2e2p = expn2.eval(tgt).item() + errors.append(np.abs(pot_actual - pot_p2e2e2p)) + + conv_factor = approx_convergence_factor(1 + np.array(ORDERS_P2E2E2P), errors) + assert conv_factor <= min(1, case.conv_factor * (1 + RTOL_P2E2E2P)), \ + (conv_factor, case.conv_factor * (1 + RTOL_P2E2E2P)) + + # You can test individual routines by typing # $ python test_misc.py 'test_p2p(cl.create_some_context)'