Steps
In the previous tutorial, we introduced the concept of operators, and briefly touched upon the concept of steps.
In this tutorial, we will explore how to write data for multiple steps, and how to read them back.
So let’s dig in!
Start editing the skeleton file ADIOS2/examples/hello/bpStepsWriteRead/bpStepsWriteRead_tutorialSkeleton.cpp.
In an MPI application first we need to always initialize MPI. We do that with the following lines:
int rank, size;
int rank, size;
int provided;
// MPI_THREAD_MULTIPLE is only required if you enable the SST MPI_DP
MPI_Init_thread(&argc, &argv, MPI_THREAD_MULTIPLE, &provided);
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
MPI_Comm_size(MPI_COMM_WORLD, &size);
This application has an optional command line argument for engine being used. If no argument is provided, the default engine is BPFile.
const std::string engine = argv[1] ? argv[1] : "BPFile";
We will define the number of steps and the size of the data that we will create.
const std::string filename = engine + "StepsWriteRead.bp";
const unsigned int nSteps = 10;
const unsigned int Nx = 60000;
Now we need to create an ADIOS2 instance.
adios2::ADIOS adios(MPI_COMM_WORLD);
Now we will populate the writer function with the following signature:
void writer(adios2::ADIOS &adios, const std::string &engine, const std::string &fname,
const size_t Nx, unsigned int nSteps, int rank, int size)
{
...
}
Let’s create some simulation data. We will create a 1D array of size Nx, and fill it with 0.block
std::vector<double> simData(Nx, 0.0);
Now we will create an IO object and set the engine type.block
adios2::IO bpIO = adios.DeclareIO("SimulationOutput");
io.SetEngine(engine);
Note
The beauty of ADIOS2 is that you write the same code for all engines. The only thing that changes is the engine name. The underlying engine handles all the intricacies of the engine’s format, and the user enjoys the API’s simplicity.
Now we will create a variable for the simulation data and the step.
const adios2::Dims shape{static_cast<size_t>(size * Nx)};
const adios2::Dims start{static_cast<size_t>(rank * Nx)};
const adios2::Dims count{Nx};
auto bpFloats = bpIO.DefineVariable<float>("bpFloats", shape, start, count);
auto bpStep = bpIO.DefineVariable<unsigned int>("bpStep");
Now we will open the file for writing.
adios2::Engine bpWriter = bpIO.Open(fname, adios2::Mode::Write);
Now we will write the data for each step.
for (unsigned int step = 0; step < nSteps; ++step)
{
const adios2::Box<adios2::Dims> sel({0}, {Nx});
bpFloats.SetSelection(sel);
bpWriter.BeginStep();
bpWriter.Put(bpFloats, simData.data());
bpWriter.Put(bpStep, step);
bpWriter.EndStep();
// Update values in the simulation data
update_array(simData, 10);
}
Now we will close the file.
bpWriter.Close();
Now we will populate the reader function with the following signature:
void reader(adios2::ADIOS &adios, const std::string &engine, const std::string &fname,
const size_t Nx, unsigned int /*nSteps*/, int rank, int /*size*/)
{
...
}
Now we will create an IO object and set the engine type.
adios2::IO bpIO = adios.DeclareIO("SimulationOutput");
io.SetEngine(engine);
Now we will open the file for reading.
adios2::Engine bpReader = bpIO.Open(fname, adios2::Mode::Read);
Now we will create a vector to store simData and a variable for the step.
std::vector<float> simData(Nx, 0);
unsigned int inStep = 0;
Now we will read the data for each step.
for (unsigned int step = 0; bpReader.BeginStep() == adios2::StepStatus::OK; ++step)
{
auto bpFloats = bpIO.InquireVariable<float>("bpFloats");
if (bpFloats)
{
const adios2::Box<adios2::Dims> sel({{Nx * rank}, {Nx}});
bpFloats.SetSelection(sel);
bpReader.Get(bpFloats, simData.data());
}
auto bpStep = bpIO.InquireVariable<unsigned int>("bpStep");
if (bpStep)
{
bpReader.Get(bpStep, &inStep);
}
bpReader.EndStep();
}
Now we will close the file.
bpReader.Close();
Now we will call the writer and reader functions:
writer(adios, engine, filename, Nx, nSteps, rank, size);
reader(adios, engine, filename, Nx, nSteps, rank, size);
Finally we need to finalize MPI.
MPI_Finalize();
The final code should look as follows (excluding try/catch and optional usage of MPI), and it was derived from the example ADIOS2/examples/hello/bpStepsWriteRead/bpStepsWriteRead.cpp.
/*
* Distributed under the OSI-approved Apache License, Version 2.0. See
* accompanying file Copyright.txt for details.
*
* bpStepsWriteRead.cpp Simple example of writing and reading data through ADIOS2 BP engine with
* multiple simulations steps for every IO step.
*
* Created on: Feb 16, 2017
* Author: William F Godoy godoywf@ornl.gov
*/
#include <algorithm> // std::for_each
#include <ios> // std::ios_base::failure
#include <iostream> // std::cout
#if ADIOS2_USE_MPI
#include <mpi.h>
#endif
#include <stdexcept> //std::invalid_argument std::exception
#include <vector>
#include <adios2.h>
void update_array(std::vector<float> &array, int val)
{
std::transform(array.begin(), array.end(), array.begin(),
[val](float v) -> float { return v + static_cast<float>(val); });
}
void writer(adios2::ADIOS &adios, const std::string &engine, const std::string &fname,
const size_t Nx, unsigned int nSteps, int rank, int size)
{
std::vector<float> simData(Nx, 0);
adios2::IO bpIO = adios.DeclareIO("WriteIO");
bpIO.SetEngine(engine);
const adios2::Dims shape{static_cast<size_t>(size * Nx)};
const adios2::Dims start{static_cast<size_t>(rank * Nx)};
const adios2::Dims count{Nx};
auto bpFloats = bpIO.DefineVariable<float>("bpFloats", shape, start, count);
auto bpStep = bpIO.DefineVariable<unsigned int>("bpStep");
adios2::Engine bpWriter = bpIO.Open(fname, adios2::Mode::Write);
for (unsigned int step = 0; step < nSteps; ++step)
{
const adios2::Box<adios2::Dims> sel({0}, {Nx});
bpFloats.SetSelection(sel);
bpWriter.BeginStep();
bpWriter.Put(bpFloats, simData.data());
bpWriter.Put(bpStep, step);
bpWriter.EndStep();
// Update values in the simulation data
update_array(simData, 10);
}
bpWriter.Close();
}
void reader(adios2::ADIOS &adios, const std::string &engine, const std::string &fname,
const size_t Nx, unsigned int /*nSteps*/, int rank, int /*size*/)
{
adios2::IO bpIO = adios.DeclareIO("ReadIO");
bpIO.SetEngine(engine);
adios2::Engine bpReader = bpIO.Open(fname, adios2::Mode::Read);
std::vector<float> simData(Nx, 0);
unsigned int inStep = 0;
for (unsigned int step = 0; bpReader.BeginStep() == adios2::StepStatus::OK; ++step)
{
auto bpFloats = bpIO.InquireVariable<float>("bpFloats");
if (bpFloats)
{
const adios2::Box<adios2::Dims> sel({{Nx * rank}, {Nx}});
bpFloats.SetSelection(sel);
bpReader.Get(bpFloats, simData.data());
}
auto bpStep = bpIO.InquireVariable<unsigned int>("bpStep");
if (bpStep)
{
bpReader.Get(bpStep, &inStep);
}
bpReader.EndStep();
if (inStep != step)
{
std::cout << "ERROR: step mismatch\n";
return;
}
}
bpReader.Close();
}
int main(int argc, char *argv[])
{
int rank, size;
#if ADIOS2_USE_MPI
int provided;
// MPI_THREAD_MULTIPLE is only required if you enable the SST MPI_DP
MPI_Init_thread(&argc, &argv, MPI_THREAD_MULTIPLE, &provided);
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
MPI_Comm_size(MPI_COMM_WORLD, &size);
#else
rank = 0;
size = 1;
#endif
const std::string engine = argv[1] ? argv[1] : "BPFile";
std::cout << "Using engine " << engine << std::endl;
const std::string filename = engine + "StepsWriteRead.bp";
const unsigned int nSteps = 10;
const unsigned int Nx = 60000;
try
{
/** ADIOS class factory of IO class objects */
#if ADIOS2_USE_MPI
adios2::ADIOS adios(MPI_COMM_WORLD);
#else
adios2::ADIOS adios;
#endif
writer(adios, engine, filename, Nx, nSteps, rank, size);
reader(adios, engine, filename, Nx, nSteps, rank, size);
}
catch (std::invalid_argument &e)
{
std::cout << "Invalid argument exception, STOPPING PROGRAM from rank " << rank << "\n";
std::cout << e.what() << "\n";
}
catch (std::ios_base::failure &e)
{
std::cout << "IO System base failure exception, STOPPING PROGRAM "
"from rank "
<< rank << "\n";
std::cout << e.what() << "\n";
}
catch (std::exception &e)
{
std::cout << "Exception, STOPPING PROGRAM from rank " << rank << "\n";
std::cout << e.what() << "\n";
}
#if ADIOS2_USE_MPI
MPI_Finalize();
#endif
return 0;
}
You can compile and run it as follows:
cd Path-To-ADIOS2/examples/hello/bpStepsWriteRead
mkdir build
cd build
cmake -DADIOS2_DIR=Path-To-ADIOS2/build/ ..
cmake --build .
mpirun -np 2 ./adios2_hello_bpStepsWriteRead_mpi
You can check the content of the output file “BPFileStepsWriteRead.bp” using bpls as follows:
Path-To-ADIOS2/build/bin/bpls ./BPFileStepsWriteRead.bp
float bpFloats 10*{120000}
uint32_t bpStep 10*scalar