Setting Up Your Own Problem
Castro problems are organized loosely into groups describing their
intent (e.g., science, hydro tests, …). These groups are
sub-directories under the Castro/Exec/
directory. Each problem is
then placed in a sub-directory of the appropriate group (for example,
Castro/Exec/hydro_tests/Sedov
holds the Sedov test problem).
To create a new problem, you will create a new directory under one of the groups and place in it the following files:
GNUmakefile
: the makefile for this problem. This will tell Castro what options to use and what network and EOS to build.
problem_initialize.H
andproblem_initialize_state_data.H
: this holds the problem initialization routines. MHD and radiation problems require an additional file.
_prob_params
(optional) : a list of runtime parameters that you problem will read. These parameters are controlled by the set in the inputs file asproblem.param
whereparam
is one of the parameters listed in_prob_params
.
Make.package
: this is a makefile fragment that is included during the build process. It tells the build system about any problem-specific files that need to be compiled.
inputs
: this is the main inputs file that controls Castro and AMReX’s behavior.
The best way to get started writing a new problem is to copy an existing problem setup into a new directory.
Runtime Parameters
The problem-specific runtime parameters are defined in _prob_params
.
This has the form:
name datatype default namelist? size
Here:
name is the name of the variable to put into
probdata_module
datatype is one of
real
,integer
,character
, orlogical
.default is the default value of the runtime parameter. It may be overridden at runtime by reading from the namelist.
namelist indicates if the variable should be able to be set as a runtime parameter. If it is empty or marked with
n
, then the variable will still be creates in theproblem
namespace, but it will not be able to be set via the commandline or inputs file. A common usage of this is to define global variables that might be set at problem initialization that are used elsewhere in Castro.size is for arrays, and gives their size. It can be any integer or variable that is known to the
probdata_module
. If you need a variable from a different module to specify the size (for example,nspec
fromnetwork
), then you express the size as a tuple:(nspec, network)
and the appropriate use statement will be added.
The variables will all be initialized for the GPU as well.
Problem Initialization
Here we describe the main problem initialization routines.
initialize_problem()
This C++ routine is primarily responsible for doing any one-time initialization for the problem (like reading in an initial model through the
model_parser.H
functionality—see thetoy_convect
problem setup for an example.This routine can also postprocess any of the parameters defined in the
_prob_params
file, which are defined inproblem
namespace.Note
many problems set the value of the
problem::center[]
array from theproblem
namespace. This is used to note the center of the problem (which does not necessarily need to be the center of the domain, e.g., for axisymmetric problems).center
is used in source terms (including rotation and gravity) and in computing some of the derived variables (like angular momentum).If you need coordinate information, it can be obtained by constructing a
Geometry
object usingDefaultGeometry()
and accessing itsProbLo()
andProbHi()
methods.problem_initialize_state_data()
:This routine will initialize the state data in a given zone. The arguments passed are:
i
,j
,k
: the index of the zone to fill the data instate
: an array containing the simulation state dataGeomData
: aGeometryData
object that can be used for obtainingdx
,problo
,probhi
, etc.
Filling data is done by simply writing to
state(i,j,k,URHO)
, etc.
Boundary Conditions
Standard boundary conditions, including outflow (zero-gradient), periodic,
and symmetry (reflect) are handled by AMReX directly. Castro has a special
hydrostatic boundary condition that can be used for the lower boundary. It
is accessed by setting the castro.lo_bc
flag to 1 in the vertical coordinate
direction, e.g., for 2-d as:
castro.lo_bc = 0 1
The flag value 1 is traditionally named “inflow” by AMReX, but generally means that the boundary implementation is left to the user. To tell Castro to use the hydrostatic boundary condition here, we set:
castro.yl_ext_bc_type = 1
castro.hse_interp_temp = 1
castro.hse_reflect_vels = 1
The first parameter tells Castro to use the HSE boundary condition for the lower
y direction.
In filling the ghost cells, hydrostatic equilibrium will be integrated
from the last interior zone into the boundary. We need one more
equation for this integration, so we either interpolate the density or
temperature into the ghost cells, depending on the value of
castro.hse_interp_temp
. Finally, castro.hse_reflect_vels
determines how we treat the velocity. The default is to give is a
zero gradient, but in tests we’ve found that reflecting the velocity
while integrating the HSE profile can be better. For modeling a
plane-parallel hydrostatic atmosphere, using the hydrostatic boundary
conditions instead of a simple symmetry boundary is essential when
using the standard CTU PPM solver.
A different special boundary condition, based on outflow, is available at
the upper boundary. This works together with the model_parser
module to fill the ghost cells at the upper boundary with the initial
model data. You set this as:
castro.hi_bc = 2 2
castro.fill_ambient_bc = 1
castro.ambient_fill_dir = 1
castro.ambient_outflow_vel = 1
where ambient_fill_dir
is the 0-based direction to fill using an
ambient state defined by the problem setup. In this example, we will
override the outflow (2) boundary condition in the y-direction. That
problem setup needs to fill the ambient_state[:]
array defined in
ambient.H
. An example of using this boundary is in the
flame_wave
problem.
The implementations of these boundary conditions is found in
Castro/Source/problems/Castro_bc_fill_nd.cpp
.
Optional Files
The follow problem-specific files are optional. There are stubs for each of these in the main source tree.
problem_checkpoint.H
,problem_restart.H
:These provides two routines, respectively
problem_checkpoint
andproblem_restart
that can be used to add information to the checkpoint files and read it in upon restart. This is useful for some global problem-specific quantities. For instance, thewdmerger
problem uses this to store center of mass position and velocity information in the checkpoint files that are used for runtime diagnostics.The name of the checkpoint directory is passed in as an argument.
problem_tagging.H
This implements problem-specific tagging for refinement, through a the function
problem_tagging
. The full hydrodynamic state (State_Type) is passed in, and the problem can mark zones for refinement by setting the tag variable for a zone to set. An example is provided by thetoy_convect
problem which refines a rectangular region (fuel layer) based on a density parameter and the H mass fraction.Problem_Derives.H
,Problem_Derive.H
, andProblem_Derive.cpp
Together, these provide a mechanism to create derived quantities that can be stored in the plotfile.
Problem_Derives.H
provides the C++ code that defines these new plot variables. It does this by adding them to thederive_lst
—a list of derived variables that Castro knows about. When adding new variables, a descriptive name, a C++ routine that does the deriving, and component ofStateData
are specified.The other two files provide the header and implementation of the function that computes the derived variable. A example is provided by the
reacting_bubble
problem, which derives several new quantities (perturbations against a background one-dimensional model, in this case).Prob.cpp
,Problem.H
These files provide problem-specific routines for computing global diagnostic information through the sum_integrated_quantities functionality that is part of the
Castro
class.An example is provided by
toy_flame
, where an estimate of the flame speed is computed by integrating the mass of fuel on the grid.
Model Parser
Many problem setups begin with a 1-d initial model that is mapped onto
the grid. The model_parser.H
provides the functions that read in
the initial model and map it on the Castro grid. To enable this, add:
USE_MODEL_PARSER = TRUE
to the problem GNUmakefile
. There are 2 other parameters that can
be set in the makefile to control the initial model storage:
MAX_NPTS_MODEL
: is the maximum number of data points in the 1-d initial model. This needs to be known at compile time so we can make the data managed for GPUs.
NUM_MODELS
: this is the number of different initial models we want to managed. Typically we only want 1, but some problems, likeflame_wave
use 2, applied to different portions of the domain.
The general form of the initial model is:
# npts = 896
# num of variables = 6
# density
# temperature
# pressure
# carbon-12
# oxygen-16
# magnesium-24
195312.5000 5437711139. 8805500.952 .4695704813E+28 0.3 0.7 0
585937.5000 5410152416. 8816689.836 0.4663923963E+28 0.3 0.7 0
The first line gives the number of points in the initial model, the next gives the number of variables, followed by the variable names (one per line), and then the data. The data begins with the coordinate and then the variables in the model, with one data point per line.
When the model is read, the variables listed in the file are matched to the ones that Castro knows about. If the variable is recognized, then it is stored in the model data, otherwise, it is ignored.
The data can then be mapped onto the grid using the interpolate()
function, e.g.,
Real dens = interpolate(height, model::idens);
This fills dens
with the density at the position height
. In
addition to density, you can specify temperature (model::itemp
),
pressure (model::ipres
), species (indexed from model::ispec
),
or an auxiliary quantity (indexed from model::iaux
).