openPMD-api
MPIBenchmark.hpp
1 /* Copyright 2018-2021 Franz Poeschel
2  *
3  * This file is part of openPMD-api.
4  *
5  * openPMD-api is free software: you can redistribute it and/or modify
6  * it under the terms of of either the GNU General Public License or
7  * the GNU Lesser General Public License as published by
8  * the Free Software Foundation, either version 3 of the License, or
9  * (at your option) any later version.
10  *
11  * openPMD-api is distributed in the hope that it will be useful,
12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  * GNU General Public License and the GNU Lesser General Public License
15  * for more details.
16  *
17  * You should have received a copy of the GNU General Public License
18  * and the GNU Lesser General Public License along with openPMD-api.
19  * If not, see <http://www.gnu.org/licenses/>.
20  */
21 
22 #pragma once
23 
24 #include "openPMD/config.hpp"
25 #if openPMD_HAVE_MPI
26 
27 #include "RandomDatasetFiller.hpp"
28 
29 #include "openPMD/DatatypeHelpers.hpp"
30 #include "openPMD/benchmark/mpi/BlockSlicer.hpp"
31 #include "openPMD/benchmark/mpi/DatasetFiller.hpp"
32 #include "openPMD/benchmark/mpi/MPIBenchmarkReport.hpp"
33 #include "openPMD/openPMD.hpp"
34 
35 #include <mpi.h>
36 
37 #include <chrono>
38 #include <exception>
39 #include <iostream>
40 #include <set>
41 #include <sstream>
42 #include <tuple>
43 #include <utility>
44 #include <vector>
45 
46 namespace openPMD
47 {
56 template <typename DatasetFillerProvider>
58 {
59 
60 public:
61  using extentT = Extent::value_type;
62  MPI_Comm communicator = MPI_COMM_WORLD;
63 
67  Extent totalExtent;
68 
69  std::shared_ptr<BlockSlicer> m_blockSlicer;
70 
71  DatasetFillerProvider m_dfp;
72 
88  std::string basePath,
89  Extent tExtent,
90  std::shared_ptr<BlockSlicer> blockSlicer,
91  DatasetFillerProvider dfp,
92  MPI_Comm comm = MPI_COMM_WORLD);
93 
106  void addConfiguration(
107  std::string compression,
108  uint8_t compressionLevel,
109  std::string backend,
110  Datatype dt,
111  typename decltype(Series::iterations)::key_type iterations,
112  int threadSize);
113 
127  void addConfiguration(
128  std::string compression,
129  uint8_t compressionLevel,
130  std::string backend,
131  Datatype dt,
132  typename decltype(Series::iterations)::key_type iterations);
133 
134  void resetConfigurations();
135 
144  template <typename Clock>
146  runBenchmark(int rootThread = 0);
147 
148 private:
149  std::string m_basePath;
150  std::vector<std::tuple<
151  std::string,
152  uint8_t,
153  std::string,
154  int,
155  Datatype,
156  typename decltype(Series::iterations)::key_type>>
157  m_configurations;
158 
159  enum Config
160  {
161  COMPRESSION = 0,
162  COMPRESSION_LEVEL,
163  BACKEND,
164  NRANKS,
165  DTYPE,
166  ITERATIONS
167  };
168 
169  std::pair<Offset, Extent> slice(int size);
170 
177  template <typename Clock>
178  struct BenchmarkExecution
179  {
181 
182  explicit BenchmarkExecution(
184  : m_benchmark{benchmark}
185  {}
186 
199  template <typename T>
200  typename Clock::duration writeBenchmark(
201  std::string const &compression,
202  uint8_t level,
203  Offset &offset,
204  Extent &extent,
205  std::string const &extension,
206  std::shared_ptr<DatasetFiller<T>> datasetFiller,
207  typename decltype(Series::iterations)::key_type iterations);
208 
218  template <typename T>
219  typename Clock::duration readBenchmark(
220  Offset &offset,
221  Extent &extent,
222  std::string extension,
223  typename decltype(Series::iterations)::key_type iterations);
224 
225  template <typename T>
226  void operator()(
228  int rootThread = 0);
229 
230  template <int n>
231  void operator()(MPIBenchmarkReport<typename Clock::duration> &, int);
232  };
233 };
234 
235 // Implementation
236 
237 template <typename DatasetFillerProvider>
238 template <typename Clock>
241 {
242  MPIBenchmarkReport<typename Clock::duration> res{this->communicator};
243  BenchmarkExecution<Clock> exec{this};
244 
245  std::set<Datatype> datatypes;
246  for (auto const &conf : m_configurations)
247  {
248  datatypes.insert(std::get<DTYPE>(conf));
249  }
250  for (Datatype dt : datatypes)
251  {
252  switchType(dt, exec, res, rootThread);
253  }
254 
255  return res;
256 }
257 
258 template <typename DatasetFillerProvider>
260  std::string basePath,
261  Extent tExtent,
262  std::shared_ptr<BlockSlicer> blockSlicer,
263  DatasetFillerProvider dfp,
264  MPI_Comm comm)
265  : communicator{comm}
266  , totalExtent{std::move(tExtent)}
267  , m_blockSlicer{std::move(blockSlicer)}
268  , m_dfp{dfp}
269  , m_basePath{std::move(basePath)}
270 {
271  if (m_blockSlicer == nullptr)
272  throw std::runtime_error("Argument blockSlicer cannot be a nullptr!");
273 }
274 
275 template <typename DatasetFillerProvider>
276 std::pair<Offset, Extent> MPIBenchmark<DatasetFillerProvider>::slice(int size)
277 {
278  int actualSize;
279  MPI_Comm_size(this->communicator, &actualSize);
280  int rank;
281  MPI_Comm_rank(this->communicator, &rank);
282  size = std::min(size, actualSize);
283  return m_blockSlicer->sliceBlock(totalExtent, size, rank);
284 }
285 
286 template <typename DatasetFillerProvider>
288  std::string compression,
289  uint8_t compressionLevel,
290  std::string backend,
291  Datatype dt,
292  typename decltype(Series::iterations)::key_type iterations,
293  int threadSize)
294 {
295  this->m_configurations.emplace_back(
296  compression, compressionLevel, backend, threadSize, dt, iterations);
297 }
298 
299 template <typename DatasetFillerProvider>
301  std::string compression,
302  uint8_t compressionLevel,
303  std::string backend,
304  Datatype dt,
305  typename decltype(Series::iterations)::key_type iterations)
306 {
307  int size;
308  MPI_Comm_size(communicator, &size);
310  compression, compressionLevel, backend, dt, iterations, size);
311 }
312 
313 template <typename DatasetFillerProvider>
315 {
316  this->m_compressions.clear();
317 }
318 
319 template <typename DatasetFillerProvider>
320 template <typename Clock>
321 template <typename T>
322 typename Clock::duration
324  std::string const &compression,
325  uint8_t level,
326  Offset &offset,
327  Extent &extent,
328  std::string const &extension,
329  std::shared_ptr<DatasetFiller<T>> datasetFiller,
330  typename decltype(Series::iterations)::key_type iterations)
331 {
332  MPI_Barrier(m_benchmark->communicator);
333  auto start = Clock::now();
334 
335  // open file for writing
336  Series series = Series(
337  m_benchmark->m_basePath + "." + extension,
339  m_benchmark->communicator);
340 
341  for (typename decltype(Series::iterations)::key_type i = 0; i < iterations;
342  i++)
343  {
344  auto writeData = datasetFiller->produceData();
345 
347  series.iterations[i].meshes["id"][MeshRecordComponent::SCALAR];
348 
349  Datatype datatype = determineDatatype(writeData);
350  Dataset dataset = Dataset(datatype, m_benchmark->totalExtent);
351  if (!compression.empty())
352  {
353  dataset.setCompression(compression, level);
354  }
355 
356  id.resetDataset(dataset);
357 
358  series.flush();
359 
360  id.storeChunk<T>(writeData, offset, extent);
361  series.flush();
362  }
363 
364  MPI_Barrier(m_benchmark->communicator);
365  auto end = Clock::now();
366 
367  // deduct the time needed for data generation
368  for (typename decltype(Series::iterations)::key_type i = 0; i < iterations;
369  i++)
370  {
371  datasetFiller->produceData();
372  }
373  auto deduct = Clock::now();
374 
375  return end - start - (deduct - end);
376 }
377 
378 template <typename DatasetFillerProvider>
379 template <typename Clock>
380 template <typename T>
381 typename Clock::duration
383  Offset &offset,
384  Extent &extent,
385  std::string extension,
386  typename decltype(Series::iterations)::key_type iterations)
387 {
388  MPI_Barrier(m_benchmark->communicator);
389  // let every thread measure time
390  auto start = Clock::now();
391 
392  Series series = Series(
393  m_benchmark->m_basePath + "." + extension,
395  m_benchmark->communicator);
396 
397  for (typename decltype(Series::iterations)::key_type i = 0; i < iterations;
398  i++)
399  {
401  series.iterations[i].meshes["id"][MeshRecordComponent::SCALAR];
402 
403  auto chunk_data = id.loadChunk<T>(offset, extent);
404  series.flush();
405  }
406 
407  MPI_Barrier(m_benchmark->communicator);
408  auto end = Clock::now();
409  return end - start;
410 }
411 
412 template <typename DatasetFillerProvider>
413 template <typename Clock>
414 template <typename T>
416  MPIBenchmarkReport<typename Clock::duration> &report, int rootThread)
417 {
418  Datatype dt = determineDatatype<T>();
419  auto dsf = std::dynamic_pointer_cast<DatasetFiller<T>>(
420  m_benchmark->m_dfp.template operator()<T>());
421  for (auto const &config : m_benchmark->m_configurations)
422  {
423  std::string compression;
424  uint8_t compressionLevel;
425  std::string backend;
426  int size;
427  Datatype dt2;
428  typename decltype(Series::iterations)::key_type iterations;
429  std::tie(
430  compression, compressionLevel, backend, size, dt2, iterations) =
431  config;
432 
433  if (dt != dt2)
434  {
435  continue;
436  }
437 
438  auto localCuboid = m_benchmark->slice(size);
439 
440  extentT blockSize = 1;
441  for (auto ext : localCuboid.second)
442  {
443  blockSize *= ext;
444  }
445  dsf->setNumberOfItems(blockSize);
446 
447  auto writeTime = writeBenchmark<T>(
448  compression,
449  compressionLevel,
450  localCuboid.first,
451  localCuboid.second,
452  backend,
453  dsf,
454  iterations);
455  auto readTime = readBenchmark<T>(
456  localCuboid.first, localCuboid.second, backend, iterations);
457  report.addReport(
458  rootThread,
459  compression,
460  compressionLevel,
461  backend,
462  size,
463  dt2,
464  iterations,
465  std::make_pair(writeTime, readTime));
466  }
467 }
468 
469 template <typename DatasetFillerProvider>
470 template <typename Clock>
471 template <int n>
474 {
475  throw std::runtime_error(
476  "Unknown/unsupported datatype requested to be benchmarked.");
477 }
478 
479 } // namespace openPMD
480 
481 #endif
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:191
Definition: Dataset.hpp:35
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
An abstract class to create one iteration of data per thread.
Definition: DatasetFiller.hpp:35
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:240
Datatype
Concrete datatype of an object available at runtime.
Definition: Datatype.hpp:45
Root level of the openPMD hierarchy.
Definition: Series.hpp:537
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:259
Public definitions of openPMD-api.
Definition: Date.cpp:28
Definition: MeshRecordComponent.hpp:29
void flush()
Execute all required remaining IO operations to write or read data.
Definition: Series.cpp:350
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:287