From 769464b3e6a41211e3070672bf9e4114d402dfa4 Mon Sep 17 00:00:00 2001 From: Matthew Emmett Date: Tue, 15 Apr 2014 14:05:16 -0700 Subject: [PATCH 01/16] Add BoxClaw. --- src/boxclaw/__init__.py | 37 ++++ src/boxclaw/cfl.py | 34 ++++ src/boxclaw/classic/__init__.py | 0 src/boxclaw/classic/setup.py | 10 ++ src/boxclaw/classic/solver.py | 31 ++++ src/boxclaw/controller.py | 24 +++ src/boxclaw/geometry.py | 108 ++++++++++++ src/boxclaw/io/__init__.py | 37 ++++ src/boxclaw/io/petsc.py | 284 ++++++++++++++++++++++++++++++ src/boxclaw/limiters/__init__.py | 16 ++ src/boxclaw/log.config | 109 ++++++++++++ src/boxclaw/plot.py | 37 ++++ src/boxclaw/setup.py | 16 ++ src/boxclaw/sharpclaw/__init__.py | 0 src/boxclaw/sharpclaw/setup.py | 12 ++ src/boxclaw/sharpclaw/solver.py | 28 +++ src/boxclaw/solution.py | 7 + src/boxclaw/state.py | 44 +++++ 18 files changed, 834 insertions(+) create mode 100644 src/boxclaw/__init__.py create mode 100644 src/boxclaw/cfl.py create mode 100644 src/boxclaw/classic/__init__.py create mode 100644 src/boxclaw/classic/setup.py create mode 100644 src/boxclaw/classic/solver.py create mode 100644 src/boxclaw/controller.py create mode 100644 src/boxclaw/geometry.py create mode 100644 src/boxclaw/io/__init__.py create mode 100644 src/boxclaw/io/petsc.py create mode 100644 src/boxclaw/limiters/__init__.py create mode 100644 src/boxclaw/log.config create mode 100644 src/boxclaw/plot.py create mode 100644 src/boxclaw/setup.py create mode 100644 src/boxclaw/sharpclaw/__init__.py create mode 100644 src/boxclaw/sharpclaw/setup.py create mode 100755 src/boxclaw/sharpclaw/solver.py create mode 100644 src/boxclaw/solution.py create mode 100644 src/boxclaw/state.py diff --git a/src/boxclaw/__init__.py b/src/boxclaw/__init__.py new file mode 100644 index 000000000..401dc5173 --- /dev/null +++ b/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') diff --git a/src/boxclaw/cfl.py b/src/boxclaw/cfl.py new file mode 100644 index 000000000..98aca84c6 --- /dev/null +++ b/src/boxclaw/cfl.py @@ -0,0 +1,34 @@ +r""" +Module for the CFL object. +""" +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) + diff --git a/src/boxclaw/classic/__init__.py b/src/boxclaw/classic/__init__.py new file mode 100644 index 000000000..e69de29bb diff --git a/src/boxclaw/classic/setup.py b/src/boxclaw/classic/setup.py new file mode 100644 index 000000000..0ad8018ba --- /dev/null +++ b/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()) diff --git a/src/boxclaw/classic/solver.py b/src/boxclaw/classic/solver.py new file mode 100644 index 000000000..34a37ae38 --- /dev/null +++ b/src/boxclaw/classic/solver.py @@ -0,0 +1,31 @@ +r""" +Module containing the PetClaw solvers + +This file currently only exists so that these solvers have a different +__module__ property, used by pyclaw.solver.Solver.__init__ to +determine the containing claw_package to use. +""" + +from __future__ import absolute_import +from clawpack import pyclaw + +class ClawSolver1D(pyclaw.ClawSolver1D): + r""" + Parallel solver for 1D problems using classic Clawpack algorithms. + """ + + __doc__ += pyclaw.util.add_parent_doc(pyclaw.ClawSolver1D) + +class ClawSolver2D(pyclaw.ClawSolver2D): + r""" + Parallel solver for 2D problems using classic Clawpack algorithms. + """ + + __doc__ += pyclaw.util.add_parent_doc(pyclaw.ClawSolver2D) + +class ClawSolver3D(pyclaw.ClawSolver3D): + r""" + Parallel solver for 3D problems using classic Clawpack algorithms. + """ + + __doc__ += pyclaw.util.add_parent_doc(pyclaw.ClawSolver3D) diff --git a/src/boxclaw/controller.py b/src/boxclaw/controller.py new file mode 100644 index 000000000..837fd55fb --- /dev/null +++ b/src/boxclaw/controller.py @@ -0,0 +1,24 @@ +""" +Module for BoxClaw controller class. +""" + +from clawpack import pyclaw + +class Controller(pyclaw.controller.Controller): + """Parallel Controller Class""" + + __doc__ += pyclaw.util.add_parent_doc(pyclaw.controller.Controller) + + def __init__(self): + super(Controller,self).__init__() + + def is_proc_0(self): + import boxlib + return boxlib.rank() == 0 + + def log_info(self, str): + import logging + if self.is_proc_0(): + logging.info(str) + else: + pass diff --git a/src/boxclaw/geometry.py b/src/boxclaw/geometry.py new file mode 100644 index 000000000..5723ff2c4 --- /dev/null +++ b/src/boxclaw/geometry.py @@ -0,0 +1,108 @@ +#!/usr/bin/env python +# encoding: utf-8 +r""" +Module containing boxclaw.geometry. +""" + +from clawpack import pyclaw +from clawpack.pyclaw import geometry as pyclaw_geometry + +class Patch(pyclaw_geometry.Patch): + """Parallel Patch class.""" + + __doc__ += pyclaw.util.add_parent_doc(pyclaw_geometry.Patch) + + def __init__(self,dimensions): + + import boxlib + import numpy as np + + super(Patch,self).__init__(dimensions) + + self._ba, box = self._create_boxarray() + + is_per = np.asarray(self.num_dim * [ 1 ], np.int32) + rb = boxlib.RealBox(boxlib.lo(box), boxlib.hi(box)) + self._geom = boxlib.bl[self.num_dim].Geometry(box, rb, 0, is_per) + + # XXX: create a multifab from the boxarray to get geometry information + tmp = boxlib.MultiFab(self._ba) + + fab = None + for i in range(tmp.size()): + fab = tmp[i] + if fab is not None: + break + + assert(fab is not None) + + lo = boxlib.lo(fab.box()) + hi = boxlib.hi(fab.box()) + grid_dimensions = [] + for i in range(self.num_dim): + lower = self.lower_global[i] + lo[i]*self.delta[i] + upper = self.lower_global[i] + (hi[i]+1)*self.delta[i] + num_cells = hi[i]-lo[i]+1 + + grid_dimensions.append(pyclaw_geometry.Dimension(lower,upper, + num_cells,name=dimensions[i].name)) + + if lower == self.lower_global[i]: + grid_dimensions[-1].on_lower_boundary = True + else: + grid_dimensions[-1].on_lower_boundary = False + + if upper == self.upper_global[i]: + grid_dimensions[-1].on_upper_boundary = True + else: + grid_dimensions[-1].on_upper_boundary = False + + self.grid = pyclaw_geometry.Grid(grid_dimensions) + + del tmp + + + def _create_boxarray(self): + """Returns a BoxLib BoxArray.""" + + import boxlib + import numpy as np + + dx = self.delta + lg = self.lower_global + nc = self.num_cells_global + + lo = [ int(lg[d]/dx[d]) for d in range(self.num_dim) ] + hi = [ lo[d]+nc[d]-1 for d in range(self.num_dim) ] + + box = boxlib.Box(lo, hi) + ba = boxlib.BoxArray([box]) + + # XXX: breaking up the box array to distribute it should be + # done more intelligently + + # boxlib supports more than one box per processor; however, + # for pyclaw we will enforce one box per processor + nproc_per_dim = int(float(boxlib.size())**(1./self.num_dim)) + + max_sizes = [ (hi[d] - lo[d]) / nproc_per_dim + 1 for d in range(self.num_dim) ] + max_sizes = boxlib.bl[self.num_dim].IntVect(*max_sizes) + ba.maxSize(max_sizes) + + return ba, box + + +class Domain(pyclaw_geometry.Domain): + """Parallel Domain Class""" + + __doc__ += pyclaw.util.add_parent_doc(pyclaw.ClawSolver2D) + + def __init__(self,geom): + if not isinstance(geom,list): + geom = [geom] + if isinstance(geom[0],Patch): + self.patches = geom + elif isinstance(geom[0],pyclaw_geometry.Dimension): + self.patches = [Patch(geom)] + + diff --git a/src/boxclaw/io/__init__.py b/src/boxclaw/io/__init__.py new file mode 100644 index 000000000..55e18df1c --- /dev/null +++ b/src/boxclaw/io/__init__.py @@ -0,0 +1,37 @@ +#!/usr/bin/env python +# encoding: utf-8 +# ====================================================================== +# Package: pyclaw.io +# File: __init__.py +# Created: Feb 10, 2008 +# Author: Kyle Mandli +# ====================================================================== +"""Output package for Pyclaw""" + +import logging +from clawpack.pyclaw.io import ascii +__all__ = ['ascii.read','ascii.write'] + +# Check for HDF 5 support +try: + import h5py + from clawpack.pyclaw.io import hdf5 + __all__ += ['hdf5.read','hdf5.write'] +except: + logging.debug("No hdf5 support found.") + +# Check for netcdf support +try: + import netCDF4 + from clawpack.pyclaw.io import netcdf + __all__ += ['netcdf.read','netcdf.write'] +except(ImportError): + logging.debug("No netcdf4 support found.") + +# Check for petsc4py support +try: + import petsc + __all__ += ['petsc.read','petsc.write'] +except(ImportError): + logging.debug("No petsc support found.") + diff --git a/src/boxclaw/io/petsc.py b/src/boxclaw/io/petsc.py new file mode 100644 index 000000000..acd6ef3e6 --- /dev/null +++ b/src/boxclaw/io/petsc.py @@ -0,0 +1,284 @@ +#!/usr/bin/env python +# encoding: utf-8 +r""" +Routines for reading and writing a petsc-style output file. + +These routines preserve petclaw/pyclaw syntax for i/o while taking advantage of +PETSc's parallel i/o capabilities to allow for parallel reads and writes of +frame data. +""" + +from petsc4py import PETSc +import pickle + + +def write(solution,frame,path='./',file_prefix='claw',write_aux=False, + options={},write_p=False): + r""" + Write out pickle and PETSc data files representing the + solution. Common data is written from process 0 in pickle + files. Shared data is written from all processes into PETSc + data files. + + :Input: + - *solution* - (:class:`~pyclaw.solution.Solution`) pyclaw + object to be output + - *frame* - (int) Frame number + - *path* - (string) Root path + - *file_prefix* - (string) Prefix for the file name. ``default = + 'claw'`` + - *write_aux* - (bool) Boolean controlling whether the associated + auxiliary array should be written out. ``default = False`` + - *options* - (dict) Optional argument dictionary, see + `PETScIO Option Table`_ + + .. _`PETScIO Option Table`: + + format : one of 'ascii' or 'binary' + clobber : if True (Default), files will be overwritten + """ + import os + + # Option parsing + option_defaults = {'format':'binary','clobber':True} + + for k in option_defaults.iterkeys(): + if options.has_key(k): + pass + else: + options[k] = option_defaults[k] + + clobber = options['clobber'] + file_format = options['format'] + + pickle_filename = os.path.join(path, '%s.pkl' % file_prefix) + str(frame).zfill(4) + if file_format == 'vtk': + viewer_filename = os.path.join(path, file_prefix+str(frame).zfill(4)+'.vtk') + else: + viewer_filename = os.path.join(path, '%s.ptc' % file_prefix) + str(frame).zfill(4) + + if solution.num_aux == 0: + write_aux = False + if write_aux: + aux_filename = os.path.join(path, '%s_aux.ptc' % file_prefix) + str(frame).zfill(4) + + if not clobber: + for f in (pickle_filename, viewer_filename, aux_filename): + if os.path.exists(f): + raise IOError('Cowardly refusing to clobber %s!' % f) + + rank = PETSc.Comm.getRank(PETSc.COMM_WORLD) + if rank==0: + pickle_file = open(pickle_filename,'wb') + # explicitly dumping a dictionary here to help out anybody trying to read the pickle file + sol_dict = {'t':solution.t,'num_eqn':solution.num_eqn,'nstates':len(solution.states), + 'num_aux':solution.num_aux,'num_dim':solution.domain.num_dim, + 'write_aux':write_aux, + 'problem_data' : solution.problem_data, + 'mapc2p': solution.state.grid.mapc2p} + if write_p: + sol_dict['num_eqn'] = solution.mp + + pickle.dump(sol_dict, pickle_file) + + # now set up the PETSc viewers + if file_format == 'ascii': + viewer = PETSc.Viewer().createASCII(viewer_filename, PETSc.Viewer.Mode.WRITE) + if write_aux: + aux_viewer = PETSc.Viewer().createASCII(aux_filename, PETSc.Viewer.Mode.WRITE) + elif file_format == 'binary': + if hasattr(PETSc.Viewer,'createMPIIO'): + viewer = PETSc.Viewer().createMPIIO(viewer_filename, PETSc.Viewer.Mode.WRITE) + else: + viewer = PETSc.Viewer().createBinary(viewer_filename, PETSc.Viewer.Mode.WRITE) + if write_aux: + if hasattr(PETSc.Viewer,'createMPIIO'): + aux_viewer = PETSc.Viewer().createMPIIO(aux_filename, PETSc.Viewer.Mode.WRITE) + else: + aux_viewer = PETSc.Viewer().createBinary(aux_filename, PETSc.Viewer.Mode.WRITE) + elif file_format == 'vtk': + viewer = PETSc.Viewer().createASCII(viewer_filename, PETSc.Viewer.Mode.WRITE, format=PETSc.Viewer.Format.ASCII_VTK) + if write_aux: + aux_viewer = PETSc.Viewer().createASCII(aux_filename, PETSc.Viewer.Mode.WRITE) + else: + raise IOError('format type %s not supported' % file_format) + + for state in solution.states: + patch = state.patch + if rank==0: + pickle.dump({'level':patch.level, + 'names':patch.name,'lower':patch.lower_global, + 'num_cells':patch.num_cells_global,'delta':patch.delta}, pickle_file) +# we will reenable this bad boy when we switch over to petsc-dev +# state.q_da.view(viewer) + if write_p: + state.gpVec.view(viewer) + else: + state.gqVec.view(viewer) + + if write_aux: + state.gauxVec.view(aux_viewer) + + viewer.flush() + viewer.destroy() + if rank==0: + pickle_file.close() + if write_aux: + aux_viewer.flush() + aux_viewer.destroy() + + +def read(solution,frame,path='./',file_prefix='claw',read_aux=False,options={}): + r""" + Read in pickles and PETSc data files representing the solution + + :Input: + - *solution* - (:class:`~pyclaw.solution.Solution`) Solution object to + read the data into. + - *frame* - (int) Frame number to be read in + - *path* - (string) Path to the current directory of the file + - *file_prefix* - (string) Prefix of the files to be read in. + ``default = 'fort'`` + - *read_aux* (bool) Whether or not an auxiliary file will try to be read + in. ``default = False`` + - *options* - (dict) Optional argument dictionary, see + `PETScIO Option Table`_ + + .. _`PETScIO Option Table`: + + format : one of 'ascii' or 'binary' + + """ + import os + + if options.has_key('format'): + file_format = options['format'] + else: + file_format = 'binary' + + pickle_filename = os.path.join(path, '%s.pkl' % file_prefix) + str(frame).zfill(4) + viewer_filename = os.path.join(path, '%s.ptc' % file_prefix) + str(frame).zfill(4) + aux_viewer_filename1 = os.path.join(path, '%s_aux.ptc' % file_prefix) + str(frame).zfill(4) + aux_viewer_filename2 = os.path.join(path, '%s_aux.ptc' % file_prefix) + str(0).zfill(4) + if os.path.exists(aux_viewer_filename1): + aux_viewer_filename = aux_viewer_filename1 + else: + aux_viewer_filename = aux_viewer_filename2 + + try: + pickle_file = open(pickle_filename,'rb') + except IOError: + print "Error: file " + pickle_filename + " does not exist or is unreadable." + raise + + # this dictionary is mostly holding debugging information, only nstates is needed + # most of this information is explicitly saved in the individual patches + value_dict = pickle.load(pickle_file) + nstates = value_dict['nstates'] + num_dim = value_dict['num_dim'] + num_aux = value_dict['num_aux'] + num_eqn = value_dict['num_eqn'] + + # now set up the PETSc viewer + if file_format == 'ascii': + viewer = PETSc.Viewer().createASCII(viewer_filename, PETSc.Viewer.Mode.READ) + if read_aux: + aux_viewer = PETSc.Viewer().createASCII(aux_viewer_filename, PETSc.Viewer.Mode.READ) + elif file_format == 'binary': + if hasattr(PETSc.Viewer,'createMPIIO'): + viewer = PETSc.Viewer().createMPIIO(viewer_filename, PETSc.Viewer.Mode.READ) + else: + viewer = PETSc.Viewer().createBinary(viewer_filename, PETSc.Viewer.Mode.READ) + if read_aux: + if os.path.exists(aux_viewer_filename): + if hasattr(PETSc.Viewer,'createMPIIO'): + aux_viewer = PETSc.Viewer().createMPIIO(aux_viewer_filename, PETSc.Viewer.Mode.READ) + else: + aux_viewer = PETSc.Viewer().createBinary(aux_viewer_filename, PETSc.Viewer.Mode.READ) + else: + from warnings import warn + aux_file_path = os.path.join(path,aux_viewer_filename) + warn('read_aux=True but aux file %s does not exist' % aux_file_path) + read_aux=False + else: + raise IOError('format type %s not supported' % file_format) + + patches = [] + for m in xrange(nstates): + patch_dict = pickle.load(pickle_file) + + level = patch_dict['level'] + names = patch_dict['names'] + lower = patch_dict['lower'] + n = patch_dict['num_cells'] + d = patch_dict['delta'] + + from clawpack import petclaw + dimensions = [] + for i in xrange(num_dim): + dimensions.append( + petclaw.Dimension(names[i],lower[i],lower[i] + n[i]*d[i],n[i])) + patch = petclaw.Patch(dimensions) + patch.level = level + state = petclaw.State(patch,num_eqn,num_aux) + state.t = value_dict['t'] + state.problem_data = value_dict['problem_data'] + state.grid.mapc2p = value_dict['mapc2p'] + +# DA View/Load is broken in Petsc-3.1.8, we can load/view the DA if needed in petsc-3.2 +# state.q_da.load(viewer) + state.gqVec.load(viewer) + + if read_aux: + state.gauxVec.load(aux_viewer) + + solution.states.append(state) + patches.append(state.patch) + solution.domain = petclaw.geometry.Domain(patches) + + pickle_file.close() + viewer.destroy() + if read_aux: + aux_viewer.destroy() + +def read_t(frame,path='./',file_prefix='claw'): + r"""Read only the petsc.pkl file and return the data + + :Input: + - *frame* - (int) Frame number to be read in + - *path* - (string) Path to the current directory of the file + - *file_prefix* - (string) Prefix of the files to be read in. + ``default = 'claw'`` + + :Output: + - (list) List of output variables + - *t* - (int) Time of frame + - *num_eqn* - (int) Number of equations in the frame + - *npatches* - (int) Number of patches + - *num_aux* - (int) Auxillary value in the frame + - *num_dim* - (int) Number of dimensions in q and aux + + """ + import os + import logging + logger = logging.getLogger('io') + + base_path = os.path.join(path,) + path = os.path.join(base_path, '%s.pkl' % file_prefix) + str(frame).zfill(4) + try: + f = open(path,'rb') + except IOError: + print "Error: file " + path + " does not exist or is unreadable." + raise + logger.debug("Opening %s file." % path) + patch_dict = pickle.load(f) + + t = patch_dict['t'] + num_eqn = patch_dict['num_eqn'] + nstates = patch_dict['nstates'] + num_aux = patch_dict['num_aux'] + num_dim = patch_dict['num_dim'] + + f.close() + + return t,num_eqn,nstates,num_aux,num_dim diff --git a/src/boxclaw/limiters/__init__.py b/src/boxclaw/limiters/__init__.py new file mode 100644 index 000000000..0b07c2f9d --- /dev/null +++ b/src/boxclaw/limiters/__init__.py @@ -0,0 +1,16 @@ +#!/usr/bin/env python +# encoding: utf-8 +""" +limiters + +The BoxClaw limiters all inherit from PyClaw limiters. +""" + +# This __init__ script only imports common utilities, most of the import +# should be done depending on the solver needed + +from clawpack.riemann import * + +__all__ = ['tvd'] +from clawpack.pyclaw.limiters import tvd + diff --git a/src/boxclaw/log.config b/src/boxclaw/log.config new file mode 100644 index 000000000..812dadbaf --- /dev/null +++ b/src/boxclaw/log.config @@ -0,0 +1,109 @@ +# boxclaw Default Logging Configuration + +# ---------------------------------------------------------------------------- +# This section defines the names of the loggers, handlers and formatters +# + +# These are the names of the different loggers +[loggers] +keys=root,io,evolve,f2py,data + +# These are the names of the different handlers that we will setup later +[handlers] +keys=file,console,syslog + +# These are the formatters used for the formatters, here we only define one +# but multiple may be defined for different tasks +[formatters] +keys=default,detailed + +# ---------------------------------------------------------------------------- +# Logger settings +# +[logger_root] +level=INFO +handlers=file,console + +[logger_io] +level=INFO +propagate=0 +qualname=io +handlers=file +channel=io +parent=(root) + +[logger_solution] +level=INFO +propagate=0 +qualname=solution +handlers=file +channel=solution +partent=(root) + +[logger_plot] +level=INFO +propagate=0 +qualname=plot +handlers=file +channel=plot +parent=(root) + +[logger_evolve] +level=INFO +propagate=0 +qualname=evolve +handlers=file,console +channel=evolve +parent=(root) + +[logger_f2py] +level=INFO +propagate=0 +qualname=f2py +handlers=file,console +channel=f2py +parent=(root) + +[logger_data] +level=INFO +propagate=0 +qualname=data +handlers=file,console +channel=data +parent=(root) +# ---------------------------------------------------------------------------- + +# ---------------------------------------------------------------------------- +# Handlers +# +[handler_file] +class=FileHandler +level=NOTSET +formatter=detailed +args=('boxclaw.log','w') + +[handler_console] +class=StreamHandler +level=INFO +formatter=default +args=(sys.stdout,) + +[handler_syslog] +class=handlers.SysLogHandler +level=NOTSET +formatter=detailed +args=(('localhost',handlers.SYSLOG_UDP_PORT), handlers.SysLogHandler.LOG_USER) +# ---------------------------------------------------------------------------- + +# ---------------------------------------------------------------------------- +# Formatters +# +[formatter_default] +format=%(asctime)s %(levelname)s CLAW: %(message)s +datefmt= + +[formatter_detailed] +format=%(asctime)s %(name)s %(levelname)s CLAW: %(lineno)d - %(message)s +datefmt= + +# ---------------------------------------------------------------------------- diff --git a/src/boxclaw/plot.py b/src/boxclaw/plot.py new file mode 100644 index 000000000..e3fb14a9f --- /dev/null +++ b/src/boxclaw/plot.py @@ -0,0 +1,37 @@ +def interactive_plot(outdir='./_output',file_format='ascii',setplot=None): + """ + Convenience function for launching an interactive plotting session. + """ + from clawpack.pyclaw.plot import plot + plot(setplot,outdir=outdir,file_format=file_format,iplot=True,htmlplot=False) + +def html_plot(outdir='./_output',file_format='ascii',setplot=None): + """ + Convenience function for creating html page with plots. + """ + from clawpack.pyclaw.plot import plot + plot(setplot,outdir=outdir,file_format=file_format,htmlplot=True,iplot=False) + +# def plotPetsc(clawobj,delay=1): +# """ +# Takes either a controller or solution object and prints each frame +# using PETSc.Viewer. +# """ +# from petsc4py import PETSc +# from clawpack.pyclaw import controller, solution + +# if isinstance(clawobj,controller.Controller): +# for n in xrange(0,clawobj.num_output_times): +# sol = clawobj.frames[n] +# viewer = PETSc.Viewer.DRAW(sol.patch.gqVec.comm) +# OptDB = PETSc.Options() +# OptDB['draw_pause'] = delay +# viewer(sol.patch.gqVec) + +# elif isinstance(clawobj,solution.Solution): +# viewer = PETSc.Viewer.DRAW(clawobj.patch.gqVec.comm) +# OptDB = PETSc.Options() +# OptDB['draw_pause'] = -1 +# viewer(clawobj.patch.gqVec) + + diff --git a/src/boxclaw/setup.py b/src/boxclaw/setup.py new file mode 100644 index 000000000..20b54c43f --- /dev/null +++ b/src/boxclaw/setup.py @@ -0,0 +1,16 @@ +#!/usr/bin/env python + +def configuration(parent_package='',top_path=None): + from numpy.distutils.misc_util import Configuration + + config = Configuration('boxclaw', parent_package, top_path) + config.add_data_files('log.config') + config.add_subpackage('classic') + config.add_subpackage('sharpclaw') + config.add_subpackage('limiters') + config.add_subpackage('io') + return config + +if __name__ == '__main__': + from numpy.distutils.core import setup + setup(**configuration(top_path='').todict()) diff --git a/src/boxclaw/sharpclaw/__init__.py b/src/boxclaw/sharpclaw/__init__.py new file mode 100644 index 000000000..e69de29bb diff --git a/src/boxclaw/sharpclaw/setup.py b/src/boxclaw/sharpclaw/setup.py new file mode 100644 index 000000000..ddb8d2028 --- /dev/null +++ b/src/boxclaw/sharpclaw/setup.py @@ -0,0 +1,12 @@ +#!/usr/bin/env python + +import os + +def configuration(parent_package='',top_path=None): + from numpy.distutils.misc_util import Configuration + config = Configuration('sharpclaw', parent_package, top_path) + return config + +if __name__ == '__main__': + from numpy.distutils.core import setup + setup(**configuration(top_path='').todict()) diff --git a/src/boxclaw/sharpclaw/solver.py b/src/boxclaw/sharpclaw/solver.py new file mode 100755 index 000000000..04305280f --- /dev/null +++ b/src/boxclaw/sharpclaw/solver.py @@ -0,0 +1,28 @@ +#!/usr/bin/env python +# encoding: utf-8 +r""" +Module containing SharpClaw solvers for BoxClaw +""" + +from __future__ import absolute_import +from clawpack import pyclaw + +class SharpClawSolver1D(pyclaw.SharpClawSolver1D): + """1D parallel SharpClaw solver. + """ + + __doc__ += pyclaw.util.add_parent_doc(pyclaw.SharpClawSolver2D) + + +class SharpClawSolver2D(pyclaw.SharpClawSolver2D): + """2D parallel SharpClaw solver. + """ + + __doc__ += pyclaw.util.add_parent_doc(pyclaw.SharpClawSolver2D) + + +class SharpClawSolver3D(pyclaw.SharpClawSolver3D): + """3D parallel SharpClaw solver. + """ + + __doc__ += pyclaw.util.add_parent_doc(pyclaw.SharpClawSolver3D) diff --git a/src/boxclaw/solution.py b/src/boxclaw/solution.py new file mode 100644 index 000000000..038433796 --- /dev/null +++ b/src/boxclaw/solution.py @@ -0,0 +1,7 @@ +from clawpack import pyclaw +from clawpack.pyclaw.solution import Solution + +class Solution(Solution): + """ Parallel Solution class. + """ + __doc__ += pyclaw.util.add_parent_doc(pyclaw.Solution) diff --git a/src/boxclaw/state.py b/src/boxclaw/state.py new file mode 100644 index 000000000..b051b56cd --- /dev/null +++ b/src/boxclaw/state.py @@ -0,0 +1,44 @@ + +import clawpack.pyclaw +import numpy as np + +class State(clawpack.pyclaw.State): + """Parallel State class""" + + def get_qbc_from_q(self,num_ghost,qbc): + """ + XXX: Fills in the interior of qbc by copying q to it. + """ + + num_dim = self.patch.num_dim + + mf = self._create_multifab(self.num_eqn, num_ghost) + a = self._get_array(mf) + + if num_dim == 1: + a[num_ghost:-num_ghost,:] = np.rollaxis(self.q, self.q.ndim-1) + elif num_dim == 2: + a[num_ghost:-num_ghost,num_ghost:-num_ghost,:] = np.rollaxis(self.q, self.q.ndim-1) + elif num_dim == 3: + a[num_ghost:-num_ghost,num_ghost:-num_ghost,num_ghost:-num_ghost,:] = np.rollaxis(self.q, self.q.ndim-1) + else: + raise Exception("Assumption (1 <= num_dim <= 3) violated.") + + mf.FillBoundary() + self.patch._geom.FillPeriodicBoundary(mf) + + qbc[...] = np.rollaxis(a, self.q.ndim-1) + + return qbc + + + def _create_multifab(self, ncomp, nghost): + import boxlib + return boxlib.MultiFab(self.patch._ba, ncomp, nghost) + + + def _get_array(self, mf): + for i in range(mf.size()): + if mf[i] is not None: + return mf[i].get_array() + return None From b2de46630217801106e755af328bdd773ab9798d Mon Sep 17 00:00:00 2001 From: Matthew Emmett Date: Wed, 16 Apr 2014 10:50:48 -0700 Subject: [PATCH 02/16] BoxClaw: Fix rollaxis for multi-dimensions, tidy and add get_auxbc_from_aux. --- src/boxclaw/state.py | 55 +++++++++++++++++++++++++++++--------------- 1 file changed, 36 insertions(+), 19 deletions(-) diff --git a/src/boxclaw/state.py b/src/boxclaw/state.py index b051b56cd..2ef3c9ddc 100644 --- a/src/boxclaw/state.py +++ b/src/boxclaw/state.py @@ -9,36 +9,53 @@ def get_qbc_from_q(self,num_ghost,qbc): """ XXX: Fills in the interior of qbc by copying q to it. """ + self._fill_boundary(self.num_eqn, num_ghost, self.q, qbc) + return qbc - num_dim = self.patch.num_dim - - mf = self._create_multifab(self.num_eqn, num_ghost) - a = self._get_array(mf) - - if num_dim == 1: - a[num_ghost:-num_ghost,:] = np.rollaxis(self.q, self.q.ndim-1) - elif num_dim == 2: - a[num_ghost:-num_ghost,num_ghost:-num_ghost,:] = np.rollaxis(self.q, self.q.ndim-1) - elif num_dim == 3: - a[num_ghost:-num_ghost,num_ghost:-num_ghost,num_ghost:-num_ghost,:] = np.rollaxis(self.q, self.q.ndim-1) - else: - raise Exception("Assumption (1 <= num_dim <= 3) violated.") + def get_auxbc_from_aux(self,num_ghost,auxbc): + self._fill_boundary(self.num_aux, num_ghost, self.aux, auxbc) + return auxbc + def _fill_boundary(self, ncomp, nghost, q, qbc): + mf = self._create_multifab(ncomp, nghost) + self._copy_into_multifab(mf, nghost, q) mf.FillBoundary() self.patch._geom.FillPeriodicBoundary(mf) - - qbc[...] = np.rollaxis(a, self.q.ndim-1) - - return qbc - + self._copy_outof_multifab(mf, qbc) def _create_multifab(self, ncomp, nghost): import boxlib return boxlib.MultiFab(self.patch._ba, ncomp, nghost) - def _get_array(self, mf): for i in range(mf.size()): if mf[i] is not None: return mf[i].get_array() return None + + def _copy_into_multifab(self, mf, ng, q): + + num_dim = self.patch.num_dim + fab = self._get_array(mf) + + if ng == 0: + fab[...] = np.rollaxis(q, 0, q.ndim) + elif num_dim == 1: + fab[ng:-ng,:] = np.rollaxis(q, 0, q.ndim) + elif num_dim == 2: + fab[ng:-ng,ng:-ng,:] = np.rollaxis(q, 0, q.ndim) + elif num_dim == 3: + fab[ng:-ng,ng:-ng,ng:-ng,:] = np.rollaxis(q, 0, q.ndim) + else: + raise Exception("Assumption (1 <= num_dim <= 3) violated.") + + def _copy_outof_multifab(self, mf, q): + num_dim = self.patch.num_dim + fab = self._get_array(mf) + q[...] = np.rollaxis(fab, self.q.ndim-1) + + + + + + From 91217e6e3bef4086caa085d12516769ca6793114 Mon Sep 17 00:00:00 2001 From: Matthew Emmett Date: Wed, 16 Apr 2014 10:51:21 -0700 Subject: [PATCH 03/16] BoxClaw: Add 'multifab' write routine. --- src/boxclaw/controller.py | 1 + src/boxclaw/geometry.py | 6 +- src/boxclaw/io/__init__.py | 13 +- src/boxclaw/io/multifab.py | 138 ++++++++++++++++++ src/boxclaw/io/petsc.py | 284 ------------------------------------- 5 files changed, 148 insertions(+), 294 deletions(-) create mode 100644 src/boxclaw/io/multifab.py delete mode 100644 src/boxclaw/io/petsc.py diff --git a/src/boxclaw/controller.py b/src/boxclaw/controller.py index 837fd55fb..546099e50 100644 --- a/src/boxclaw/controller.py +++ b/src/boxclaw/controller.py @@ -11,6 +11,7 @@ class Controller(pyclaw.controller.Controller): def __init__(self): super(Controller,self).__init__() + self.output_format = 'multifab' def is_proc_0(self): import boxlib diff --git a/src/boxclaw/geometry.py b/src/boxclaw/geometry.py index 5723ff2c4..2c0632e1e 100644 --- a/src/boxclaw/geometry.py +++ b/src/boxclaw/geometry.py @@ -19,11 +19,11 @@ def __init__(self,dimensions): super(Patch,self).__init__(dimensions) - self._ba, box = self._create_boxarray() + self._ba, self._gbox = self._create_boxarray() is_per = np.asarray(self.num_dim * [ 1 ], np.int32) - rb = boxlib.RealBox(boxlib.lo(box), boxlib.hi(box)) - self._geom = boxlib.bl[self.num_dim].Geometry(box, rb, 0, is_per) + rb = boxlib.RealBox(boxlib.lo(self._gbox), boxlib.hi(self._gbox)) + self._geom = boxlib.bl[self.num_dim].Geometry(self._gbox, rb, 0, is_per) # XXX: create a multifab from the boxarray to get geometry information tmp = boxlib.MultiFab(self._ba) diff --git a/src/boxclaw/io/__init__.py b/src/boxclaw/io/__init__.py index 55e18df1c..6ec95568f 100644 --- a/src/boxclaw/io/__init__.py +++ b/src/boxclaw/io/__init__.py @@ -19,19 +19,18 @@ __all__ += ['hdf5.read','hdf5.write'] except: logging.debug("No hdf5 support found.") - + # Check for netcdf support try: import netCDF4 - from clawpack.pyclaw.io import netcdf + from clawpack.pyclaw.io import netcdf __all__ += ['netcdf.read','netcdf.write'] except(ImportError): logging.debug("No netcdf4 support found.") -# Check for petsc4py support try: - import petsc - __all__ += ['petsc.read','petsc.write'] + import multifab + __all__ += ['multifab.write'] except(ImportError): - logging.debug("No petsc support found.") - + logging.debug("No boxlib support found.") + diff --git a/src/boxclaw/io/multifab.py b/src/boxclaw/io/multifab.py new file mode 100644 index 000000000..6d3934634 --- /dev/null +++ b/src/boxclaw/io/multifab.py @@ -0,0 +1,138 @@ +#!/usr/bin/env python +# encoding: utf-8 +r""" +Routines for reading and writing multifabs. + +These routines preserve petclaw/pyclaw syntax for i/o while taking advantage of +PETSc's parallel i/o capabilities to allow for parallel reads and writes of +frame data. +""" + +import boxlib +import pickle +import numpy as np + +def write(solution,frame,path='./',file_prefix='claw',write_aux=False, + options={},write_p=False): + r""" + Write out pickle and multifab data files representing the + solution. Common data is written from process 0 in pickle + files. Shared data is written from all processes into BoxLib + plot files. + + :Input: + - *solution* - (:class:`~pyclaw.solution.Solution`) pyclaw + object to be output + - *frame* - (int) Frame number + - *path* - (string) Root path + - *file_prefix* - (string) Prefix for the file name. ``default = + 'claw'`` + - *write_aux* - (bool) Boolean controlling whether the associated + auxiliary array should be written out. ``default = False`` + - *options* - (dict) Optional argument dictionary. + """ + + import os + + # Option parsing + # option = {'format':'binary','clobber':True}.update(options) + # clobber = options['clobber'] + + clobber = True + + pickle_filename = os.path.join(path, '%s.pkl' % file_prefix) + str(frame).zfill(4) + plt_filename = os.path.join(path, '%s.plt' % file_prefix) + str(frame).zfill(4) + + # if solution.num_aux == 0: + # write_aux = False + # if write_aux: + # aux_filename = os.path.join(path, '%s_aux.ptc' % file_prefix) + str(frame).zfill(4) + + if not clobber: + for f in (pickle_filename, viewer_filename, aux_filename): + if os.path.exists(f): + raise IOError('Cowardly refusing to clobber %s!' % f) + + rank = boxlib.rank() + if rank==0: + pickle_file = open(pickle_filename,'wb') + # explicitly dumping a dictionary here to help out anybody trying to read the pickle file + sol_dict = {'t':solution.t,'num_eqn':solution.num_eqn,'nstates':len(solution.states), + 'num_aux':solution.num_aux,'num_dim':solution.domain.num_dim, + 'write_aux':write_aux, + 'problem_data' : solution.problem_data, + 'mapc2p': solution.state.grid.mapc2p} + if write_p: + sol_dict['num_eqn'] = solution.mp + + pickle.dump(sol_dict, pickle_file) + + for state in solution.states: + patch = state.patch + if rank==0: + pickle.dump({'level':patch.level, + 'names':patch.name,'lower':patch.lower_global, + 'num_cells':patch.num_cells_global,'delta':patch.delta}, pickle_file) + + mf = solution.state._create_multifab(solution.num_eqn, 0) + q = solution.state._get_array(mf) + # print q.shape + # print solution.state.q.shape + # print np.rollaxis(solution.state.q, 0, q.ndim).shape + q[...] = np.rollaxis(solution.state.q, 0, q.ndim) + + try: + if not os.path.exists(plt_filename): + os.mkdir(plt_filename) + if not os.path.exists(plt_filename + '/Level_0'): + os.mkdir(plt_filename + '/Level_0') + except: + pass + + # build multifab + mf = solution.state._create_multifab(solution.state.num_eqn, 0) + solution.state._copy_into_multifab(mf, 0, solution.state.q) + + # write header + if rank==0: + with open(os.path.join(plt_filename, 'Header'), 'w') as f: + write = lambda s: f.write(str(s) + '\n') + + plo = solution.domain.patch.lower_global + phi = solution.domain.patch.upper_global + dx = solution.domain.patch.delta + + write("HyperCLaw-V1.1") # yeah, i have no idea why... + write(solution.state.num_eqn) # number of variables + for i in range(solution.state.num_eqn): # variable names + write('q'+str(i)) + write(solution.state.num_dim) # number of dimensions + write(solution.state.t) # time + write('0') # finest amr level number + write(' '.join(map(str, plo))) # lower corner of problem domain + write(' '.join(map(str, phi))) # upper corner of problem domain + write('1') # refinement ratio + write(solution.state.patch._gbox) # grid domain + write(frame) # time step number + write(' '.join(map(str, dx))) # dx + write('0') # cartesian coords + write('0') # boundary data? (0=no) + + # now write info for each amr level (only one in this case) + write('0 %d %f' % (mf.size(), solution.state.t)) + write(frame) + for i in range(mf.size()): + bx = solution.state.patch._ba.get(i) + lo = bx.smallEnd() + hi = bx.bigEnd() + for d in range(solution.state.num_dim): + write("%lf %lf" % (plo[d] + lo[d]*dx[d], plo[d] + hi[d]*dx[d])) + write("Level_0/Cell") + + # write cell data + mf.writeOut(plt_filename + '/Level_0/Cell') + + if rank==0: + pickle_file.close() + + diff --git a/src/boxclaw/io/petsc.py b/src/boxclaw/io/petsc.py deleted file mode 100644 index acd6ef3e6..000000000 --- a/src/boxclaw/io/petsc.py +++ /dev/null @@ -1,284 +0,0 @@ -#!/usr/bin/env python -# encoding: utf-8 -r""" -Routines for reading and writing a petsc-style output file. - -These routines preserve petclaw/pyclaw syntax for i/o while taking advantage of -PETSc's parallel i/o capabilities to allow for parallel reads and writes of -frame data. -""" - -from petsc4py import PETSc -import pickle - - -def write(solution,frame,path='./',file_prefix='claw',write_aux=False, - options={},write_p=False): - r""" - Write out pickle and PETSc data files representing the - solution. Common data is written from process 0 in pickle - files. Shared data is written from all processes into PETSc - data files. - - :Input: - - *solution* - (:class:`~pyclaw.solution.Solution`) pyclaw - object to be output - - *frame* - (int) Frame number - - *path* - (string) Root path - - *file_prefix* - (string) Prefix for the file name. ``default = - 'claw'`` - - *write_aux* - (bool) Boolean controlling whether the associated - auxiliary array should be written out. ``default = False`` - - *options* - (dict) Optional argument dictionary, see - `PETScIO Option Table`_ - - .. _`PETScIO Option Table`: - - format : one of 'ascii' or 'binary' - clobber : if True (Default), files will be overwritten - """ - import os - - # Option parsing - option_defaults = {'format':'binary','clobber':True} - - for k in option_defaults.iterkeys(): - if options.has_key(k): - pass - else: - options[k] = option_defaults[k] - - clobber = options['clobber'] - file_format = options['format'] - - pickle_filename = os.path.join(path, '%s.pkl' % file_prefix) + str(frame).zfill(4) - if file_format == 'vtk': - viewer_filename = os.path.join(path, file_prefix+str(frame).zfill(4)+'.vtk') - else: - viewer_filename = os.path.join(path, '%s.ptc' % file_prefix) + str(frame).zfill(4) - - if solution.num_aux == 0: - write_aux = False - if write_aux: - aux_filename = os.path.join(path, '%s_aux.ptc' % file_prefix) + str(frame).zfill(4) - - if not clobber: - for f in (pickle_filename, viewer_filename, aux_filename): - if os.path.exists(f): - raise IOError('Cowardly refusing to clobber %s!' % f) - - rank = PETSc.Comm.getRank(PETSc.COMM_WORLD) - if rank==0: - pickle_file = open(pickle_filename,'wb') - # explicitly dumping a dictionary here to help out anybody trying to read the pickle file - sol_dict = {'t':solution.t,'num_eqn':solution.num_eqn,'nstates':len(solution.states), - 'num_aux':solution.num_aux,'num_dim':solution.domain.num_dim, - 'write_aux':write_aux, - 'problem_data' : solution.problem_data, - 'mapc2p': solution.state.grid.mapc2p} - if write_p: - sol_dict['num_eqn'] = solution.mp - - pickle.dump(sol_dict, pickle_file) - - # now set up the PETSc viewers - if file_format == 'ascii': - viewer = PETSc.Viewer().createASCII(viewer_filename, PETSc.Viewer.Mode.WRITE) - if write_aux: - aux_viewer = PETSc.Viewer().createASCII(aux_filename, PETSc.Viewer.Mode.WRITE) - elif file_format == 'binary': - if hasattr(PETSc.Viewer,'createMPIIO'): - viewer = PETSc.Viewer().createMPIIO(viewer_filename, PETSc.Viewer.Mode.WRITE) - else: - viewer = PETSc.Viewer().createBinary(viewer_filename, PETSc.Viewer.Mode.WRITE) - if write_aux: - if hasattr(PETSc.Viewer,'createMPIIO'): - aux_viewer = PETSc.Viewer().createMPIIO(aux_filename, PETSc.Viewer.Mode.WRITE) - else: - aux_viewer = PETSc.Viewer().createBinary(aux_filename, PETSc.Viewer.Mode.WRITE) - elif file_format == 'vtk': - viewer = PETSc.Viewer().createASCII(viewer_filename, PETSc.Viewer.Mode.WRITE, format=PETSc.Viewer.Format.ASCII_VTK) - if write_aux: - aux_viewer = PETSc.Viewer().createASCII(aux_filename, PETSc.Viewer.Mode.WRITE) - else: - raise IOError('format type %s not supported' % file_format) - - for state in solution.states: - patch = state.patch - if rank==0: - pickle.dump({'level':patch.level, - 'names':patch.name,'lower':patch.lower_global, - 'num_cells':patch.num_cells_global,'delta':patch.delta}, pickle_file) -# we will reenable this bad boy when we switch over to petsc-dev -# state.q_da.view(viewer) - if write_p: - state.gpVec.view(viewer) - else: - state.gqVec.view(viewer) - - if write_aux: - state.gauxVec.view(aux_viewer) - - viewer.flush() - viewer.destroy() - if rank==0: - pickle_file.close() - if write_aux: - aux_viewer.flush() - aux_viewer.destroy() - - -def read(solution,frame,path='./',file_prefix='claw',read_aux=False,options={}): - r""" - Read in pickles and PETSc data files representing the solution - - :Input: - - *solution* - (:class:`~pyclaw.solution.Solution`) Solution object to - read the data into. - - *frame* - (int) Frame number to be read in - - *path* - (string) Path to the current directory of the file - - *file_prefix* - (string) Prefix of the files to be read in. - ``default = 'fort'`` - - *read_aux* (bool) Whether or not an auxiliary file will try to be read - in. ``default = False`` - - *options* - (dict) Optional argument dictionary, see - `PETScIO Option Table`_ - - .. _`PETScIO Option Table`: - - format : one of 'ascii' or 'binary' - - """ - import os - - if options.has_key('format'): - file_format = options['format'] - else: - file_format = 'binary' - - pickle_filename = os.path.join(path, '%s.pkl' % file_prefix) + str(frame).zfill(4) - viewer_filename = os.path.join(path, '%s.ptc' % file_prefix) + str(frame).zfill(4) - aux_viewer_filename1 = os.path.join(path, '%s_aux.ptc' % file_prefix) + str(frame).zfill(4) - aux_viewer_filename2 = os.path.join(path, '%s_aux.ptc' % file_prefix) + str(0).zfill(4) - if os.path.exists(aux_viewer_filename1): - aux_viewer_filename = aux_viewer_filename1 - else: - aux_viewer_filename = aux_viewer_filename2 - - try: - pickle_file = open(pickle_filename,'rb') - except IOError: - print "Error: file " + pickle_filename + " does not exist or is unreadable." - raise - - # this dictionary is mostly holding debugging information, only nstates is needed - # most of this information is explicitly saved in the individual patches - value_dict = pickle.load(pickle_file) - nstates = value_dict['nstates'] - num_dim = value_dict['num_dim'] - num_aux = value_dict['num_aux'] - num_eqn = value_dict['num_eqn'] - - # now set up the PETSc viewer - if file_format == 'ascii': - viewer = PETSc.Viewer().createASCII(viewer_filename, PETSc.Viewer.Mode.READ) - if read_aux: - aux_viewer = PETSc.Viewer().createASCII(aux_viewer_filename, PETSc.Viewer.Mode.READ) - elif file_format == 'binary': - if hasattr(PETSc.Viewer,'createMPIIO'): - viewer = PETSc.Viewer().createMPIIO(viewer_filename, PETSc.Viewer.Mode.READ) - else: - viewer = PETSc.Viewer().createBinary(viewer_filename, PETSc.Viewer.Mode.READ) - if read_aux: - if os.path.exists(aux_viewer_filename): - if hasattr(PETSc.Viewer,'createMPIIO'): - aux_viewer = PETSc.Viewer().createMPIIO(aux_viewer_filename, PETSc.Viewer.Mode.READ) - else: - aux_viewer = PETSc.Viewer().createBinary(aux_viewer_filename, PETSc.Viewer.Mode.READ) - else: - from warnings import warn - aux_file_path = os.path.join(path,aux_viewer_filename) - warn('read_aux=True but aux file %s does not exist' % aux_file_path) - read_aux=False - else: - raise IOError('format type %s not supported' % file_format) - - patches = [] - for m in xrange(nstates): - patch_dict = pickle.load(pickle_file) - - level = patch_dict['level'] - names = patch_dict['names'] - lower = patch_dict['lower'] - n = patch_dict['num_cells'] - d = patch_dict['delta'] - - from clawpack import petclaw - dimensions = [] - for i in xrange(num_dim): - dimensions.append( - petclaw.Dimension(names[i],lower[i],lower[i] + n[i]*d[i],n[i])) - patch = petclaw.Patch(dimensions) - patch.level = level - state = petclaw.State(patch,num_eqn,num_aux) - state.t = value_dict['t'] - state.problem_data = value_dict['problem_data'] - state.grid.mapc2p = value_dict['mapc2p'] - -# DA View/Load is broken in Petsc-3.1.8, we can load/view the DA if needed in petsc-3.2 -# state.q_da.load(viewer) - state.gqVec.load(viewer) - - if read_aux: - state.gauxVec.load(aux_viewer) - - solution.states.append(state) - patches.append(state.patch) - solution.domain = petclaw.geometry.Domain(patches) - - pickle_file.close() - viewer.destroy() - if read_aux: - aux_viewer.destroy() - -def read_t(frame,path='./',file_prefix='claw'): - r"""Read only the petsc.pkl file and return the data - - :Input: - - *frame* - (int) Frame number to be read in - - *path* - (string) Path to the current directory of the file - - *file_prefix* - (string) Prefix of the files to be read in. - ``default = 'claw'`` - - :Output: - - (list) List of output variables - - *t* - (int) Time of frame - - *num_eqn* - (int) Number of equations in the frame - - *npatches* - (int) Number of patches - - *num_aux* - (int) Auxillary value in the frame - - *num_dim* - (int) Number of dimensions in q and aux - - """ - import os - import logging - logger = logging.getLogger('io') - - base_path = os.path.join(path,) - path = os.path.join(base_path, '%s.pkl' % file_prefix) + str(frame).zfill(4) - try: - f = open(path,'rb') - except IOError: - print "Error: file " + path + " does not exist or is unreadable." - raise - logger.debug("Opening %s file." % path) - patch_dict = pickle.load(f) - - t = patch_dict['t'] - num_eqn = patch_dict['num_eqn'] - nstates = patch_dict['nstates'] - num_aux = patch_dict['num_aux'] - num_dim = patch_dict['num_dim'] - - f.close() - - return t,num_eqn,nstates,num_aux,num_dim From d2aa5addbe2753825e3f1223cffda3c9775041f3 Mon Sep 17 00:00:00 2001 From: Matthew Emmett Date: Wed, 16 Apr 2014 10:59:55 -0700 Subject: [PATCH 04/16] BoxClaw: Add 'use_boxlib' arguments. --- examples/acoustics_1d_homogeneous/acoustics_1d.py | 6 ++++-- examples/advection_1d/advection_1d.py | 8 +++++--- examples/advection_2d/advection_2d.py | 6 ++++-- 3 files changed, 13 insertions(+), 7 deletions(-) diff --git a/examples/acoustics_1d_homogeneous/acoustics_1d.py b/examples/acoustics_1d_homogeneous/acoustics_1d.py index fda562552..9d375e401 100644 --- a/examples/acoustics_1d_homogeneous/acoustics_1d.py +++ b/examples/acoustics_1d_homogeneous/acoustics_1d.py @@ -18,8 +18,8 @@ The final solution is identical to the initial data because both waves have crossed the domain exactly once. """ - -def setup(use_petsc=False,kernel_language='Fortran',solver_type='classic',outdir='./_output',weno_order=5, disable_output=False): + +def setup(use_petsc=False,use_boxlib=False,kernel_language='Fortran',solver_type='classic',outdir='./_output',weno_order=5, disable_output=False): """ This example solves the 1-dimensional acoustics equations in a homogeneous medium. @@ -32,6 +32,8 @@ def setup(use_petsc=False,kernel_language='Fortran',solver_type='classic',outdir #================================================================= if use_petsc: import clawpack.petclaw as pyclaw + elif use_boxlib: + import clawpack.boxclaw as pyclaw else: from clawpack import pyclaw diff --git a/examples/advection_1d/advection_1d.py b/examples/advection_1d/advection_1d.py index f6682183f..27a674fdc 100755 --- a/examples/advection_1d/advection_1d.py +++ b/examples/advection_1d/advection_1d.py @@ -17,25 +17,27 @@ crossed the domain exactly once. """ -def setup(nx=100,kernel_language='Python', use_petsc=False, +def setup(nx=100,kernel_language='Python', use_petsc=False, use_boxlib=False, solver_type='classic', weno_order=5, outdir='./_output'): import numpy as np from clawpack import riemann if use_petsc: import clawpack.petclaw as pyclaw + elif use_boxlib: + import clawpack.boxclaw as pyclaw else: from clawpack import pyclaw if solver_type=='classic': if kernel_language == 'Fortran': solver = pyclaw.ClawSolver1D(riemann.advection_1D) - elif kernel_language=='Python': + elif kernel_language=='Python': solver = pyclaw.ClawSolver1D(riemann.advection_1D_py.advection_1D) elif solver_type=='sharpclaw': if kernel_language == 'Fortran': solver = pyclaw.SharpClawSolver1D(riemann.advection_1D) - elif kernel_language=='Python': + elif kernel_language=='Python': solver = pyclaw.SharpClawSolver1D(riemann.advection_1D_py.advection_1D) solver.weno_order=weno_order else: raise Exception('Unrecognized value of solver_type.') diff --git a/examples/advection_2d/advection_2d.py b/examples/advection_2d/advection_2d.py index bffd225af..73f34740e 100755 --- a/examples/advection_2d/advection_2d.py +++ b/examples/advection_2d/advection_2d.py @@ -29,8 +29,8 @@ def qinit(state): state.q[:,i,j] = 1.0 else: state.q[:,i,j] = 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'): """ Example python script for solving the 2d advection equation. """ @@ -38,6 +38,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 From 2348d04e16d251e66a3e30135e6f85ac8437aa09 Mon Sep 17 00:00:00 2001 From: Matthew Emmett Date: Wed, 16 Apr 2014 13:50:52 -0700 Subject: [PATCH 05/16] Check for 'multifab' in format_list in Solution.write. --- src/boxclaw/io/__init__.py | 13 ++----------- src/boxclaw/io/multifab.py | 3 +++ src/pyclaw/solution.py | 3 +++ 3 files changed, 8 insertions(+), 11 deletions(-) diff --git a/src/boxclaw/io/__init__.py b/src/boxclaw/io/__init__.py index 6ec95568f..19209e7e0 100644 --- a/src/boxclaw/io/__init__.py +++ b/src/boxclaw/io/__init__.py @@ -1,18 +1,10 @@ -#!/usr/bin/env python -# encoding: utf-8 -# ====================================================================== -# Package: pyclaw.io -# File: __init__.py -# Created: Feb 10, 2008 -# Author: Kyle Mandli -# ====================================================================== """Output package for Pyclaw""" import logging from clawpack.pyclaw.io import ascii + __all__ = ['ascii.read','ascii.write'] -# Check for HDF 5 support try: import h5py from clawpack.pyclaw.io import hdf5 @@ -20,7 +12,6 @@ except: logging.debug("No hdf5 support found.") -# Check for netcdf support try: import netCDF4 from clawpack.pyclaw.io import netcdf @@ -32,5 +23,5 @@ import multifab __all__ += ['multifab.write'] except(ImportError): - logging.debug("No boxlib support found.") + logging.debug("No multifab io support found.") diff --git a/src/boxclaw/io/multifab.py b/src/boxclaw/io/multifab.py index 6d3934634..c97423396 100644 --- a/src/boxclaw/io/multifab.py +++ b/src/boxclaw/io/multifab.py @@ -12,6 +12,9 @@ import pickle import numpy as np +import logging +logger = logging.getLogger('io') + def write(solution,frame,path='./',file_prefix='claw',write_aux=False, options={},write_p=False): r""" diff --git a/src/pyclaw/solution.py b/src/pyclaw/solution.py index 269e85608..f3889fa4d 100755 --- a/src/pyclaw/solution.py +++ b/src/pyclaw/solution.py @@ -286,6 +286,9 @@ def write(self,frame,path='./',file_format='ascii',file_prefix=None, if 'petsc' in format_list: from clawpack.petclaw import io write_func = io.petsc.write + elif 'multifab' in format_list: + from clawpack.boxclaw import io + write_func = io.multifab.write else: from clawpack.pyclaw import io write_func = io.ascii.write From aa7aaf693bc92a32772d1380cfe4e2bb99190c41 Mon Sep 17 00:00:00 2001 From: Matthew Emmett Date: Wed, 16 Apr 2014 13:51:53 -0700 Subject: [PATCH 06/16] Add 'use_boxlib' to acoustics_3d example. --- examples/acoustics_3d_variable/acoustics_3d_interface.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/examples/acoustics_3d_variable/acoustics_3d_interface.py b/examples/acoustics_3d_variable/acoustics_3d_interface.py index 361f6f518..e78b8e8ea 100755 --- a/examples/acoustics_3d_variable/acoustics_3d_interface.py +++ b/examples/acoustics_3d_variable/acoustics_3d_interface.py @@ -3,7 +3,7 @@ import numpy as np -def setup(use_petsc=False,outdir='./_output',solver_type='classic',disable_output=False,**kwargs): +def setup(use_petsc=False,use_boxlib=False,outdir='./_output',solver_type='classic',disable_output=False,**kwargs): """ Example python script for solving the 3d acoustics equations. """ @@ -11,6 +11,8 @@ def setup(use_petsc=False,outdir='./_output',solver_type='classic',disable_outpu if use_petsc: import clawpack.petclaw as pyclaw + elif use_boxlib: + import clawpack.boxclaw as pyclaw else: from clawpack import pyclaw @@ -126,6 +128,5 @@ def setup(use_petsc=False,outdir='./_output',solver_type='classic',disable_outpu if __name__=="__main__": - import sys from clawpack.pyclaw.util import run_app_from_main output = run_app_from_main(setup) From 31699ae742ee2e7f82ae740f790673c888ef70bd Mon Sep 17 00:00:00 2001 From: Matthew Emmett Date: Fri, 18 Apr 2014 13:27:24 -0700 Subject: [PATCH 07/16] BoxClaw: Fix lower/upper bounds of patches. --- src/boxclaw/geometry.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/boxclaw/geometry.py b/src/boxclaw/geometry.py index 2c0632e1e..7462e3952 100644 --- a/src/boxclaw/geometry.py +++ b/src/boxclaw/geometry.py @@ -40,8 +40,8 @@ def __init__(self,dimensions): hi = boxlib.hi(fab.box()) grid_dimensions = [] for i in range(self.num_dim): - lower = self.lower_global[i] + lo[i]*self.delta[i] - upper = self.lower_global[i] + (hi[i]+1)*self.delta[i] + lower = (lo[i]+0) * self.delta[i] + upper = (hi[i]+1) * self.delta[i] num_cells = hi[i]-lo[i]+1 grid_dimensions.append(pyclaw_geometry.Dimension(lower,upper, From 37a7d5970b2c5f45c6bca124652e9d3fcafb8e4c Mon Sep 17 00:00:00 2001 From: Matthew Emmett Date: Fri, 18 Apr 2014 13:28:46 -0700 Subject: [PATCH 08/16] BoxClaw: Add logging to geometry module. Don't let patches be deepcopied. Fix domain decomposition. --- src/boxclaw/geometry.py | 41 ++++++++++++++++++++++++++++++----------- 1 file changed, 30 insertions(+), 11 deletions(-) diff --git a/src/boxclaw/geometry.py b/src/boxclaw/geometry.py index 7462e3952..a8654a68f 100644 --- a/src/boxclaw/geometry.py +++ b/src/boxclaw/geometry.py @@ -7,11 +7,18 @@ from clawpack import pyclaw from clawpack.pyclaw import geometry as pyclaw_geometry +import logging +logger = logging.getLogger('root') + class Patch(pyclaw_geometry.Patch): """Parallel Patch class.""" __doc__ += pyclaw.util.add_parent_doc(pyclaw_geometry.Patch) + def __deepcopy__(self,memo): + # don't recreate boxarrays + return self + def __init__(self,dimensions): import boxlib @@ -19,10 +26,10 @@ def __init__(self,dimensions): super(Patch,self).__init__(dimensions) - self._ba, self._gbox = self._create_boxarray() + self._ba, self._gbox = self._create_boxarray() is_per = np.asarray(self.num_dim * [ 1 ], np.int32) - rb = boxlib.RealBox(boxlib.lo(self._gbox), boxlib.hi(self._gbox)) + rb = boxlib.RealBox(boxlib.lo(self._gbox), boxlib.hi(self._gbox)) self._geom = boxlib.bl[self.num_dim].Geometry(self._gbox, rb, 0, is_per) # XXX: create a multifab from the boxarray to get geometry information @@ -65,7 +72,10 @@ def __init__(self,dimensions): def _create_boxarray(self): """Returns a BoxLib BoxArray.""" + # note that boxlib supports more than one box per processor + import boxlib + import math import numpy as np dx = self.delta @@ -78,17 +88,28 @@ def _create_boxarray(self): box = boxlib.Box(lo, hi) ba = boxlib.BoxArray([box]) - # XXX: breaking up the box array to distribute it should be - # done more intelligently - - # boxlib supports more than one box per processor; however, - # for pyclaw we will enforce one box per processor - nproc_per_dim = int(float(boxlib.size())**(1./self.num_dim)) + max_sizes = self.num_dim * [ 1 ] + + nprocs = boxlib.size() + if nprocs % 2**self.num_dim == 0: + # divide domain into cubes + nproc_per_dim = nprocs / 2**self.num_dim + 1 + print nproc_per_dim + for d in range(self.num_dim): + max_sizes[d] = int(math.ceil(float(hi[d] - lo[d] + 1) / nproc_per_dim)) + else: + # divide domain into slabs + max_sizes[0] = int(math.ceil(float(hi[0] - lo[0] + 1) / nprocs)) + for d in range(1, self.num_dim): + max_sizes[d] = hi[d] - lo[d] + 1 - max_sizes = [ (hi[d] - lo[d]) / nproc_per_dim + 1 for d in range(self.num_dim) ] max_sizes = boxlib.bl[self.num_dim].IntVect(*max_sizes) ba.maxSize(max_sizes) + if boxlib.rank() == 0: + logger.info("max_sizes: " + str(max_sizes)) + logger.info("boxarray: " + str(ba)) + return ba, box @@ -104,5 +125,3 @@ def __init__(self,geom): self.patches = geom elif isinstance(geom[0],pyclaw_geometry.Dimension): self.patches = [Patch(geom)] - - From d558d19bd8579558442e0f367f4375b0d324c2f6 Mon Sep 17 00:00:00 2001 From: Matthew Emmett Date: Mon, 28 Apr 2014 13:09:01 -0700 Subject: [PATCH 09/16] BoxClaw: Fix Header in boxclaw.io.write so that VisIt works. --- src/boxclaw/io/multifab.py | 12 ++++++++---- src/boxclaw/state.py | 3 ++- 2 files changed, 10 insertions(+), 5 deletions(-) diff --git a/src/boxclaw/io/multifab.py b/src/boxclaw/io/multifab.py index c97423396..327a3c1ea 100644 --- a/src/boxclaw/io/multifab.py +++ b/src/boxclaw/io/multifab.py @@ -96,6 +96,8 @@ def write(solution,frame,path='./',file_prefix='claw',write_aux=False, mf = solution.state._create_multifab(solution.state.num_eqn, 0) solution.state._copy_into_multifab(mf, 0, solution.state.q) + finest_level = 0 + # write header if rank==0: with open(os.path.join(plt_filename, 'Header'), 'w') as f: @@ -105,16 +107,18 @@ def write(solution,frame,path='./',file_prefix='claw',write_aux=False, phi = solution.domain.patch.upper_global dx = solution.domain.patch.delta - write("HyperCLaw-V1.1") # yeah, i have no idea why... + write("BoxCLaw") write(solution.state.num_eqn) # number of variables for i in range(solution.state.num_eqn): # variable names write('q'+str(i)) write(solution.state.num_dim) # number of dimensions write(solution.state.t) # time - write('0') # finest amr level number + write(finest_level) # finest amr level number write(' '.join(map(str, plo))) # lower corner of problem domain write(' '.join(map(str, phi))) # upper corner of problem domain - write('1') # refinement ratio + + if finest_level == 0: + write('') write(solution.state.patch._gbox) # grid domain write(frame) # time step number write(' '.join(map(str, dx))) # dx @@ -129,7 +133,7 @@ def write(solution,frame,path='./',file_prefix='claw',write_aux=False, lo = bx.smallEnd() hi = bx.bigEnd() for d in range(solution.state.num_dim): - write("%lf %lf" % (plo[d] + lo[d]*dx[d], plo[d] + hi[d]*dx[d])) + write("%lf %lf" % (lo[d]*dx[d], (hi[d]+1)*dx[d])) write("Level_0/Cell") # write cell data diff --git a/src/boxclaw/state.py b/src/boxclaw/state.py index 2ef3c9ddc..02c47d0d5 100644 --- a/src/boxclaw/state.py +++ b/src/boxclaw/state.py @@ -7,7 +7,7 @@ class State(clawpack.pyclaw.State): def get_qbc_from_q(self,num_ghost,qbc): """ - XXX: Fills in the interior of qbc by copying q to it. + Fills in the interior of qbc by copying q to it. """ self._fill_boundary(self.num_eqn, num_ghost, self.q, qbc) return qbc @@ -22,6 +22,7 @@ def _fill_boundary(self, ncomp, nghost, q, qbc): mf.FillBoundary() self.patch._geom.FillPeriodicBoundary(mf) self._copy_outof_multifab(mf, qbc) + del mf def _create_multifab(self, ncomp, nghost): import boxlib From da0687cd0ca0eada28a51f01b268d10da490fa65 Mon Sep 17 00:00:00 2001 From: Matthew Emmett Date: Mon, 28 Apr 2014 13:21:31 -0700 Subject: [PATCH 10/16] BoxClaw: Add README and tidy up a bit. --- src/boxclaw/README | 17 +++++++++++++++++ src/boxclaw/cfl.py | 6 ++---- src/boxclaw/geometry.py | 2 +- 3 files changed, 20 insertions(+), 5 deletions(-) create mode 100644 src/boxclaw/README diff --git a/src/boxclaw/README b/src/boxclaw/README new file mode 100644 index 000000000..6a39547c1 --- /dev/null +++ b/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 + diff --git a/src/boxclaw/cfl.py b/src/boxclaw/cfl.py index 98aca84c6..9142492f6 100644 --- a/src/boxclaw/cfl.py +++ b/src/boxclaw/cfl.py @@ -1,8 +1,6 @@ -r""" -Module for the CFL object. -""" + class CFL(object): - """ Parallel CFL object, responsible for computing the + """Parallel CFL object, responsible for computing the Courant-Friedrichs-Lewy condition across all processes. """ diff --git a/src/boxclaw/geometry.py b/src/boxclaw/geometry.py index a8654a68f..b7fa1dcc5 100644 --- a/src/boxclaw/geometry.py +++ b/src/boxclaw/geometry.py @@ -1,7 +1,7 @@ #!/usr/bin/env python # encoding: utf-8 r""" -Module containing boxclaw.geometry. +Module containing BoxClaw geometry. """ from clawpack import pyclaw From d64ec40b95c4f14ffbcbc6b74f1155f6f302e058 Mon Sep 17 00:00:00 2001 From: Matthew Emmett Date: Tue, 29 Apr 2014 09:54:37 -0700 Subject: [PATCH 11/16] BoxClaw: Tidy (remove commented block). --- src/boxclaw/plot.py | 24 ------------------------ 1 file changed, 24 deletions(-) diff --git a/src/boxclaw/plot.py b/src/boxclaw/plot.py index e3fb14a9f..97322ef77 100644 --- a/src/boxclaw/plot.py +++ b/src/boxclaw/plot.py @@ -11,27 +11,3 @@ def html_plot(outdir='./_output',file_format='ascii',setplot=None): """ from clawpack.pyclaw.plot import plot plot(setplot,outdir=outdir,file_format=file_format,htmlplot=True,iplot=False) - -# def plotPetsc(clawobj,delay=1): -# """ -# Takes either a controller or solution object and prints each frame -# using PETSc.Viewer. -# """ -# from petsc4py import PETSc -# from clawpack.pyclaw import controller, solution - -# if isinstance(clawobj,controller.Controller): -# for n in xrange(0,clawobj.num_output_times): -# sol = clawobj.frames[n] -# viewer = PETSc.Viewer.DRAW(sol.patch.gqVec.comm) -# OptDB = PETSc.Options() -# OptDB['draw_pause'] = delay -# viewer(sol.patch.gqVec) - -# elif isinstance(clawobj,solution.Solution): -# viewer = PETSc.Viewer.DRAW(clawobj.patch.gqVec.comm) -# OptDB = PETSc.Options() -# OptDB['draw_pause'] = -1 -# viewer(clawobj.patch.gqVec) - - From aa50fcd0fd3690ea22d26c427e35a1a1e81000d4 Mon Sep 17 00:00:00 2001 From: Matthew Emmett Date: Wed, 18 Jun 2014 14:47:47 -0700 Subject: [PATCH 12/16] BoxClaw: Fix domain decomposition logic in 1D. --- src/boxclaw/geometry.py | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/src/boxclaw/geometry.py b/src/boxclaw/geometry.py index b7fa1dcc5..26658c6a9 100644 --- a/src/boxclaw/geometry.py +++ b/src/boxclaw/geometry.py @@ -91,10 +91,9 @@ def _create_boxarray(self): max_sizes = self.num_dim * [ 1 ] nprocs = boxlib.size() - if nprocs % 2**self.num_dim == 0: + if (self.num_dim > 1) and (nprocs % 2**self.num_dim == 0): # divide domain into cubes nproc_per_dim = nprocs / 2**self.num_dim + 1 - print nproc_per_dim for d in range(self.num_dim): max_sizes[d] = int(math.ceil(float(hi[d] - lo[d] + 1) / nproc_per_dim)) else: From 309fe7801b927dd5f9327bb4d2f0d9d660fa59e2 Mon Sep 17 00:00:00 2001 From: Matthew Emmett Date: Mon, 14 Jul 2014 13:13:35 -0700 Subject: [PATCH 13/16] boxlib: Revert whitespace changes and rename README to README.rst. --- examples/acoustics_1d_homogeneous/acoustics_1d.py | 14 +++++++------- .../acoustics_3d_interface.py | 15 ++++++++------- examples/advection_1d/advection_1d.py | 12 ++++++------ examples/advection_2d/advection_2d.py | 12 ++++++------ src/boxclaw/{README => README.rst} | 4 ++-- 5 files changed, 29 insertions(+), 28 deletions(-) rename src/boxclaw/{README => README.rst} (78%) diff --git a/examples/acoustics_1d_homogeneous/acoustics_1d.py b/examples/acoustics_1d_homogeneous/acoustics_1d.py index 5c589fbe4..adcf010df 100644 --- a/examples/acoustics_1d_homogeneous/acoustics_1d.py +++ b/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, @@ -63,7 +63,7 @@ def setup(use_petsc=False,use_boxlib=False,kernel_language='Fortran',solver_type 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)) @@ -86,11 +86,11 @@ def setup(use_petsc=False,use_boxlib=False,kernel_language='Fortran',solver_type 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 @@ -108,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)' @@ -122,7 +122,7 @@ def setplot(plotdata): plotitem.plotstyle = '-' plotitem.color = 'b' plotitem.kwargs = {'linewidth':3,'markersize':5} - + return plotdata diff --git a/examples/acoustics_3d_variable/acoustics_3d_interface.py b/examples/acoustics_3d_variable/acoustics_3d_interface.py index 729e30d7a..1886c2128 100755 --- a/examples/acoustics_3d_variable/acoustics_3d_interface.py +++ b/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 \\ @@ -40,7 +40,7 @@ def setup(use_petsc=False,use_boxlib=False,outdir='./_output',solver_type='class solver.limiters = pyclaw.limiters.tvd.MC elif solver_type=='sharpclaw': solver = pyclaw.SharpClawSolver3D(riemann.vc_acoustics_3D) - + else: raise Exception('Unrecognized solver_type.') @@ -69,7 +69,7 @@ def setup(use_petsc=False,use_boxlib=False,outdir='./_output',solver_type='class solver.lim_type = 1 solver.limiters = [4] - + mx=mx; my=my; mz=mz # Grid resolution zr = 1.0 # Impedance in right half @@ -78,7 +78,7 @@ def setup(use_petsc=False,use_boxlib=False,outdir='./_output',solver_type='class 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 @@ -118,9 +118,9 @@ def setup(use_petsc=False,use_boxlib=False,outdir='./_output',solver_type='class 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. @@ -139,5 +139,6 @@ def setup(use_petsc=False,use_boxlib=False,outdir='./_output',solver_type='class if __name__=="__main__": + import sys from clawpack.pyclaw.util import run_app_from_main output = run_app_from_main(setup) diff --git a/examples/advection_1d/advection_1d.py b/examples/advection_1d/advection_1d.py index 746e2ea4b..d9d6b93d9 100755 --- a/examples/advection_1d/advection_1d.py +++ b/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. @@ -33,7 +33,7 @@ def setup(nx=100, kernel_language='Python', use_petsc=False, use_boxlib=False, s 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': @@ -74,9 +74,9 @@ def setup(nx=100, kernel_language='Python', use_petsc=False, use_boxlib=False, s 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) @@ -92,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) diff --git a/examples/advection_2d/advection_2d.py b/examples/advection_2d/advection_2d.py index 3d33c90dd..b7b867258 100755 --- a/examples/advection_2d/advection_2d.py +++ b/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. @@ -76,9 +76,9 @@ def setup(use_petsc=False,use_boxlib=False,outdir='./_output',solver_type='class def setplot(plotdata): - """ + """ Plot solution using VisClaw. - """ + """ from clawpack.visclaw import colormaps plotdata.clearfigures() # clear any old figures,axes,items data @@ -98,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) @@ -114,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) diff --git a/src/boxclaw/README b/src/boxclaw/README.rst similarity index 78% rename from src/boxclaw/README rename to src/boxclaw/README.rst index 6a39547c1..75c2a6b8e 100644 --- a/src/boxclaw/README +++ b/src/boxclaw/README.rst @@ -7,9 +7,9 @@ install BoxLib, `download BoxLib`_ and .. code-block:: sh + git clone https://ccse.lbl.gov/pub/Downloads/BoxLib.git cd BoxLib/Src/Python - python setup.py build - python setup.py install + python setup.py install --user .. _BoxLib: https://ccse.lbl.gov/BoxLib/index.html From d13de03c73830899b9165d3aa5f771eb9d6b749a Mon Sep 17 00:00:00 2001 From: Matthew Emmett Date: Mon, 14 Jul 2014 13:20:43 -0700 Subject: [PATCH 14/16] boxlib: Revert whitespace changes in solution.py. --- src/pyclaw/solution.py | 122 ++++++++++++++++++++--------------------- 1 file changed, 61 insertions(+), 61 deletions(-) diff --git a/src/pyclaw/solution.py b/src/pyclaw/solution.py index c5e5ed2cb..866918ca4 100755 --- a/src/pyclaw/solution.py +++ b/src/pyclaw/solution.py @@ -15,18 +15,18 @@ class Solution(object): r""" Pyclaw patch container class - + :Input and Output: - + Input and output of solution objects is handle via the io package. Solution contains the generic methods :meth:`write`, :meth:`read` and :meth:`plot` which then figure out the correct method to call. Please - see the io package for the particulars of each format and method and + see the io package for the particulars of each format and method and the methods in this class for general input and output information. - + :Properties: - - If there is only one state and patch belonging to this solution, + + If there is only one state and patch belonging to this solution, the solution will appear to have many of the attributes assigned to its one state and patch. Some parameters that have in the past been parameters for all patch,s are also reachable although Solution does not @@ -36,27 +36,27 @@ class Solution(object): 'dimensions' State Attributes: 't','num_eqn','q','aux','capa','problem_data' - - + + :Initialization: - + The initialization of a Solution can happen one of these ways - + 1. `args` is empty and an empty Solution is created 2. `args` is an integer (the number of components of q), a single State, or a list of States and is followed by the appropriate :ref:`geometry ` object which can be one of: - + - (:class:`~pyclaw.geometry.Domain`) - (:class:`~pyclaw.geometry.Patch`) - A domain is created with the patch or list of patches provided. - - (:class:`~pyclaw.geometry.Dimension`) - A domain and - patch is created with the dimensions or list of + - (:class:`~pyclaw.geometry.Dimension`) - A domain and + patch is created with the dimensions or list of dimensions provided. - 3. `args` is a variable number of arguments that describes the + 3. `args` is a variable number of arguments that describes the location of a file to be read in to initialize the object - + :Examples: >>> import clawpack.pyclaw as pyclaw @@ -82,7 +82,7 @@ def __setattr__(self, key, value): self.__dict__[key] = value # ========== Attributes ================================================== - + # ========== Properties ================================================== @property def state(self): @@ -99,18 +99,18 @@ def start_frame(self): object is initialized by loading frame from file""" return self._start_frame _start_frame = 0 - + # ========== Class Methods =============================================== def __init__(self,*arg,**kargs): r"""Solution Initiatlization Routine - + See :class:`Solution` for more info. """ # select package to build solver objects from, by default this will be # the package that contains the module implementing the derived class - # for example, if Solution is implemented in 'clawpack.petclaw.solution', then + # for example, if Solution is implemented in 'clawpack.petclaw.solution', then # the computed claw_package will be 'clawpack.petclaw' import sys @@ -179,16 +179,16 @@ def get_clawpack_dot_xxx(modname): return modname.rpartition('.')[0] self.states.append(State(self.domain,arg[0])) if self.states == [] or self.domain is None: raise Exception("Invalid argument list") - - + + def is_valid(self): r""" Checks to see if this solution is valid - + The Solution checks to make sure it is valid by checking each of its states. If an invalid state is found, a message is logged what specifically made this solution invalid. - + :Output: - (bool) - True if valid, false otherwise """ @@ -201,51 +201,51 @@ def __str__(self): for state in self.states: output = output + str(state) return str(output) - - + + def set_all_states(self,attr,value,overwrite=True): r""" Sets all member states attribute 'attr' to value - + :Input: - *attr* - (string) Attribute name to be set - *value* - (id) Value for attribute - - *overwrite* - (bool) Whether to overwrite the attribute if it + - *overwrite* - (bool) Whether to overwrite the attribute if it already exists. ``default = True`` """ for state in self.states: if getattr(state,attr) is None or overwrite: - setattr(state,attr,value) - - + setattr(state,attr,value) + + def _get_base_state_attribute(self, name): r""" Return base state attribute - + :Output: - (id) - Value of attribute from ``states[0]`` """ return getattr(self.states[0],name) - - + + def __copy__(self): return self.__class__(self) - - + + def __deepcopy__(self,memo={}): import copy # Create basic container result = self.__class__() result.__init__() - + # Populate the states for state in self.states: result.states.append(copy.deepcopy(state)) result.domain = copy.deepcopy(self.domain) - + return result - - + + # ========== IO Functions ================================================ def write(self,frame,path='./',file_format='ascii',file_prefix=None, write_aux=False,options={},write_p=False): @@ -258,15 +258,15 @@ def write(self,frame,path='./',file_format='ascii',file_prefix=None, :Input: - *frame* - (int) Frame number to append to the file output - - *path* - (string) Root path, will try and create the path if it + - *path* - (string) Root path, will try and create the path if it does not already exist. ``default = './'`` - - *format* - (string or list of strings) a string or list of strings + - *format* - (string or list of strings) a string or list of strings containing the desired output formats. ``default = 'ascii'`` - *file_prefix* - (string) Prefix for the file name. Defaults to the particular io modules default. - - *write_aux* - (book) Write the auxillary array out as well if + - *write_aux* - (book) Write the auxillary array out as well if present. ``default = False`` - - *options* - (dict) Dictionary of optional arguments dependent on + - *options* - (dict) Dictionary of optional arguments dependent on which format is being used. ``default = {}`` """ # Determine if we need to create the path @@ -275,7 +275,7 @@ def write(self,frame,path='./',file_format='ascii',file_prefix=None, try: os.makedirs(path) except OSError: - print "directory already exists, ignoring" + print "directory already exists, ignoring" # Call the correct write function based on the output format if isinstance(file_format,str): @@ -306,46 +306,46 @@ def write(self,frame,path='./',file_format='ascii',file_prefix=None, msg = "Wrote out solution in format %s for time t=%s" % (form,self.t) logging.getLogger('pyclaw.io').info(msg) - + def read(self,frame,path='./_output',file_format='ascii',file_prefix=None, read_aux=True,options={}, **kargs): r""" Reads in a Solution object from a file - - Reads in and initializes this Solution with the data specified. This - function will raise an IOError if it was unsuccessful. + + Reads in and initializes this Solution with the data specified. This + function will raise an IOError if it was unsuccessful. Any format must conform to the following call signiture and return True if the file has been successfully read into the given solution or False otherwise. Options is a dictionary of parameters that each format can specify. See the ascii module for an example.:: - + read_(solution,path,frame,file_prefix,options={}) - + ```` is the name of the format in question. - + :Input: - *frame* - (int) Frame number to be read in - - *path* - (string) Base path to the files to be read. + - *path* - (string) Base path to the files to be read. ``default = './'`` - - *file_format* - (string) Format of the file, should match on of the + - *file_format* - (string) Format of the file, should match on of the modules inside of the io package. ``default = 'ascii'`` - - *file_prefix* - (string) Name prefix in front of all the files, + - *file_prefix* - (string) Name prefix in front of all the files, defaults to whatever the format defaults to, e.g. fort for ascii - - *options* - (dict) Dictionary of optional arguments dependent on + - *options* - (dict) Dictionary of optional arguments dependent on the format being read in. ``default = {}`` - + :Output: - (bool) - True if read was successful, False otherwise """ - + if file_format=='petsc': from clawpack.petclaw import io read_func = io.petsc.read elif file_format == 'binary': - from clawpack.pyclaw import io + from clawpack.pyclaw import io read_func = io.binary.read - elif file_format=='ascii': + elif file_format=='ascii': from clawpack.pyclaw import io read_func = io.ascii.read @@ -356,8 +356,8 @@ def read(self,frame,path='./_output',file_format='ascii',file_prefix=None, read_func(self,frame,path,file_prefix=file_prefix, read_aux=read_aux,options=options) logging.getLogger('pyclaw.io').info("Read in solution for time t=%s" % self.t) - - + + def plot(self): r""" Plot the solution From 2fbf81e63e51a9cb87984cbb302b407928922658 Mon Sep 17 00:00:00 2001 From: Matthew Emmett Date: Mon, 14 Jul 2014 13:56:28 -0700 Subject: [PATCH 15/16] boxlib: Add BoxClaw regressions tests for acoustics_1d_homogeneous and advection_1d examples. --- .../acoustics_1d_homogeneous/test_acoustics.py | 13 +++++++++++++ examples/advection_1d/test_advection.py | 13 ++++++++++++- src/boxclaw/util.py | 18 ++++++++++++++++++ 3 files changed, 43 insertions(+), 1 deletion(-) create mode 100644 src/boxclaw/util.py diff --git a/examples/acoustics_1d_homogeneous/test_acoustics.py b/examples/acoustics_1d_homogeneous/test_acoustics.py index 78752aeec..7b17386fb 100644 --- a/examples/acoustics_1d_homogeneous/test_acoustics.py +++ b/examples/acoustics_1d_homogeneous/test_acoustics.py @@ -46,3 +46,16 @@ def acoustics_verify(claw): from itertools import chain for test in chain(classic_tests, sharp_tests_rk, sharp_tests_lmm, weno_tests): yield test + + +if __name__=='__main__': + for test in test_1d_acoustics(): + test() + + # monkey patch util.build_variant_arg_dicts and rerun tests + import clawpack.pyclaw.util + import clawpack.boxclaw.util + clawpack.pyclaw.util.build_variant_arg_dicts = clawpack.boxclaw.util.boxlib_build_variant_arg_dicts + + for test in test_1d_acoustics(): + test() diff --git a/examples/advection_1d/test_advection.py b/examples/advection_1d/test_advection.py index a1f0e2ae1..820f7be4a 100644 --- a/examples/advection_1d/test_advection.py +++ b/examples/advection_1d/test_advection.py @@ -45,5 +45,16 @@ def advection_verify(claw): for test in chain(classic_tests, sharp_tests_rk, sharp_tests_lmm, weno_tests): yield test + if __name__=='__main__': - test_1d_advection() + for test in test_1d_advection(): + test() + + # monkey patch util.build_variant_arg_dicts and rerun tests + import clawpack.pyclaw.util + import clawpack.boxclaw.util + clawpack.pyclaw.util.build_variant_arg_dicts = clawpack.boxclaw.util.boxlib_build_variant_arg_dicts + + for test in test_1d_advection(): + test() + diff --git a/src/boxclaw/util.py b/src/boxclaw/util.py new file mode 100644 index 000000000..d8317dfc0 --- /dev/null +++ b/src/boxclaw/util.py @@ -0,0 +1,18 @@ +""" +BoxClaw related utilities. +""" + +def boxlib_build_variant_arg_dicts(kernel_languages=('Fortran',)): + import itertools + + # test petsc4py only if it is available + try: + import boxlib + except ImportError: + return [] + + opt_names = 'use_boxlib','kernel_language' + opt_product = itertools.product((True,),kernel_languages) + arg_dicts = [dict(zip(opt_names,argset)) for argset in opt_product] + + return arg_dicts From 085c8ddcd3b41156a2aed4d97d0eed11af9626cd Mon Sep 17 00:00:00 2001 From: Matthew Emmett Date: Fri, 14 Nov 2014 03:59:24 +0300 Subject: [PATCH 16/16] BoxClaw: Finish move to Fortran based BoxLib. --- src/boxclaw/__init__.py | 3 ++ src/boxclaw/cfl.py | 12 ++++---- src/boxclaw/controller.py | 4 +-- src/boxclaw/geometry.py | 59 +++++++++----------------------------- src/boxclaw/io/multifab.py | 20 ++++++------- src/boxclaw/state.py | 21 ++++---------- src/boxclaw/util.py | 4 +-- 7 files changed, 41 insertions(+), 82 deletions(-) diff --git a/src/boxclaw/__init__.py b/src/boxclaw/__init__.py index 401dc5173..3e346d3c2 100644 --- a/src/boxclaw/__init__.py +++ b/src/boxclaw/__init__.py @@ -3,6 +3,9 @@ import os import logging, logging.config +import fboxlib +fboxlib.open() + # Default logging configuration file _DEFAULT_LOG_CONFIG_PATH = os.path.join(os.path.dirname(__file__),'log.config') del os diff --git a/src/boxclaw/cfl.py b/src/boxclaw/cfl.py index 9142492f6..5a0878eb5 100644 --- a/src/boxclaw/cfl.py +++ b/src/boxclaw/cfl.py @@ -1,4 +1,6 @@ +import fboxlib + class CFL(object): """Parallel CFL object, responsible for computing the Courant-Friedrichs-Lewy condition across all processes. @@ -10,14 +12,14 @@ def __init__(self, global_max): def get_global_max(self): r""" - Compute the maximum CFL number over all processes for the current step. + 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) + self._global_max = fboxlib.reduce_max(self._local_max) return self._global_max def get_cached_max(self): @@ -27,6 +29,4 @@ 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) - + self._global_max = fboxlib.reduce_max(new_local_max) diff --git a/src/boxclaw/controller.py b/src/boxclaw/controller.py index 546099e50..ae0565c8b 100644 --- a/src/boxclaw/controller.py +++ b/src/boxclaw/controller.py @@ -2,6 +2,7 @@ Module for BoxClaw controller class. """ +import fboxlib from clawpack import pyclaw class Controller(pyclaw.controller.Controller): @@ -14,8 +15,7 @@ def __init__(self): self.output_format = 'multifab' def is_proc_0(self): - import boxlib - return boxlib.rank() == 0 + return fboxlib.mpi_rank() == 0 def log_info(self, str): import logging diff --git a/src/boxclaw/geometry.py b/src/boxclaw/geometry.py index 26658c6a9..73978385f 100644 --- a/src/boxclaw/geometry.py +++ b/src/boxclaw/geometry.py @@ -4,6 +4,9 @@ Module containing BoxClaw geometry. """ +import fboxlib +import math + from clawpack import pyclaw from clawpack.pyclaw import geometry as pyclaw_geometry @@ -16,35 +19,17 @@ class Patch(pyclaw_geometry.Patch): __doc__ += pyclaw.util.add_parent_doc(pyclaw_geometry.Patch) def __deepcopy__(self,memo): - # don't recreate boxarrays return self def __init__(self,dimensions): - import boxlib - import numpy as np - super(Patch,self).__init__(dimensions) - self._ba, self._gbox = self._create_boxarray() - - is_per = np.asarray(self.num_dim * [ 1 ], np.int32) - rb = boxlib.RealBox(boxlib.lo(self._gbox), boxlib.hi(self._gbox)) - self._geom = boxlib.bl[self.num_dim].Geometry(self._gbox, rb, 0, is_per) + self._la = self._create_layout() - # XXX: create a multifab from the boxarray to get geometry information - tmp = boxlib.MultiFab(self._ba) + bxs = self._la.local_boxes + lo, hi = self._la.get_box(bxs[0]) - fab = None - for i in range(tmp.size()): - fab = tmp[i] - if fab is not None: - break - - assert(fab is not None) - - lo = boxlib.lo(fab.box()) - hi = boxlib.hi(fab.box()) grid_dimensions = [] for i in range(self.num_dim): lower = (lo[i]+0) * self.delta[i] @@ -52,7 +37,7 @@ def __init__(self,dimensions): num_cells = hi[i]-lo[i]+1 grid_dimensions.append(pyclaw_geometry.Dimension(lower,upper, - num_cells,name=dimensions[i].name)) + num_cells,name=dimensions[i].name)) if lower == self.lower_global[i]: grid_dimensions[-1].on_lower_boundary = True @@ -66,17 +51,9 @@ def __init__(self,dimensions): self.grid = pyclaw_geometry.Grid(grid_dimensions) - del tmp - - def _create_boxarray(self): - """Returns a BoxLib BoxArray.""" - - # note that boxlib supports more than one box per processor - - import boxlib - import math - import numpy as np + def _create_layout(self): + """Returns a FBoxLib layout.""" dx = self.delta lg = self.lower_global @@ -85,12 +62,8 @@ def _create_boxarray(self): lo = [ int(lg[d]/dx[d]) for d in range(self.num_dim) ] hi = [ lo[d]+nc[d]-1 for d in range(self.num_dim) ] - box = boxlib.Box(lo, hi) - ba = boxlib.BoxArray([box]) - max_sizes = self.num_dim * [ 1 ] - - nprocs = boxlib.size() + nprocs = fboxlib.mpi_size() if (self.num_dim > 1) and (nprocs % 2**self.num_dim == 0): # divide domain into cubes nproc_per_dim = nprocs / 2**self.num_dim + 1 @@ -102,14 +75,10 @@ def _create_boxarray(self): for d in range(1, self.num_dim): max_sizes[d] = hi[d] - lo[d] + 1 - max_sizes = boxlib.bl[self.num_dim].IntVect(*max_sizes) - ba.maxSize(max_sizes) - - if boxlib.rank() == 0: - logger.info("max_sizes: " + str(max_sizes)) - logger.info("boxarray: " + str(ba)) - - return ba, box + ba = fboxlib.boxarray(boxes=[[lo, hi]]) + ba.maxsize(max_sizes) + la = fboxlib.layout(ba) + return la class Domain(pyclaw_geometry.Domain): diff --git a/src/boxclaw/io/multifab.py b/src/boxclaw/io/multifab.py index 327a3c1ea..ed588afcf 100644 --- a/src/boxclaw/io/multifab.py +++ b/src/boxclaw/io/multifab.py @@ -8,7 +8,7 @@ frame data. """ -import boxlib +import fboxlib import pickle import numpy as np @@ -56,7 +56,7 @@ def write(solution,frame,path='./',file_prefix='claw',write_aux=False, if os.path.exists(f): raise IOError('Cowardly refusing to clobber %s!' % f) - rank = boxlib.rank() + rank = fboxlib.mpi_rank() if rank==0: pickle_file = open(pickle_filename,'wb') # explicitly dumping a dictionary here to help out anybody trying to read the pickle file @@ -119,27 +119,25 @@ def write(solution,frame,path='./',file_prefix='claw',write_aux=False, if finest_level == 0: write('') - write(solution.state.patch._gbox) # grid domain + lo, hi = solution.domain.patch._la.get_box(1) + write(str(lo) + " " + str(hi)) # grid domain XXX + write('') write(frame) # time step number write(' '.join(map(str, dx))) # dx write('0') # cartesian coords write('0') # boundary data? (0=no) # now write info for each amr level (only one in this case) - write('0 %d %f' % (mf.size(), solution.state.t)) + write('0 %d %f' % (mf.nboxes, solution.state.t)) write(frame) - for i in range(mf.size()): - bx = solution.state.patch._ba.get(i) - lo = bx.smallEnd() - hi = bx.bigEnd() + for i in range(1, mf.nboxes+1): + lo, hi = solution.state.patch._la.get_box(i) for d in range(solution.state.num_dim): write("%lf %lf" % (lo[d]*dx[d], (hi[d]+1)*dx[d])) write("Level_0/Cell") # write cell data - mf.writeOut(plt_filename + '/Level_0/Cell') + mf.write(plt_filename, 'Level_0/Cell') if rank==0: pickle_file.close() - - diff --git a/src/boxclaw/state.py b/src/boxclaw/state.py index 02c47d0d5..8cb1faeaa 100644 --- a/src/boxclaw/state.py +++ b/src/boxclaw/state.py @@ -1,4 +1,5 @@ +import fboxlib import clawpack.pyclaw import numpy as np @@ -19,20 +20,15 @@ def get_auxbc_from_aux(self,num_ghost,auxbc): def _fill_boundary(self, ncomp, nghost, q, qbc): mf = self._create_multifab(ncomp, nghost) self._copy_into_multifab(mf, nghost, q) - mf.FillBoundary() - self.patch._geom.FillPeriodicBoundary(mf) + mf.fill_boundary() self._copy_outof_multifab(mf, qbc) del mf def _create_multifab(self, ncomp, nghost): - import boxlib - return boxlib.MultiFab(self.patch._ba, ncomp, nghost) + return fboxlib.multifab(self.patch._la, ncomp, nghost) def _get_array(self, mf): - for i in range(mf.size()): - if mf[i] is not None: - return mf[i].get_array() - return None + return mf.fab(1).array def _copy_into_multifab(self, mf, ng, q): @@ -51,12 +47,5 @@ def _copy_into_multifab(self, mf, ng, q): raise Exception("Assumption (1 <= num_dim <= 3) violated.") def _copy_outof_multifab(self, mf, q): - num_dim = self.patch.num_dim - fab = self._get_array(mf) + fab = self._get_array(mf) q[...] = np.rollaxis(fab, self.q.ndim-1) - - - - - - diff --git a/src/boxclaw/util.py b/src/boxclaw/util.py index d8317dfc0..d9e278d83 100644 --- a/src/boxclaw/util.py +++ b/src/boxclaw/util.py @@ -5,9 +5,9 @@ def boxlib_build_variant_arg_dicts(kernel_languages=('Fortran',)): import itertools - # test petsc4py only if it is available + # test fboxlib only if it is available try: - import boxlib + import fboxlib except ImportError: return []