From 4e64202bafb801aa57a62562469f9695d562e5ca Mon Sep 17 00:00:00 2001 From: Andreas Kloeckner <inform@tiker.net> Date: Thu, 5 Oct 2017 10:44:56 -0500 Subject: [PATCH] Test, fix order finder --- sumpy/expansion/level_to_order.py | 19 +++++++++++++++---- test/test_misc.py | 23 ++++++++++++++++++++++- 2 files changed, 37 insertions(+), 5 deletions(-) diff --git a/sumpy/expansion/level_to_order.py b/sumpy/expansion/level_to_order.py index 93a29268..4dc22f25 100644 --- a/sumpy/expansion/level_to_order.py +++ b/sumpy/expansion/level_to_order.py @@ -140,13 +140,24 @@ class SimpleExpansionOrderFinder(object): tree.stick_out_factor * tree.root_extent / (1 << level)) + factor = box_lengthscale * helmholtz_k / (2*float(np.pi)) + from math import factorial helm_order = 1 + helm_error = self.err_const_helmholtz * factor while True: - helm_error = ( - 1/factorial(helm_order+1) - * self.err_const_helmholtz - * (box_lengthscale * helmholtz_k)**(helm_order+1)) + helm_error = helm_error * factor / (helm_order+1) + + if helm_order < 4: + # this may overflow for large orders + helm_error_direct = ( + 1/factorial(helm_order+1) + * self.err_const_helmholtz + * factor**(helm_order+1)) + print(helm_error, helm_error_direct) + assert (abs(helm_error - helm_error_direct) + < 1e-13 * abs(helm_error_direct)) + if helm_error < self.tol: break diff --git a/test/test_misc.py b/test/test_misc.py index 9edd1966..97feecd2 100644 --- a/test/test_misc.py +++ b/test/test_misc.py @@ -32,7 +32,8 @@ import pyopencl as cl # noqa: F401 from pyopencl.tools import ( # noqa pytest_generate_tests_for_pyopencl as pytest_generate_tests) -from sumpy.kernel import BiharmonicKernel, YukawaKernel +from sumpy.kernel import (LaplaceKernel, HelmholtzKernel, + BiharmonicKernel, YukawaKernel) # {{{ pde check for kernels @@ -114,6 +115,26 @@ def test_pde_check(dim, order=4): assert eoc_rec.order_estimate() > order-2-0.1 +class FakeTree: + def __init__(self, dimensions, root_extent, stick_out_factor): + self.dimensions = dimensions + self.root_extent = root_extent + self.stick_out_factor = stick_out_factor + + +@pytest.mark.parametrize("knl", [LaplaceKernel(2), HelmholtzKernel(2)]) +def test_order_finder(knl): + from sumpy.expansion.level_to_order import SimpleExpansionOrderFinder + + ofind = SimpleExpansionOrderFinder(1e-5) + + tree = FakeTree(knl.dim, 200, 0.5) + orders = [ + ofind(knl, frozenset([("k", 5)]), tree, level) + for level in range(30)] + print(orders) + + # You can test individual routines by typing # $ python test_misc.py 'test_p2p(cl.create_some_context)' -- GitLab