Plotfile

AMReX has its own native plotfile format. Many visualization tools are available for AMReX plotfiles (see the chapter on Visualization). AMReX provides the following two functions for writing a generic AMReX plotfile. Many AMReX application codes may have their own plotfile routines that store additional information such as compiler options, git hashes of the source codes and ParmParse runtime parameters.

void WriteSingleLevelPlotfile (const std::string &plotfilename,
                               const MultiFab &mf,
                               const Vector<std::string> &varnames,
                               const Geometry &geom,
                               Real time,
                               int level_step);

void WriteMultiLevelPlotfile (const std::string &plotfilename,
                              int nlevels,
                              const Vector<const MultiFab*> &mf,
                              const Vector<std::string> &varnames,
                              const Vector<Geometry> &geom,
                              Real time,
                              const Vector<int> &level_steps,
                              const Vector<IntVect> &ref_ratio);

WriteSingleLevelPlotfile is for single level runs and WriteMultiLevelPlotfile is for multiple levels. The name of the plotfile is specified by the plotfilename argument. This is the top level directory name for the plotfile. In AMReX convention, the plotfile name consist of letters followed by numbers (e.g., plt00258). amrex::Concatenate is a useful helper function for making such strings.

int istep = 258;
const std::string& pfname = amrex::Concatenate("plt",istep); // plt00258

// By default there are 5 digits, but we can change it to say 4.
const std::string& pfname2 = amrex::Concatenate("plt",istep,4); // plt0258

istep =1234567;  // Having more than 5 digits is OK.
const std::string& pfname3 = amrex::Concatenate("plt",istep); // plt1234567

The argument mf above (MultiFab for single level and Vector<const MultiFab*> for multi-level) is the data to be written to the disk. Note that many visualization tools expect this to be cell-centered data. So for nodal data, we need to convert them to cell-centered data through some kind of averaging. Also note that if you have data at each AMR level in several MultiFabs, you need to build a new MultiFab at each level to hold all the data on that level. This involves local data copy in memory and is not expected to significantly increase the total wall time for writing plotfiles. For the multi-level version, the function expects Vector<const MultiFab*>, whereas the multi-level data are often stored as Vector<std::unique_ptr<MultiFab>>. AMReX has a helper function for this and one can use it as follows,

WriteMultiLevelPlotfile(......, amrex::GetVecOfConstPtrs(mf), ......);

The argument varnames has the names for each component of the MultiFab data. The size of the Vector should be equal to the number of components. The argument geom is for passing Geometry objects that contain the physical domain information. The argument time is for the time associated with the data. The argument level_step is for the current time step associated with the data. For multi-level plotfiles, the argument nlevels is the total number of levels, and we also need to provide the refinement ratio via an Vector of size nlevels-1.

We note that AMReX does not overwrite old plotfiles if the new plotfile has the same name. The old plotfiles will be renamed to new directories named like plt00350.old.46576787980.

Async Output

AMReX provides the ability to print MultiFabs, plotfiles and particle data asynchronously. Asynchronous output works by creating a copy of the data at the time of the call, which is written to disk by a persistent thread created during AMReX’s initialization. This allows the calculation to continue immediately, which can drastically reduce walltime spent writing to disk.

If the number of output files is less than the number of MPI ranks, AMReX’s async output requires MPI to be initialized with THREAD_MULTIPLE support. THREAD_MULTIPLE support allows multiple unique threads to run unique MPI calls simultaneously. This support is required to allow AMReX applications to perform MPI work while the Async Output concurrently pings ranks to signal that they can safely begin writing to their assigned files. However, THREAD_MULTIPLE can introduce additional overhead as each threads’ MPI operations must be scheduled safely around each other. Therefore, AMReX uses a lower level of support, SERIALIZED, by default and applications have to turn on THREAD_MULTIPLE support.

To turn on Async Output, use the input flag amrex.async_out=1. The number of output files can also be set, using amrex.async_out_nfiles. The default number of files is 64. If the number of ranks is larger than the number of files, THREAD_MULTIPLE must be turned on by adding MPI_THREAD_MULTIPLE=TRUE to the GNUMakefile. Otherwise, AMReX will throw an error.

Async Output works for a wide range of AMReX calls, including:

