Add preliminary BoxLib based backend.

examples/acoustics_1d_homogeneous/acoustics_1d.py

@@ -7,8 +7,8 @@
 
 Solve the (linear) acoustics equations:
 
-.. math:: 
-    p_t + K u_x & = 0 \\ 
+.. math::
+    p_t + K u_x & = 0 \\
     u_t + p_x / \rho & = 0.
 
 Here p is the pressure, u is the velocity, K is the bulk modulus,
@@ -18,20 +18,23 @@
 The final solution is identical to the initial data because both waves have
 crossed the domain exactly once.
 """
+
 from numpy import sqrt, exp, cos
 from clawpack import riemann
-    
-def setup(use_petsc=False,kernel_language='Fortran',solver_type='classic',
+
+def setup(use_petsc=False,use_boxlib=False,kernel_language='Fortran',solver_type='classic',
           outdir='./_output',weno_order=5,time_integrator='SSP104', disable_output=False):
 
     if use_petsc:
         import clawpack.petclaw as pyclaw
+    elif use_boxlib:
+        import clawpack.boxclaw as pyclaw
     else:
         from clawpack import pyclaw
 
     if kernel_language == 'Fortran':
         riemann_solver = riemann.acoustics_1D
-    elif kernel_language=='Python': 
+    elif kernel_language=='Python':
         riemann_solver = riemann.acoustics_1D_py.acoustics_1D
 
     if solver_type=='classic':
@@ -60,7 +63,7 @@ def setup(use_petsc=False,kernel_language='Fortran',solver_type='classic',
     state.problem_data['bulk']=bulk
     state.problem_data['zz']=sqrt(rho*bulk) # Impedance
     state.problem_data['cc']=sqrt(bulk/rho) # Sound speed
- 
+
     xc=domain.grid.x.centers
     beta=100; gamma=0; x0=0.75
     state.q[0,:] = exp(-beta * (xc-x0)**2) * cos(gamma * (xc - x0))
@@ -83,11 +86,11 @@ def setup(use_petsc=False,kernel_language='Fortran',solver_type='classic',
 
 
 def setplot(plotdata):
-    """ 
+    """
     Specify what is to be plotted at each frame.
     Input:  plotdata, an instance of visclaw.data.ClawPlotData.
     Output: a modified version of plotdata.
-    """ 
+    """
     plotdata.clearfigures()  # clear any old figures,axes,items data
 
     # Figure for pressure
@@ -105,7 +108,7 @@ def setplot(plotdata):
     plotitem.plotstyle = '-o'
     plotitem.color = 'b'
     plotitem.kwargs = {'linewidth':2,'markersize':5}
-    
+
     # Set up for axes in this figure:
     plotaxes = plotfigure.new_plotaxes()
     plotaxes.axescmd = 'subplot(212)'
@@ -119,7 +122,7 @@ def setplot(plotdata):
     plotitem.plotstyle = '-'
     plotitem.color = 'b'
     plotitem.kwargs = {'linewidth':3,'markersize':5}
-    
+
     return plotdata
 
 

examples/acoustics_3d_variable/acoustics_3d_interface.py

@@ -6,8 +6,8 @@
 
 Solve the variable-coefficient acoustics equations in 3D:
 
-.. math:: 
-    p_t + K(x,y,z) (u_x + v_y + w_z) & = 0 \\ 
+.. math::
+    p_t + K(x,y,z) (u_x + v_y + w_z) & = 0 \\
     u_t + p_x / \rho(x,y,z) & = 0 \\
     v_t + p_y / \rho(x,y,z) & = 0 \\
     w_t + p_z / \rho(x,y,z) & = 0 \\
@@ -18,10 +18,10 @@
 This example shows how to solve a problem with variable coefficients.
 The left and right halves of the domain consist of different materials.
 """
- 
+
 import numpy as np
 
