Merge pull request #1 from rjleveque/absorbing_layer

Absorbing layer updates

Author: cjvogl

examples/acoustics_2d_radial/README.rst

@@ -10,3 +10,21 @@ compute the "true solution" and then setplot.py contains code to produce a
 scatter plot of the computed 2d pressure vs. distance from origin and compare
 with the 1d solution.
 
+Extrapolation BCs
+------------------
+
+The code is set up to use extrapolation boundary conditions at all
+boundaries. This does a reasonably good job of providing non-reflecting 
+boundaries, but there are some artifacts visible at later times.
+
+Absorbing boundary layer
+------------------------
+
+New cababilities have been added to Clawpack 5.5.0 to allow extending the
+computational domain with an aborbing boundary layer that does a better job
+of eliminating artificial reflections.  [Add more discussion.]
+
+To try this version::
+
+    make all -f Makefile_abl
+

examples/acoustics_2d_radial/setplot_abl.py

@@ -41,10 +41,10 @@ def setplot(plotdata=None):
     abldata.read(plotdata.outdir + '/abl.data', force=True)
     clawdata = ClawData()
     clawdata.read(plotdata.outdir + '/claw.data', force=True)
-    clawdata.lower[0] += abldata.depth_lower[0]
-    clawdata.upper[0] -= abldata.depth_upper[0]
-    clawdata.lower[1] += abldata.depth_lower[1]
-    clawdata.upper[1] -= abldata.depth_upper[1]
+    x1 = clawdata.lower[0] + abldata.depth_lower[0]
+    x2 = clawdata.upper[0] - abldata.depth_upper[0]
+    y1 = clawdata.lower[1] + abldata.depth_lower[1]
+    y2 = clawdata.upper[1] - abldata.depth_upper[1]
 
     plotdata.clearfigures()  # clear any old figures,axes,items data
 
@@ -70,13 +70,12 @@ def setplot(plotdata=None):
     plotitem.add_colorbar = True
 
     def plot_original_domain(current_data):
-        from matplotlib.pyplot import gca
+        from matplotlib.pyplot import gca, text
         ax = gca()
-        x = [clawdata.lower[0], clawdata.upper[0], clawdata.upper[0], \
-             clawdata.lower[0], clawdata.lower[0]]
-        y = [clawdata.lower[1], clawdata.lower[1], clawdata.upper[1], \
-             clawdata.upper[1], clawdata.lower[1]]
+        x = [x1,x2,x2,x1,x1]
+        y = [y1,y1,y2,y2,y1]
         ax.plot(x, y, '--k')
+        text(-0.6,1.05,'Absorbing Boundary Layer')
 
     plotaxes.afteraxes = plot_original_domain
 
@@ -102,10 +101,10 @@ def p_vs_r(current_data):
         x = current_data.x
         y = current_data.y
         r = sqrt(x**2 + y**2)
-        r = masked_where(x < clawdata.lower[0], r)
-        r = masked_where(x > clawdata.upper[0], r)
-        r = masked_where(y < clawdata.lower[1], r)
-        r = masked_where(y > clawdata.upper[1], r)
+        r = masked_where(x < x1, r)
+        r = masked_where(x > x2, r)
+        r = masked_where(y < y1, r)
+        r = masked_where(y > y2, r)
         q = current_data.q
         p = MaskedArray(q[0,:,:], mask=r.mask)
         return r,p
@@ -127,7 +126,7 @@ def p_vs_r(current_data):
 
     def make_legend(current_data):
         import matplotlib.pyplot as plt
-        plt.legend(('2d data', '1d reference solution'))
+        plt.legend(('2d data (interior only)', '1d reference solution'))
 
     plotaxes.afteraxes = make_legend