24 #include "openPMD/config.hpp" 27 #include "RandomDatasetFiller.hpp" 29 #include "openPMD/openPMD.hpp" 30 #include "openPMD/DatatypeHelpers.hpp" 31 #include "openPMD/benchmark/mpi/MPIBenchmarkReport.hpp" 32 #include "openPMD/benchmark/mpi/DatasetFiller.hpp" 33 #include "openPMD/benchmark/mpi/BlockSlicer.hpp" 56 template<
typename DatasetFillerProv
ider >
61 using extentT = Extent::value_type;
62 MPI_Comm communicator = MPI_COMM_WORLD;
69 std::shared_ptr< BlockSlicer > m_blockSlicer;
71 DatasetFillerProvider m_dfp;
91 std::shared_ptr< BlockSlicer > blockSlicer,
92 DatasetFillerProvider dfp,
93 MPI_Comm comm = MPI_COMM_WORLD
107 std::string compression,
108 uint8_t compressionLevel,
111 typename decltype( Series::iterations )::key_type iterations,
127 std::string compression,
128 uint8_t compressionLevel,
131 typename decltype( Series::iterations)::key_type iterations
134 void resetConfigurations( );
145 template<
typename Clock >
151 std::string m_basePath;
159 typename decltype( Series::iterations)::key_type>>
183 template<
typename Clock >
184 struct BenchmarkExecution
190 m_benchmark { benchmark }
209 typename Clock::duration writeBenchmark(
210 std::string
const & compression,
214 std::string
const & extension,
216 typename decltype( Series::iterations)::key_type iterations
231 typename Clock::duration readBenchmark(
234 std::string extension,
235 typename decltype( Series::iterations)::key_type iterations
238 template<
typename T >
260 template<
typename DatasetFillerProv
ider >
261 template<
typename Clock >
268 BenchmarkExecution< Clock > exec {
this };
270 std::set< Datatype > datatypes;
271 for(
auto const & conf: m_configurations )
273 datatypes.insert( std::get< DTYPE >( conf ) );
275 for( Datatype dt: datatypes )
289 template<
typename DatasetFillerProv
ider >
291 std::string basePath,
293 std::shared_ptr< BlockSlicer > blockSlicer,
294 DatasetFillerProvider dfp,
297 communicator { comm },
298 totalExtent { std::move( tExtent ) },
299 m_blockSlicer { std::move( blockSlicer ) },
301 m_basePath { std::move( basePath ) }
303 if( m_blockSlicer ==
nullptr )
304 throw std::runtime_error(
"Argument blockSlicer cannot be a nullptr!");
308 template<
typename DatasetFillerProv
ider >
328 return m_blockSlicer->sliceBlock(
336 template<
typename DatasetFillerProv
ider >
338 std::string compression,
339 uint8_t compressionLevel,
342 typename decltype( Series::iterations)::key_type iterations,
346 this->m_configurations
358 template<
typename DatasetFillerProv
ider >
360 std::string compression,
361 uint8_t compressionLevel,
364 typename decltype( Series::iterations)::key_type iterations
383 template<
typename DatasetFillerProv
ider >
391 template<
typename DatasetFillerProv
ider >
392 template<
typename Clock >
393 template<
typename T >
394 typename Clock::duration
396 std::string
const & compression,
400 std::string
const & extension,
402 typename decltype( Series::iterations)::key_type iterations
405 MPI_Barrier( m_benchmark->communicator );
406 auto start = Clock::now( );
410 m_benchmark->m_basePath +
"." + extension,
412 m_benchmark->communicator
415 for(
typename decltype( Series::iterations)::key_type i = 0;
419 auto writeData = datasetFiller->produceData( );
424 series.iterations[i].meshes[
"id"][MeshRecordComponent::SCALAR];
426 Datatype datatype = determineDatatype( writeData );
431 if( !compression.empty( ) )
433 dataset.setCompression(
439 id.resetDataset( dataset );
451 MPI_Barrier( m_benchmark->communicator );
452 auto end = Clock::now( );
455 for(
typename decltype( Series::iterations)::key_type i = 0;
459 datasetFiller->produceData( );
461 auto deduct = Clock::now( );
463 return end - start - ( deduct - end );
467 template<
typename DatasetFillerProv
ider >
468 template<
typename Clock >
469 template<
typename T >
470 typename Clock::duration
474 std::string extension,
475 typename decltype( Series::iterations)::key_type iterations
478 MPI_Barrier( m_benchmark->communicator );
480 auto start = Clock::now( );
483 m_benchmark->m_basePath +
"." + extension,
485 m_benchmark->communicator
488 for(
typename decltype( Series::iterations)::key_type i = 0;
494 series.iterations[i].meshes[
"id"][MeshRecordComponent::SCALAR];
497 auto chunk_data =
id.loadChunk<
T >(
504 MPI_Barrier( m_benchmark->communicator );
505 auto end = Clock::now( );
510 template<
typename DatasetFillerProv
ider >
511 template<
typename Clock >
512 template<
typename T >
519 Datatype dt = determineDatatype< T >( );
525 for(
auto const & config: m_benchmark->m_configurations )
527 std::string compression;
528 uint8_t compressionLevel;
532 typename decltype( Series::iterations)::key_type iterations;
547 auto localCuboid = m_benchmark->slice( size );
549 extentT blockSize = 1;
550 for(
auto ext: localCuboid.second )
554 dsf->setNumberOfItems( blockSize );
556 auto writeTime = writeBenchmark< T >(
565 auto readTime = readBenchmark< T >(
589 template<
typename DatasetFillerProv
ider >
590 template<
typename Clock >
598 throw std::runtime_error(
"Unknown/unsupported datatype requested to be benchmarked." );
void addReport(int rootThread, std::string compression, uint8_t level, std::string extension, int threadSize, Datatype dt, typename decltype(Series::iterations)::key_type iterations, std::pair< Duration, Duration > const &report)
Add results for a certain compression strategy and level.
Definition: MPIBenchmarkReport.hpp:245
Definition: Dataset.hpp:36
An abstract class to create one iteration of data per thread.
Definition: DatasetFiller.hpp:37
auto switchType(Datatype dt, Action action, Args &&... args) -> decltype(action. template operator()< char >(std::forward< Args >(args)...))
Generalizes switching over an openPMD datatype.
Definition: DatatypeHelpers.hpp:131
Extent totalExtent
Total extent of the hypercuboid used in the benchmark.
Definition: MPIBenchmark.hpp:67
create new series and truncate existing (files)
MPIBenchmarkReport< typename Clock::duration > runBenchmark(int rootThread=0)
Main function for running a benchmark.
Definition: MPIBenchmark.hpp:263
Datatype
Concrete datatype of an object available at runtime.
Definition: Datatype.hpp:45
Root level of the openPMD hierarchy.
Definition: Series.hpp:476
Class representing a benchmark.
Definition: MPIBenchmark.hpp:57
The report for a single benchmark produced by <openPMD/benchmark/mpi/MPIBenchmark>.
Definition: MPIBenchmarkReport.hpp:44
MPIBenchmark(std::string basePath, Extent tExtent, std::shared_ptr< BlockSlicer > blockSlicer, DatasetFillerProvider dfp, MPI_Comm comm=MPI_COMM_WORLD)
Construct an MPI benchmark manually.
Definition: MPIBenchmark.hpp:290
Public definitions of openPMD-api.
Definition: Date.cpp:29
Definition: MeshRecordComponent.hpp:30
void flush()
Execute all required remaining IO operations to write or read data.
Definition: Series.cpp:348
open series as read-only, fails if series is not found
void addConfiguration(std::string compression, uint8_t compressionLevel, std::string backend, Datatype dt, typename decltype(Series::iterations)::key_type iterations, int threadSize)
Definition: MPIBenchmark.hpp:337