amrex::WriteSingleLevelPlotfile()
amrex::WriteMultiLevelPlotfile()
amrex::WriteMLMF()
VisMF::AsyncWrite()
ParticleContainer::Checkpoint()
ParticleContainer::WritePlotFile()
Amr::writePlotFile()
Amr::writeSmallPlotFile()
Amr::checkpoint()
AmrLevel::writePlotFile()
StateData::checkPoint()
FabSet::write()

Be aware: when using Async Output, a thread is spawned and exclusively used to perform output throughout the runtime. As such, you may oversubscribe resources if you launch an AMReX application that assigns all available hardware threads in another way, such as OpenMP. If you see any degradation when using Async Output and OpenMP, try using one less thread in OMP_NUM_THREADS to prevent oversubscription and get more consistent results.

HDF5 Plotfile

Besides AMReX’s native plotfile, applications can also write plotfile in the HDF5 format, which is a cross-platform, self-describing file format. The HDF5 plotfiles store the same information as the native format, and has the additional compression capability that can reduce the file size. Currently supported compression libraries include SZ and ZFP.

To enable HDF5 output, AMReX must be compiled and linked to an HDF5 library with parallel I/O support, by adding USE_HDF5=TRUE and HDF5_HOME=/path/to/hdf5/install/dir to the GNUMakefile. Many HPC systems have an HDF5 module available that can be loaded with module load hdf5 or module load cray-hdf5-parallel. To download and compile HDF5 from source code, please go to HDF5 Download webpage and follow the instructions (latest version is recommended and remember to turn on parallel I/O).

Following are two functions for writing a generic AMReX plotfile in HDF5 format, which are very similar to the AMReX native write functions.

void WriteSingleLevelPlotfileHDF5 (const std::string &plotfilename,
                                   const MultiFab &mf,
                                   const Vector<std::string> &varnames,
                                   const Geometry &geom,
                                   Real t,
                                   int level_step,
                                   const std::string &compression);

void WriteMultiLevelPlotfileHDF5 (const std::string &plotfilename,
                                  int nlevels,
                                  const Vector<const MultiFab*> &mf,
                                  const Vector<std::string> &varnames,
                                  const Vector<Geometry> &geom,
                                  Real time,
                                  const Vector<int> &level_steps,
                                  const Vector<IntVect> &ref_ratio,
                                  const std::string &compression);

WriteSingleLevelPlotfileHDF5 is for single level runs and WriteMultiLevelPlotfileHDF5 is for multiple levels. Their arguments are the same as the native ones except the last one, which optional, and specifies the compression parameters. These two functions write plotfiles with a Chombo-compatible HDF5 file schema, which can be read by visualization tools such as VisIt and ParaView using their built-in Chombo reader plugin (see the chapter on Visualization)

HDF5 Plotfile Compression

To enable SZ or ZFP data compression on the HDF5 datasets, the corresponding compression library and its HDF5 plugin must be available. To compile SZ or ZFP plugin, please refer to their documentation: H5Z-SZ and H5Z-ZFP, and adding USE_HDF5_SZ=TRUE, SZ_HOME=, or USE_HDF5_ZFP=TRUE, ZFP_HOME=, H5Z_HOME= to the GNUMakefile.

ZLIB compression is available without external libraries or other make flags. Different compression levels (at the cost of read/write time) can be used, just like GZIP.

The string argument compression in the above two functions controls whether to enable data compression and its parameters. Currently supported options include:

No compression
- None@0
ZLIB compression
- ZLIB@compression_level
SZ compression
- SZ@/path/to/sz.config
ZFP compression
- ZFP_RATE@rate
- ZFP_PRECISION@precision
- ZFP_ACCURACY@accuracy
- ZFP_REVERSIBLE@reversible

Using compression requires data to be stored in a chunked format. The size of these chunks can (and generally should) be configured by changing the HDF5_CHUNK_SIZE environment variable, with a default value of 1024 elements provided.

HDF5 Asynchronous Output

The HDF5 output also comes with its own asynchronous I/O support, which is different from the native async output mentioned in the previous section. To use the HDF5 asynchronous I/O VOL connector, download and compile by following the instructions at vol-async.

Since the HDF5 asynchronous I/O in AMReX does not use double buffering, vol-async must be compiled with -DENABLE_WRITE_MEMCPY=1 added to CFLAGS. When compiling AMReX, add USE_HDF5_ASYNC = TRUE, ABT_HOME=, ASYNC_HOME=, and MPI_THREAD_MULTIPLE=TRUE to the GNUMakefile. Refer to amrex/Tests/HDF5Benchmark/GNUmakefile for the example usage.