-def setup(use_petsc=False,outdir='./_output',solver_type='classic',
+def setup(use_petsc=False,use_boxlib=False,outdir='./_output',solver_type='classic',
           mx=30,my=30,mz=30,disable_output=False,problem='heterogeneous',**kwargs):
     """
     Example python script for solving the 3d acoustics equations.
@@ -30,6 +30,8 @@ def setup(use_petsc=False,outdir='./_output',solver_type='classic',
 
     if use_petsc:
         import clawpack.petclaw as pyclaw
+    elif use_boxlib:
+        import clawpack.boxclaw as pyclaw
     else:
         from clawpack import pyclaw
 
@@ -38,7 +40,7 @@ def setup(use_petsc=False,outdir='./_output',solver_type='classic',
         solver.limiters = pyclaw.limiters.tvd.MC
     elif solver_type=='sharpclaw':
         solver = pyclaw.SharpClawSolver3D(riemann.vc_acoustics_3D)
-        
+
     else:
         raise Exception('Unrecognized solver_type.')
 
@@ -67,7 +69,7 @@ def setup(use_petsc=False,outdir='./_output',solver_type='classic',
             solver.lim_type = 1
 
         solver.limiters = [4]
-        
+
         mx=mx; my=my; mz=mz # Grid resolution
 
         zr = 1.0  # Impedance in right half
@@ -76,7 +78,7 @@ def setup(use_petsc=False,outdir='./_output',solver_type='classic',
     if problem == 'heterogeneous':
         if solver_type=='classic':
             solver.dimensional_split=False
-        
+
         solver.bc_lower[0]    =pyclaw.BC.wall
         solver.bc_lower[1]    =pyclaw.BC.wall
         solver.bc_lower[2]    =pyclaw.BC.wall
@@ -116,9 +118,9 @@ def setup(use_petsc=False,outdir='./_output',solver_type='classic',
         r = np.sqrt((X-x0)**2 + (Y-y0)**2 + (Z-z0)**2)
         width=0.1
         state.q[0,:,:,:] = (np.abs(r-0.3)<=width)*(1.+np.cos(np.pi*(r-0.3)/width))
-    else: 
+    else:
         raise Exception('Unrecognized problem name')
-        
+
     # Set initial velocities to zero
     state.q[1,:,:,:] = 0.
     state.q[2,:,:,:] = 0.
@@ -137,6 +139,5 @@ def setup(use_petsc=False,outdir='./_output',solver_type='classic',
 
 
 if __name__=="__main__":
-    import sys
     from clawpack.pyclaw.util import run_app_from_main
     output = run_app_from_main(setup)

examples/advection_1d/advection_1d.py

@@ -7,7 +7,7 @@
 
 Solve the linear advection equation:
 
-.. math:: 
+.. math::
     q_t + u q_x & = 0.
 
 Here q is the density of some conserved quantity and u is the velocity.
@@ -19,19 +19,21 @@
 import numpy as np
 from clawpack import riemann
 
-def setup(nx=100, kernel_language='Python', use_petsc=False, solver_type='classic', weno_order=5, 
+def setup(nx=100, kernel_language='Python', use_petsc=False, use_boxlib=False, solver_type='classic', weno_order=5,
           time_integrator='SSP104', outdir='./_output'):
 
     if use_petsc:
         import clawpack.petclaw as pyclaw
+    elif use_boxlib:
+        import clawpack.boxclaw as pyclaw
     else:
         from clawpack import pyclaw
 
     if kernel_language == 'Fortran':
         riemann_solver = riemann.advection_1D
     elif kernel_language == 'Python':
         riemann_solver = riemann.advection_1D_py.advection_1D
-            
+
     if solver_type=='classic':
         solver = pyclaw.ClawSolver1D(riemann_solver)
     elif solver_type=='sharpclaw':
@@ -72,9 +74,9 @@ def setup(nx=100, kernel_language='Python', use_petsc=False, solver_type='classi
     return claw
 
 def setplot(plotdata):
-    """ 
+    """
     Plot solution using VisClaw.
-    """ 
+    """
     plotdata.clearfigures()  # clear any old figures,axes,items data
 
     plotfigure = plotdata.new_plotfigure(name='q', figno=1)
@@ -90,10 +92,10 @@ def setplot(plotdata):
     plotitem.plotstyle = '-o'
     plotitem.color = 'b'
     plotitem.kwargs = {'linewidth':2,'markersize':5}
-    
+
     return plotdata
 
- 
+
 if __name__=="__main__":
     from clawpack.pyclaw.util import run_app_from_main
     output = run_app_from_main(setup,setplot)

examples/advection_2d/advection_2d.py

@@ -6,7 +6,7 @@
 
 Solve the two-dimensional advection equation
 
-.. math:: 
+.. math::
     q_t + u q_x + v q_y & = 0
 
 Here q is a conserved quantity, and (u,v) is the velocity vector.
@@ -24,12 +24,14 @@ def qinit(state):
     """
     X, Y = state.grid.p_centers
     state.q[0,:,:] = 0.9*(0.1<X)*(X<0.6)*(0.1<Y)*(Y<0.6) + 0.1
-
 
-def setup(use_petsc=False,outdir='./_output',solver_type='classic'):
+
+def setup(use_petsc=False,use_boxlib=False,outdir='./_output',solver_type='classic'):
 
     if use_petsc:
         import clawpack.petclaw as pyclaw
+    elif use_boxlib:
+        import clawpack.boxclaw as pyclaw
     else:
         from clawpack import pyclaw
 
@@ -74,9 +76,9 @@ def setup(use_petsc=False,outdir='./_output',solver_type='classic'):
 
 
 def setplot(plotdata):
-    """ 
+    """
     Plot solution using VisClaw.
-    """ 
+    """
     from clawpack.visclaw import colormaps
 
     plotdata.clearfigures()  # clear any old figures,axes,items data
@@ -96,7 +98,7 @@ def setplot(plotdata):
     plotitem.pcolor_cmin = 0.0
     plotitem.pcolor_cmax = 1.0
     plotitem.add_colorbar = True
-    
+
     # Figure for contour plot
     plotfigure = plotdata.new_plotfigure(name='contour', figno=1)
 
@@ -112,10 +114,10 @@ def setplot(plotdata):
     plotitem.contour_min = 0.01
     plotitem.contour_max = 0.99
     plotitem.amr_contour_colors = ['b','k','r']
-    
+
     return plotdata
 
-    
+
 if __name__=="__main__":
     from clawpack.pyclaw.util import run_app_from_main
     output = run_app_from_main(setup,setplot)

src/boxclaw/README

@@ -0,0 +1,17 @@
+
+BoxLib backend for PyClaw
+=========================
+
+This backend uses the C++ version of the CCSE BoxLib_ library.  To
+install BoxLib, `download BoxLib`_ and
+
+.. code-block:: sh
+
+   cd BoxLib/Src/Python
+   python setup.py build
+   python setup.py install
+
+
+.. _BoxLib: https://ccse.lbl.gov/BoxLib/index.html
+.. _`download BoxLib`: https://ccse.lbl.gov/Downloads/downloadBoxLib.html
+

src/boxclaw/__init__.py

@@ -0,0 +1,37 @@
+"""Main boxclaw package"""
+
+import os
+import logging, logging.config
+
+# Default logging configuration file
+_DEFAULT_LOG_CONFIG_PATH = os.path.join(os.path.dirname(__file__),'log.config')
+del os
+
+# Setup loggers
+logging.config.fileConfig(_DEFAULT_LOG_CONFIG_PATH)
+
+__all__ = []
+
+# Module imports
+__all__.extend(['Controller','Dimension','Patch','Domain','Solution','State','CFL','riemann'])
+from .controller import Controller
+from clawpack.boxclaw.geometry import Patch, Domain
+from clawpack.pyclaw.geometry import Dimension
+from .solution import Solution
+from .state import State
+from .cfl import CFL
+
+__all__.extend(['ClawSolver1D','ClawSolver2D','ClawSolver3D','SharpClawSolver1D','SharpClawSolver2D','SharpClawSolver3D'])
+from .classic.solver import ClawSolver1D,ClawSolver2D,ClawSolver3D
+from .sharpclaw.solver import SharpClawSolver1D,SharpClawSolver2D,SharpClawSolver3D
+
+__all__.append('BC')
+from clawpack.pyclaw.solver import BC
+
+# Sub-packages
+import limiters
+from limiters import *
+__all__.extend(limiters.__all__)
+
+import plot
+__all__.append('plot')

src/boxclaw/cfl.py

@@ -0,0 +1,32 @@
+
+class CFL(object):
+    """Parallel CFL object, responsible for computing the
+    Courant-Friedrichs-Lewy condition across all processes.
+    """
+
+    def __init__(self, global_max):
+        self._local_max = global_max
+        self._global_max = global_max
+
+    def get_global_max(self):
+        r"""
+        Compute the maximum CFL number over all processes for the current step.
+
+        This is used to determine whether the CFL condition was
+        violated and adjust the timestep.
+        """
+        import boxlib
+        self._reduce_vec.array = self._local_max
+        self._global_max = boxlib.bl[0].ReduceRealMax(self._local_max)
+        return self._global_max
+
+    def get_cached_max(self):
+        return self._global_max
+
+    def set_local_max(self,new_local_max):
+        self._local_max = new_local_max
+
+    def update_global_max(self,new_local_max):
+        import boxlib
+        self._global_max = boxlib.bl[0].ReduceRealMax(new_local_max)
+

src/boxclaw/classic/setup.py

@@ -0,0 +1,10 @@
+#!/usr/bin/env python
+
+def configuration(parent_package='',top_path=None):
+    from numpy.distutils.misc_util import Configuration
+    config = Configuration('classic', parent_package, top_path)
+    return config
+
+if __name__ == '__main__':
+    from numpy.distutils.core import setup
+    setup(**configuration(top_path='').todict())