Skip to content

Commit 4f83a7e

Browse files
ADIOS2: Override access mode (#1638)
* Main implementation * Add test * Documentation * Guard against ADIOS2 v2.7 * Extend the test to parallel
1 parent 74fdc47 commit 4f83a7e

File tree

5 files changed

+214
-0
lines changed

5 files changed

+214
-0
lines changed

docs/source/details/backendconfig.rst

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -122,6 +122,9 @@ Explanation of the single keys:
122122

123123
* ``adios2.engine.type``: A string that is passed directly to ``adios2::IO:::SetEngine`` for choosing the ADIOS2 engine to be used.
124124
Please refer to the `official ADIOS2 documentation <https://adios2.readthedocs.io/en/latest/engines/engines.html>`_ for a list of available engines.
125+
* ``adios2.engine.access_mode``: One of ``"Write", "Read", "Append", "ReadRandomAccess"``.
126+
Only needed in specific use cases, the access mode is usually determined from the specified ``openPMD::Access``.
127+
Useful for finetuning the backend-specific behavior of ADIOS2 when overwriting existing Iterations in file-based Append mode.
125128
* ``adios2.engine.parameters``: An associative array of string-formatted engine parameters, passed directly through to ``adios2::IO::SetParameters``.
126129
Please refer to the `official ADIOS2 documentation <https://adios2.readthedocs.io/en/latest/engines/engines.html>`_ for the available engine parameters.
127130
The openPMD-api does not interpret these values and instead simply forwards them to ADIOS2.

docs/source/usage/workflow.rst

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -102,6 +102,8 @@ The openPMD-api distinguishes between a number of different access modes:
102102
We suggest to fully define iterations when using Append mode (i.e. as if using Create mode) to avoid implementation-specific behavior.
103103
Appending to an openPMD Series is only supported on a per-iteration level.
104104

105+
**Tip:** Use the ``adios2.engine.access_mode`` :ref:`backend key <backendconfig>` of the :ref:`ADIOS2 backend <backends-adios2>` to finetune the backend-specific behavior of Append mode for niche use cases.
106+
105107
**Warning:** There is no reading involved in using Append mode.
106108
It is a user's responsibility to ensure that the appended dataset and the appended-to dataset are compatible with each other.
107109
The results of using incompatible backend configurations are undefined.

src/IO/ADIOS/ADIOS2IOHandler.cpp

Lines changed: 42 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -30,6 +30,7 @@
3030
#include "openPMD/IterationEncoding.hpp"
3131
#include "openPMD/auxiliary/Environment.hpp"
3232
#include "openPMD/auxiliary/Filesystem.hpp"
33+
#include "openPMD/auxiliary/JSON_internal.hpp"
3334
#include "openPMD/auxiliary/Mpi.hpp"
3435
#include "openPMD/auxiliary/StringManip.hpp"
3536
#include "openPMD/auxiliary/TypeTraits.hpp"
@@ -40,6 +41,7 @@
4041
#include <iterator>
4142
#include <memory>
4243
#include <set>
44+
#include <sstream>
4345
#include <stdexcept>
4446
#include <string>
4547
#include <type_traits>
@@ -1496,6 +1498,46 @@ void ADIOS2IOHandlerImpl::touch(
14961498

14971499
adios2::Mode ADIOS2IOHandlerImpl::adios2AccessMode(std::string const &fullPath)
14981500
{
1501+
if (m_config.json().contains("engine") &&
1502+
m_config["engine"].json().contains("access_mode"))
1503+
{
1504+
auto const &access_mode_json = m_config["engine"]["access_mode"].json();
1505+
auto maybe_access_mode_string =
1506+
json::asLowerCaseStringDynamic(access_mode_json);
1507+
if (!maybe_access_mode_string.has_value())
1508+
{
1509+
throw error::BackendConfigSchema(
1510+
{"adios2", "engine", "access_mode"}, "Must be of string type.");
1511+
}
1512+
auto access_mode_string = *maybe_access_mode_string;
1513+
using pair_t = std::pair<char const *, adios2::Mode>;
1514+
constexpr std::array<pair_t, 4> modeNames{
1515+
pair_t{"write", adios2::Mode::Write},
1516+
pair_t{"read", adios2::Mode::Read},
1517+
pair_t{"append", adios2::Mode::Append}
1518+
#if openPMD_HAS_ADIOS_2_8
1519+
,
1520+
pair_t{"readrandomaccess", adios2::Mode::ReadRandomAccess}
1521+
#endif
1522+
};
1523+
for (auto const &[name, mode] : modeNames)
1524+
{
1525+
if (name == access_mode_string)
1526+
{
1527+
return mode;
1528+
}
1529+
}
1530+
std::stringstream error;
1531+
error << "Unsupported value '" << access_mode_string
1532+
<< "'. Must be one of:";
1533+
for (auto const &pair : modeNames)
1534+
{
1535+
error << " '" << pair.first << "'";
1536+
}
1537+
error << '.';
1538+
throw error::BackendConfigSchema(
1539+
{"adios2", "engine", "access_mode"}, error.str());
1540+
}
14991541
switch (m_handler->m_backendAccess)
15001542
{
15011543
case Access::CREATE:

test/ParallelIOTest.cpp

Lines changed: 115 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -2082,4 +2082,119 @@ TEST_CASE("joined_dim", "[parallel]")
20822082
}
20832083
}
20842084

2085+
#if openPMD_HAVE_ADIOS2_BP5
2086+
// Parallel version of the same test from SerialIOTest.cpp
2087+
TEST_CASE("adios2_flush_via_step")
2088+
{
2089+
int size_i(0), rank_i(0);
2090+
MPI_Comm_rank(MPI_COMM_WORLD, &rank_i);
2091+
MPI_Comm_size(MPI_COMM_WORLD, &size_i);
2092+
Extent::value_type const size(size_i), rank(rank_i);
2093+
2094+
Series write(
2095+
"../samples/adios2_flush_via_step_parallel/simData_%T.bp5",
2096+
Access::CREATE,
2097+
MPI_COMM_WORLD,
2098+
R"(adios2.engine.parameters.FlattenSteps = "on")");
2099+
std::vector<float> data(10);
2100+
for (Iteration::IterationIndex_t i = 0; i < 5; ++i)
2101+
{
2102+
Iteration it = write.writeIterations()[i];
2103+
auto E_x = it.meshes["E"]["x"];
2104+
E_x.resetDataset({Datatype::FLOAT, {size, 10, 10}});
2105+
for (Extent::value_type j = 0; j < 10; ++j)
2106+
{
2107+
std::iota(
2108+
data.begin(), data.end(), i * 100 * size + rank * 100 + j * 10);
2109+
E_x.storeChunk(data, {rank, j, 0}, {1, 1, 10});
2110+
write.flush(R"(adios2.engine.preferred_flush_target = "new_step")");
2111+
}
2112+
it.close();
2113+
}
2114+
2115+
#if openPMD_HAS_ADIOS_2_10_1
2116+
for (auto access : {Access::READ_RANDOM_ACCESS, Access::READ_LINEAR})
2117+
{
2118+
Series read(
2119+
"../samples/adios2_flush_via_step_parallel/simData_%T.%E",
2120+
access,
2121+
MPI_COMM_WORLD);
2122+
std::vector<float> load_data(100 * size);
2123+
data.resize(100 * size);
2124+
for (auto iteration : read.readIterations())
2125+
{
2126+
std::iota(
2127+
data.begin(),
2128+
data.end(),
2129+
iteration.iterationIndex * size * 100);
2130+
iteration.meshes["E"]["x"].loadChunkRaw(
2131+
load_data.data(), {0, 0, 0}, {size, 10, 10});
2132+
iteration.close();
2133+
REQUIRE(load_data == data);
2134+
}
2135+
}
2136+
#endif
2137+
2138+
/*
2139+
* Now emulate restarting from a checkpoint after a crash and continuing to
2140+
* write to the output Series. The semantics of openPMD::Access::APPEND
2141+
* don't fully fit here since that mode is for adding new Iterations to an
2142+
* existing Series. What we truly want to do is to continue writing to an
2143+
* Iteration without replacing it with a new one. So we must use the option
2144+
* adios2.engine.access_mode = "append" to tell the ADIOS2 backend that new
2145+
* steps should be added to an existing Iteration file.
2146+
*/
2147+
2148+
write = Series(
2149+
"../samples/adios2_flush_via_step_parallel/simData_%T.bp5",
2150+
Access::APPEND,
2151+
MPI_COMM_WORLD,
2152+
R"(
2153+
[adios2.engine]
2154+
access_mode = "append"
2155+
parameters.FlattenSteps = "on"
2156+
)");
2157+
for (Iteration::IterationIndex_t i = 0; i < 5; ++i)
2158+
{
2159+
Iteration it = write.writeIterations()[i];
2160+
auto E_x = it.meshes["E"]["y"];
2161+
E_x.resetDataset({Datatype::FLOAT, {size, 10, 10}});
2162+
for (Extent::value_type j = 0; j < 10; ++j)
2163+
{
2164+
std::iota(
2165+
data.begin(), data.end(), i * 100 * size + rank * 100 + j * 10);
2166+
E_x.storeChunk(data, {rank, j, 0}, {1, 1, 10});
2167+
write.flush(R"(adios2.engine.preferred_flush_target = "new_step")");
2168+
}
2169+
it.close();
2170+
}
2171+
2172+
#if openPMD_HAS_ADIOS_2_10_1
2173+
for (auto access : {Access::READ_RANDOM_ACCESS, Access::READ_LINEAR})
2174+
{
2175+
Series read(
2176+
"../samples/adios2_flush_via_step_parallel/simData_%T.%E",
2177+
access,
2178+
MPI_COMM_WORLD);
2179+
std::vector<float> load_data(100 * size);
2180+
data.resize(100 * size);
2181+
for (auto iteration : read.readIterations())
2182+
{
2183+
std::iota(
2184+
data.begin(),
2185+
data.end(),
2186+
iteration.iterationIndex * size * 100);
2187+
iteration.meshes["E"]["x"].loadChunkRaw(
2188+
load_data.data(), {0, 0, 0}, {size, 10, 10});
2189+
iteration.meshes["E"]["y"].loadChunkRaw(
2190+
load_data.data(), {0, 0, 0}, {size, 10, 10});
2191+
iteration.close();
2192+
REQUIRE(load_data == data);
2193+
REQUIRE(load_data == data);
2194+
}
2195+
}
2196+
#endif
2197+
}
2198+
#endif
2199+
20852200
#endif // openPMD_HAVE_ADIOS2 && openPMD_HAVE_MPI

test/SerialIOTest.cpp

Lines changed: 52 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -4475,6 +4475,58 @@ TEST_CASE("adios2_flush_via_step")
44754475
}
44764476
}
44774477
#endif
4478+
4479+
/*
4480+
* Now emulate restarting from a checkpoint after a crash and continuing to
4481+
* write to the output Series. The semantics of openPMD::Access::APPEND
4482+
* don't fully fit here since that mode is for adding new Iterations to an
4483+
* existing Series. What we truly want to do is to continue writing to an
4484+
* Iteration without replacing it with a new one. So we must use the option
4485+
* adios2.engine.access_mode = "append" to tell the ADIOS2 backend that new
4486+
* steps should be added to an existing Iteration file.
4487+
*/
4488+
4489+
write = Series(
4490+
"../samples/adios2_flush_via_step/simData_%T.bp5",
4491+
Access::APPEND,
4492+
R"(
4493+
[adios2.engine]
4494+
access_mode = "append"
4495+
parameters.FlattenSteps = "on"
4496+
)");
4497+
for (Iteration::IterationIndex_t i = 0; i < 5; ++i)
4498+
{
4499+
Iteration it = write.writeIterations()[i];
4500+
auto E_x = it.meshes["E"]["y"];
4501+
E_x.resetDataset({Datatype::FLOAT, {10, 10}});
4502+
for (Extent::value_type j = 0; j < 10; ++j)
4503+
{
4504+
std::iota(data.begin(), data.end(), i * 100 + j * 10);
4505+
E_x.storeChunk(data, {j, 0}, {1, 10});
4506+
write.flush(R"(adios2.engine.preferred_flush_target = "new_step")");
4507+
}
4508+
it.close();
4509+
}
4510+
4511+
#if openPMD_HAS_ADIOS_2_10_1
4512+
for (auto access : {Access::READ_RANDOM_ACCESS, Access::READ_LINEAR})
4513+
{
4514+
Series read("../samples/adios2_flush_via_step/simData_%T.%E", access);
4515+
std::vector<float> load_data(100);
4516+
data.resize(100);
4517+
for (auto iteration : read.readIterations())
4518+
{
4519+
std::iota(data.begin(), data.end(), iteration.iterationIndex * 100);
4520+
iteration.meshes["E"]["x"].loadChunkRaw(
4521+
load_data.data(), {0, 0}, {10, 10});
4522+
iteration.meshes["E"]["y"].loadChunkRaw(
4523+
load_data.data(), {0, 0}, {10, 10});
4524+
iteration.close();
4525+
REQUIRE(load_data == data);
4526+
REQUIRE(load_data == data);
4527+
}
4528+
}
4529+
#endif
44784530
}
44794531
#endif
44804532

0 commit comments

Comments
 (0)