Alternative HDF5 Plotfile Schema

WriteSingleLevelPlotfileHDF5 and WriteMultiLevelPlotfileHDF5 write HDF5 plotfiles that store all the data on an AMR level as one 1D HDF5 dataset. Each AMR box’s data is linearized and the data of different variables are concatenated, resulting in an interleaved pattern for each variable. This could be undesirable when compression is used, as it may lead to applying the compression algorithm to multiple variables with different value ranges and characteristics, and reduce the compression ratio. To overcome this issue, two additional functions are provided to write each variable into individual HDF5 datasets: WriteSingleLevelPlotfileHDF5MultiDset and WriteMultiLevelPlotfileHDF5MultiDset. They use the exact same arguments as WriteSingleLevelPlotfileHDF5 and WriteMultiLevelPlotfileHDF5. However, this alternative schema is not yet supported by the visualization tools.

Checkpoint File

Checkpoint files are used for restarting simulations from where the checkpoints are written. Each application code has its own set of data needed for restart. AMReX provides I/O functions for basic data structures like MultiFab and BoxArray. These functions can be used to build codes for reading and writing checkpoint files. Since each application code has its own requirement, there is no standard AMReX checkpoint format. However we have provided an example restart capability in the tutorial Advection AmrCore. Refer to the functions ReadCheckpointFile() and WriteCheckpointFile() in this tutorial.

A checkpoint file is actually a directory with name, e.g., chk00010 containing a Header (text) file, along with subdirectories Level_0, Level_1, etc. containing the MultiFab data at each level of refinement. The Header file contains problem-specific data (such as the finest level, simulation time, time step, etc.), along with a printout of the BoxArray at each level of refinement.

When starting a simulation from a checkpoint file, a typical sequence in the code could be:

Read in the Header file data (except for the BoxArray data).
For each level of refinement, do the following in order:

– Read in the BoxArray

– Build a DistributionMapping

– Define any MultiFab, FluxRegister, etc. objects that are built upon the BoxArray and the DistributionMapping

– Read in the MultiFab data

We do this one level at a time because when you create a distribution map, it checks how much allocated MultiFab data already exists before assigning grids to processors.

Typically a checkpoint file is a directory containing some text files and sub-directories (e.g., Level_0 and Level_1) containing various data. It is a good idea that we fist make these directories ready for subsequently writing to the disk. For example, to build directories chk00010, chk00010/Level_0, and chk00010/Level_1, you could write:

const std::string& checkpointname = amrex::Concatenate("chk",10);

amrex::Print() << "Writing checkpoint " << checkpointname << "\n";

const int nlevels = 2;

bool callBarrier = true;

// ---- prebuild a hierarchy of directories
// ---- dirName is built first.  if dirName exists, it is renamed.  then build
// ---- dirName/subDirPrefix_0 .. dirName/subDirPrefix_nlevels-1
// ---- if callBarrier is true, call ParallelDescriptor::Barrier()
// ---- after all directories are built
// ---- ParallelDescriptor::IOProcessor() creates the directories
amrex::PreBuildDirectorHierarchy(checkpointname, "Level_", nlevels, callBarrier);

A checkpoint file of AMReX application codes often has a clear text Header file that only the I/O process writes to it using std::ofstream. The Header file contains problem-dependent information such as the time, the physical domain size, grids, etc. that are necessary for restarting the simulation. To guarantee that precision is not lost for storing floating point number like time in clear text file, the file stream’s precision needs to be set properly. And a stream buffer can also be used. For example,

