diff --git a/sumpy/m2p.py b/sumpy/m2p.py index b86009ed7129f7e5cbd64b9dc9c94cd1c5fa235b..d3e2bb0641118a5dbded4b35dc012e17c80074ab 100644 --- a/sumpy/m2p.py +++ b/sumpy/m2p.py @@ -118,7 +118,7 @@ void m2p( % endfor loc_mpole_base += ${dimensions}; - #define COEFF(i) mpole_coeff_l[loc_mpole_base+${dimensions}+i] + #define COEFF(i) mpole_coeff_l[loc_mpole_base+i] % for var, expr in vars_and_exprs: % if var.startswith("output"): diff --git a/sumpy/symbolic/__init__.py b/sumpy/symbolic/__init__.py index 6376a158b5d6d45053c2f9d311e26b8eb626d834..fe51e6aaeb8b7c45ef8edfa8e076a23806bc7098 100644 --- a/sumpy/symbolic/__init__.py +++ b/sumpy/symbolic/__init__.py @@ -175,8 +175,10 @@ def make_coulomb_kernel_in(var_name, dimensions): if dimensions == 2: return sp.log(sp.sqrt((dist.T*dist)[0,0])) - else: + elif dimensions == 3: return 1/sp.sqrt((dist.T*dist)[0,0]) + else: + raise RuntimeError("unsupported dimensionality") diff --git a/test/test_kernels.py b/test/test_kernels.py index 81e56c9ec4c27df078cb3e6a23d62fbedef99246..755b19687781e433787ac7d82dc708ea9551d025 100644 --- a/test/test_kernels.py +++ b/test/test_kernels.py @@ -74,7 +74,7 @@ def test_p2m2p(ctx_getter): ctx = ctx_getter() queue = cl.CommandQueue(ctx) - res = 10 + res = 100 dimensions = 3 sources = np.random.rand(dimensions, 5).astype(np.float32) @@ -102,7 +102,7 @@ def test_p2m2p(ctx_getter): from sumpy.expansion import TaylorExpansion texp = TaylorExpansion( make_coulomb_kernel_in("b", dimensions), - order=2, dimensions=dimensions) + order=3, dimensions=dimensions) coeff_dtype = np.float32 @@ -150,17 +150,11 @@ def test_p2m2p(ctx_getter): if 0: import matplotlib.pyplot as pt - pt.imshow(potential_dev_direct.get().reshape(res, res)) + pt.imshow((potential_dev-potential_dev_direct).get().reshape(res, res)) + pt.colorbar() pt.show() - pt.imshow(potential_dev.get().reshape(res, res)) - pt.show() - - print potential_dev-potential_dev_direct - print potential_dev - print potential_dev_direct - #print mpole_coeff - # }}} + assert la.norm((potential_dev-potential_dev_direct).get())/res**2 < 1e-3