From c04082c578a74678f9e38bda3817b226e9f52f1b Mon Sep 17 00:00:00 2001 From: Matt Wala Date: Fri, 19 Jul 2019 16:41:27 -0500 Subject: [PATCH 1/4] Delete broken/obsolete h2dterms and l2dterms wrappers, and replace with FMMLib-based expansion order finder Closes #39 - [ ] Point requirements.txt back to pyfmmlib master after pyfmmlib!15 is merged --- requirements.txt | 2 +- sumpy/expansion/level_to_order.py | 90 ++++++++++++++----------------- test/test_misc.py | 20 ++++++- 3 files changed, 59 insertions(+), 53 deletions(-) diff --git a/requirements.txt b/requirements.txt index 8e48bc1c..ccf135bd 100644 --- a/requirements.txt +++ b/requirements.txt @@ -5,4 +5,4 @@ git+https://github.com/inducer/islpy git+https://github.com/inducer/pyopencl git+https://gitlab.tiker.net/inducer/boxtree git+https://github.com/inducer/loopy -git+https://github.com/inducer/pyfmmlib +git+https://gitlab.tiker.net/inducer/pyfmmlib@l3dterms-and-h3dterms diff --git a/sumpy/expansion/level_to_order.py b/sumpy/expansion/level_to_order.py index 8f012e55..5e29b861 100644 --- a/sumpy/expansion/level_to_order.py +++ b/sumpy/expansion/level_to_order.py @@ -23,73 +23,61 @@ THE SOFTWARE. """ __doc__ = """ -.. autofunction:: h2d_level_to_order_lookup -.. autofunction:: l2d_level_to_order_lookup - +.. autoclass:: FMMLibExpansionOrderFinder .. autoclass:: SimpleExpansionOrderFinder """ import numpy as np -def h2d_level_to_order_lookup(tree, helmholtz_k, epsilon): - """ - Compute a mapping from level number to expansion order, - Helmholtz 2D case. - - This wraps the function *h2dterms* from :mod:`pyfmmlib`. - - :arg tree: An instance of :class:`boxtree.Tree`. - :arg helmholtz_k: Helmholtz parameter - :arg epsilon: Precision - - :return: A :class:`numpy.array` of length `tree.nlevels` +class FMMLibExpansionOrderFinder(object): + r"""Return expansion orders that meet the tolerance for a given level + using routines wrapped from :mod:`pyfmmlib`. """ - if tree.dimensions != 2: - raise ValueError("tree must be 2d") - - orders = np.empty(tree.nlevels, dtype=int) - bbox_area = np.max( - tree.bounding_box[1] - tree.bounding_box[0]) ** 2 - - from pyfmmlib import h2dterms - for level in range(tree.nlevels): - nterms, ier = h2dterms(bbox_area / 2 ** level, helmholtz_k, epsilon) - if ier != 0: - raise RuntimeError( - "h2dterms returned error code {ier}".format(ier=ier)) - orders[level] = nterms - - return orders + def __init__(self, tol, extra_order=0): + """ + :arg tol: error tolerance + :arg extra_order: order increase to accommodate, say, the taking of + derivatives of the FMM expansions. + """ + self.tol = tol + self.extra_order = extra_order + def __call__(self, kernel, kernel_args, tree, level): + import pyfmmlib -def l2d_level_to_order_lookup(tree, epsilon): - """ - Compute a mapping from level number to expansion order, - Laplace 2D case. + from sumpy.kernel import LaplaceKernel, HelmholtzKernel - This wraps the function *l2dterms* from :mod:`pyfmmlib`. + assert isinstance(kernel, (LaplaceKernel, HelmholtzKernel)) + assert tree.dimensions in (2, 3) - :arg tree: An instance of :class:`boxtree.Tree`. - :arg epsilon: Precision + if isinstance(kernel, LaplaceKernel): + if tree.dimensions == 2: + nterms, ier = pyfmmlib.l2dterms(self.tol) + if ier: + raise RuntimeError("l2dterms returned error code '%d'" % ier) - :return: A :class:`numpy.array` of length `tree.nlevels` - """ + elif tree.dimensions == 3: + nterms, ier = pyfmmlib.l3dterms(self.tol) + if ier: + raise RuntimeError("l3dterms returned error code '%d'" % ier) - if tree.dimensions != 2: - raise ValueError("tree must be 2d") + elif isinstance(kernel, HelmholtzKernel): + helmholtz_k = dict(kernel_args)[kernel.helmholtz_k_name] + size = tree.root_extent / 2 ** level - from pyfmmlib import l2dterms - nterms, ier = l2dterms(epsilon) - if ier != 0: - raise RuntimeError( - "l2dterms returned error code {ier}".format(ier=ier)) + if tree.dimensions == 2: + nterms, ier = pyfmmlib.h2dterms(size, helmholtz_k, self.tol) + if ier: + raise RuntimeError("h2dterms returned error code '%d'" % ier) - orders = np.empty(tree.nlevels, dtype=int) - orders.fill(nterms) + elif tree.dimensions == 3: + nterms, ier = pyfmmlib.h3dterms(size, helmholtz_k, self.tol) + if ier: + raise RuntimeError("h3dterms returned error code '%d'" % ier) - return orders + return nterms + self.extra_order class SimpleExpansionOrderFinder(object): @@ -114,7 +102,7 @@ class SimpleExpansionOrderFinder(object): extra_order=1): """ :arg extra_order: order increase to accommodate, say, the taking of - oderivatives f the FMM expansions. + derivatives of the FMM expansions. """ self.tol = tol diff --git a/test/test_misc.py b/test/test_misc.py index 33e8535a..8866c289 100644 --- a/test/test_misc.py +++ b/test/test_misc.py @@ -124,7 +124,9 @@ class FakeTree: self.stick_out_factor = stick_out_factor -@pytest.mark.parametrize("knl", [LaplaceKernel(2), HelmholtzKernel(2)]) +@pytest.mark.parametrize("knl", [ + LaplaceKernel(2), HelmholtzKernel(2), + LaplaceKernel(3), HelmholtzKernel(3)]) def test_order_finder(knl): from sumpy.expansion.level_to_order import SimpleExpansionOrderFinder @@ -137,6 +139,22 @@ def test_order_finder(knl): print(orders) +@pytest.mark.parametrize("knl", [ + LaplaceKernel(2), HelmholtzKernel(2), + LaplaceKernel(3), HelmholtzKernel(3)]) +def test_fmmlib_order_finder(knl): + pytest.importorskip("pyfmmlib") + from sumpy.expansion.level_to_order import FMMLibExpansionOrderFinder + + ofind = FMMLibExpansionOrderFinder(1e-5) + + tree = FakeTree(knl.dim, 200, 0.5) + orders = [ + ofind(knl, frozenset([("k", 5)]), tree, level) + for level in range(30)] + print(orders) + + # {{{ expansion toys p2e2e2p test cases def approx_convergence_factor(orders, errors): -- GitLab From b93e14d952ca1335eda7da295db0b07117c560f5 Mon Sep 17 00:00:00 2001 From: Matt Wala Date: Fri, 19 Jul 2019 16:50:45 -0500 Subject: [PATCH 2/4] Remove obsolete test --- test/test_fmm.py | 33 --------------------------------- 1 file changed, 33 deletions(-) diff --git a/test/test_fmm.py b/test/test_fmm.py index 71e3f044..2bc6d69d 100644 --- a/test/test_fmm.py +++ b/test/test_fmm.py @@ -42,10 +42,6 @@ from sumpy.expansion.local import ( LaplaceConformingVolumeTaylorLocalExpansion, HelmholtzConformingVolumeTaylorLocalExpansion) -from sumpy.expansion.level_to_order import ( - h2d_level_to_order_lookup, - l2d_level_to_order_lookup) - import pytest import logging @@ -60,35 +56,6 @@ else: faulthandler.enable() -@pytest.mark.parametrize("lookup_func, extra_args", [ - (h2d_level_to_order_lookup, (100,)), - (l2d_level_to_order_lookup, ()), - ]) -def test_level_to_order_lookup(ctx_getter, lookup_func, extra_args): - pytest.importorskip("pyfmmlib") - ctx = ctx_getter() - queue = cl.CommandQueue(ctx) - - nsources = 500 - dtype = np.float64 - epsilon = 1e-9 - - from boxtree.tools import ( - make_normal_particle_array as p_normal) - - sources = p_normal(queue, nsources, 2, dtype, seed=15) - from boxtree import TreeBuilder - tb = TreeBuilder(ctx) - - tree, _ = tb(queue, sources, max_particles_in_box=30, debug=True) - - levels = lookup_func(tree, *(extra_args + (epsilon,))) - - assert len(levels) == tree.nlevels - # Should be monotonically nonincreasing. - assert all(levels[i] >= levels[i+1] for i in range(tree.nlevels - 1)) - - @pytest.mark.parametrize("knl, local_expn_class, mpole_expn_class", [ (LaplaceKernel(2), VolumeTaylorLocalExpansion, VolumeTaylorMultipoleExpansion), (LaplaceKernel(2), LaplaceConformingVolumeTaylorLocalExpansion, -- GitLab From c5cef4fb40db872e7a182f19148df9bc1d657a45 Mon Sep 17 00:00:00 2001 From: Matt Wala Date: Fri, 19 Jul 2019 16:53:29 -0500 Subject: [PATCH 3/4] Add a reasonable assertion for order finder tests --- test/test_misc.py | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/test/test_misc.py b/test/test_misc.py index 8866c289..a25b2e80 100644 --- a/test/test_misc.py +++ b/test/test_misc.py @@ -138,6 +138,9 @@ def test_order_finder(knl): for level in range(30)] print(orders) + # Order should not increase with level + assert (np.diff(orders) <= 0).all() + @pytest.mark.parametrize("knl", [ LaplaceKernel(2), HelmholtzKernel(2), @@ -154,6 +157,9 @@ def test_fmmlib_order_finder(knl): for level in range(30)] print(orders) + # Order should not increase with level + assert (np.diff(orders) <= 0).all() + # {{{ expansion toys p2e2e2p test cases -- GitLab From a94b31bbcb996ee062d0ae58651c3b971986d573 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Andreas=20Kl=C3=B6ckner?= Date: Sat, 20 Jul 2019 01:29:31 +0200 Subject: [PATCH 4/4] Point requirements.txt back to pyfmmlib master after pyfmmlib!15 (merged) is merged --- requirements.txt | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/requirements.txt b/requirements.txt index ccf135bd..8e48bc1c 100644 --- a/requirements.txt +++ b/requirements.txt @@ -5,4 +5,4 @@ git+https://github.com/inducer/islpy git+https://github.com/inducer/pyopencl git+https://gitlab.tiker.net/inducer/boxtree git+https://github.com/inducer/loopy -git+https://gitlab.tiker.net/inducer/pyfmmlib@l3dterms-and-h3dterms +git+https://github.com/inducer/pyfmmlib -- GitLab