// write Header file
if (ParallelDescriptor::IOProcessor()) {

    VisMF::IO_Buffer io_buffer(VisMF::IO_Buffer_Size);
    std::ofstream HeaderFile;
    HeaderFile.rdbuf()->pubsetbuf(io_buffer.dataPtr(), io_buffer.size());
    std::string HeaderFileName(checkpointname + "/Header");
    HeaderFile.open(HeaderFileName.c_str(), std::ofstream::out   |
                                            std::ofstream::trunc |
                                            std::ofstream::binary);

    if( ! HeaderFile.good()) {
        amrex::FileOpenFailed(HeaderFileName);
    }

    HeaderFile.precision(17);

    // write out title line
    HeaderFile << "Checkpoint file for AmrCoreAdv\n";

    // write out finest_level
    HeaderFile << finest_level << "\n";

    // write out array of istep
    for (int i = 0; i < istep.size(); ++i) {
        HeaderFile << istep[i] << " ";
    }
    HeaderFile << "\n";

    // write out array of dt
    for (int i = 0; i < dt.size(); ++i) {
        HeaderFile << dt[i] << " ";
    }
    HeaderFile << "\n";

    // write out array of t_new
    for (int i = 0; i < t_new.size(); ++i) {
        HeaderFile << t_new[i] << " ";
    }
    HeaderFile << "\n";

    // write the BoxArray at each level
    for (int lev = 0; lev <= finest_level; ++lev) {
        boxArray(lev).writeOn(HeaderFile);
        HeaderFile << '\n';
    }
}

amrex::VisMF is a class that can be used to perform MultiFab I/O in parallel. How many processes are allowed to perform I/O simultaneously can be set via

VisMF::SetNOutFiles(64);  // up to 64 processes, which is also the default.

The optimal number is of course system dependent. The following code shows how to write a MultiFab.

// write the MultiFab data to, e.g., chk00010/Level_0/
for (int lev = 0; lev <= finest_level; ++lev) {
    VisMF::Write(phi_new[lev],
                 amrex::MultiFabFileFullPrefix(lev, checkpointname, "Level_", "phi"));
}

It should also be noted that all the data including those in ghost cells are written/read by VisMF::Write/Read.

For reading the Header file, AMReX can have the I/O process read the file from the disk and broadcast it to others as Vector<char>. Then all processes can read the information with std::istringstream. For example,

std::string File(restart_chkfile + "/Header");

VisMF::IO_Buffer io_buffer(VisMF::GetIOBufferSize());

Vector<char> fileCharPtr;
ParallelDescriptor::ReadAndBcastFile(File, fileCharPtr);
std::string fileCharPtrString(fileCharPtr.dataPtr());
std::istringstream is(fileCharPtrString, std::istringstream::in);

std::string line, word;

// read in title line
std::getline(is, line);

// read in finest_level
is >> finest_level;
GotoNextLine(is);

// read in array of istep
std::getline(is, line);
{
    std::istringstream lis(line);
    int i = 0;
    while (lis >> word) {
        istep[i++] = std::stoi(word);
    }
}

// read in array of dt
std::getline(is, line);
{
    std::istringstream lis(line);
    int i = 0;
    while (lis >> word) {
        dt[i++] = std::stod(word);
    }
}

// read in array of t_new
std::getline(is, line);
{
    std::istringstream lis(line);
    int i = 0;
    while (lis >> word) {
        t_new[i++] = std::stod(word);
    }
}

The following code how to read in a BoxArray, create a DistributionMapping, build MultiFab and FluxRegister data, and read in a MultiFab from a checkpoint file, on a level-by-level basis:

for (int lev = 0; lev <= finest_level; ++lev) {

    // read in level 'lev' BoxArray from Header
    BoxArray ba;
    ba.readFrom(is);
    GotoNextLine(is);

    // create a distribution mapping
    DistributionMapping dm { ba, ParallelDescriptor::NProcs() };

    // set BoxArray grids and DistributionMapping dmap in AMReX_AmrMesh.H class
    SetBoxArray(lev, ba);
    SetDistributionMap(lev, dm);

    // build MultiFab and FluxRegister data
    int ncomp = 1;
    int nghost = 0;
    phi_old[lev].define(grids[lev], dmap[lev], ncomp, nghost);
    phi_new[lev].define(grids[lev], dmap[lev], ncomp, nghost);
    if (lev > 0 && do_reflux) {
        flux_reg[lev] = std::make_unique<FluxRegister>(grids[lev], dmap[lev], refRatio(lev-1), lev, ncomp);
    }
}

// read in the MultiFab data
for (int lev = 0; lev <= finest_level; ++lev) {
    VisMF::Read(phi_new[lev],
                amrex::MultiFabFileFullPrefix(lev, restart_chkfile, "Level_", "phi"));
}

It should be emphasized that calling VisMF::Read with an empty MultiFab (i.e., no memory allocated for floating point data) will result in a MultiFab with a new DistributionMapping that could be different from any other existing DistributionMapping objects and is not recommended.