From c4442e2c39edc180b6a04a564561a82d36e4a997 Mon Sep 17 00:00:00 2001
From: Andreas Kloeckner <inform@tiker.net>
Date: Fri, 26 May 2017 10:29:32 -0700
Subject: [PATCH] Expansion toys: restrict_inner, restrict_outer, l_inf

---
 examples/expansion-toys.py |  8 +++-----
 sumpy/toys.py              | 32 +++++++++++++++++++++++++++++---
 2 files changed, 32 insertions(+), 8 deletions(-)

diff --git a/examples/expansion-toys.py b/examples/expansion-toys.py
index 72be2c7e..782bbcf5 100644
--- a/examples/expansion-toys.py
+++ b/examples/expansion-toys.py
@@ -28,16 +28,14 @@ def main():
 
     diff = mexp - pt_src
     diff = mexp2 - pt_src
-    diff = lexp - pt_src
+    diff = lexp2 - pt_src
 
+    print(t.l_inf(diff, 1.2, center=lexp2.center))
     if 1:
-        t.logplot(fp, diff, cmap="jet", vmin=-3)
+        t.logplot(fp, diff, cmap="jet", vmin=-3, vmax=0)
         plt.colorbar()
         plt.show()
 
 
-
-
-
 if __name__ == "__main__":
     main()
diff --git a/sumpy/toys.py b/sumpy/toys.py
index f61a9090..64431538 100644
--- a/sumpy/toys.py
+++ b/sumpy/toys.py
@@ -266,8 +266,8 @@ class ConstantPotential(PotentialSource):
 
 class OneOnBallPotential(PotentialSource):
     def __init__(self, toy_ctx, center, radius):
-        super(PointSources, self).__init__(toy_ctx)
-        self.center = center
+        super(OneOnBallPotential, self).__init__(toy_ctx)
+        self.center = np.asarray(center)
         self.radius = radius
 
     def eval(self, targets):
@@ -299,7 +299,7 @@ class PointSources(PotentialSource):
 class ExpansionPotentialSource(PotentialSource):
     def __init__(self, toy_ctx, center, order, coeffs):
         super(ExpansionPotentialSource, self).__init__(toy_ctx)
-        self.center = center
+        self.center = np.asarray(center)
         self.order = order
         self.coeffs = coeffs
 
@@ -397,4 +397,30 @@ def logplot(fp, psource, **kwargs):
     fp.show_scalar_in_matplotlib(
             np.log10(np.abs(psource.eval(fp.points) + 1e-15)), **kwargs)
 
+
+def restrict_inner(psource, radius, center=None):
+    if center is None:
+        center = psource.center
+
+    return psource * OneOnBallPotential(psource.toy_ctx, center, radius)
+
+
+def restrict_outer(psource, radius, center=None):
+    if center is None:
+        center = psource.center
+
+    return psource * (1-OneOnBallPotential(psource.toy_ctx, center, radius))
+
+
+def l_inf(psource, radius, center=None, npoints=100):
+    if center is None:
+        center = psource.center
+
+    restr = psource * OneOnBallPotential(psource.toy_ctx, center, radius)
+
+    from sumpy.visualization import FieldPlotter
+    fp = FieldPlotter(center, extent=2*radius, npoints=npoints)
+    return np.max(np.abs(restr.eval(fp.points)))
+
+
 # vim: foldmethod=marker
-- 
GitLab