openPMD-api
MPIBenchmarkReport.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 "openPMD/Datatype.hpp"
28 #include "openPMD/Series.hpp"
29 
30 #include "string.h"
31 #include <map>
32 #include <mpi.h>
33 #include <tuple>
34 #include <vector>
35 
36 namespace openPMD
37 {
43 template <typename Duration>
45 {
46  MPI_Comm communicator;
47 
48  MPIBenchmarkReport(MPI_Comm);
49 
53  std::map<
54  std::tuple<
55  int, // rank
56  std::string, // compression
57  uint8_t, // compression level
58  std::string, // extension
59  int, // thread size
60  Datatype,
61  typename decltype(Series::iterations)::key_type>,
62  std::pair<Duration, Duration> >
64 
65  enum Selector
66  {
67  RANK = 0,
68  COMPRESSION,
69  COMPRESSION_LEVEL,
70  BACKEND,
71  NRANKS,
72  DTYPE,
73  ITERATIONS
74  };
75 
88  void addReport(
89  int rootThread,
90  std::string compression,
91  uint8_t level,
92  std::string extension,
93  int threadSize,
94  Datatype dt,
95  typename decltype(Series::iterations)::key_type iterations,
96  std::pair<Duration, Duration> const &report);
97 
109  std::pair<Duration, Duration> getReport(
110  int rank,
111  std::string compression,
112  uint8_t level,
113  std::string extension,
114  int threadSize,
115  Datatype dt,
116  typename decltype(Series::iterations)::key_type iterations);
117 
118 private:
119  template <typename D, typename Dummy = D>
120  struct MPIDatatype
121  {};
122 
123  template <typename Dummy>
124  struct MPIDatatype<char, Dummy>
125  {
126  MPI_Datatype dt = MPI_CHAR;
127  };
128  template <typename Dummy>
129  struct MPIDatatype<unsigned char, Dummy>
130  {
131  MPI_Datatype dt = MPI_UNSIGNED_CHAR;
132  };
133  template <typename Dummy>
134  struct MPIDatatype<short, Dummy>
135  {
136  MPI_Datatype dt = MPI_SHORT;
137  };
138  template <typename Dummy>
139  struct MPIDatatype<int, Dummy>
140  {
141  MPI_Datatype dt = MPI_INT;
142  };
143  template <typename Dummy>
144  struct MPIDatatype<long, Dummy>
145  {
146  MPI_Datatype dt = MPI_LONG;
147  };
148  template <typename Dummy>
149  struct MPIDatatype<float, Dummy>
150  {
151  MPI_Datatype dt = MPI_FLOAT;
152  };
153  template <typename Dummy>
154  struct MPIDatatype<double, Dummy>
155  {
156  MPI_Datatype dt = MPI_DOUBLE;
157  };
158  template <typename Dummy>
159  struct MPIDatatype<unsigned short, Dummy>
160  {
161  MPI_Datatype dt = MPI_UNSIGNED_SHORT;
162  };
163  template <typename Dummy>
164  struct MPIDatatype<unsigned int, Dummy>
165  {
166  MPI_Datatype dt = MPI_UNSIGNED;
167  };
168  template <typename Dummy>
169  struct MPIDatatype<unsigned long, Dummy>
170  {
171  MPI_Datatype dt = MPI_UNSIGNED_LONG;
172  };
173  template <typename Dummy>
174  struct MPIDatatype<long double, Dummy>
175  {
176  MPI_Datatype dt = MPI_LONG_DOUBLE;
177  };
178  template <typename Dummy>
179  struct MPIDatatype<long long, Dummy>
180  {
181  MPI_Datatype dt = MPI_LONG_LONG_INT;
182  };
183 
184  MPIDatatype<typename Duration::rep> m_mpiDatatype;
185  MPI_Datatype mpiType = m_mpiDatatype.dt;
186 };
187 
188 // implementation
189 
190 template <typename Duration>
192  int rootThread,
193  std::string compression,
194  uint8_t level,
195  std::string extension,
196  int threadSize,
197  Datatype dt,
198  typename decltype(Series::iterations)::key_type iterations,
199  std::pair<Duration, Duration> const &report)
200 {
201  using rep = typename Duration::rep;
202  // auto mpi_dt = MPIDatatype<rep>::dt;
203  int rank;
204  MPI_Comm_rank(communicator, &rank);
205  int size;
206  MPI_Comm_size(communicator, &size);
207  MPI_Comm restricted;
208  MPI_Comm_split(
209  communicator, rank < threadSize ? 0 : MPI_UNDEFINED, rank, &restricted);
210  rep readWrite[2];
211  if (rank < threadSize)
212  {
213  readWrite[0] = report.first.count();
214  readWrite[1] = report.second.count();
215  }
216  rep *recv = nullptr;
217  if (rank == rootThread)
218  {
219  recv = new rep[2 * threadSize];
220  }
221 
222  if (restricted != MPI_COMM_NULL)
223  {
224  MPI_Gather(
225  readWrite,
226  2, // should be 2 but doesnt work then..
227  this->mpiType,
228  recv,
229  2,
230  this->mpiType,
231  rootThread,
232  restricted);
233  }
234 
235  if (rank == rootThread)
236  {
237  for (int i = 0; i < threadSize; i++)
238  {
239  Duration dWrite{recv[2 * i]};
240  Duration dRead{recv[2 * i + 1]};
241  this->durations.emplace(
242  std::make_tuple(
243  i,
244  compression,
245  level,
246  extension,
247  threadSize,
248  dt,
249  iterations),
250  std::make_pair(dWrite, dRead));
251  }
252  delete[] recv;
253  }
254  if (restricted != MPI_COMM_NULL)
255  {
256  MPI_Comm_free(&restricted);
257  }
258 }
259 
260 template <typename Duration>
262  : communicator{comm}
263 {}
264 
265 template <typename Duration>
266 std::pair<Duration, Duration> MPIBenchmarkReport<Duration>::getReport(
267  int rank,
268  std::string compression,
269  uint8_t level,
270  std::string extension,
271  int threadSize,
272  Datatype dt,
273  typename decltype(Series::iterations)::key_type iterations)
274 {
275  auto it = this->durations.find(std::make_tuple(
276  rank, compression, level, extension, threadSize, dt, iterations));
277  if (it == this->durations.end())
278  {
279  throw std::runtime_error(
280  "Requested report not found. (Reports are available on the root "
281  "thread only)");
282  }
283  else
284  {
285  return it->second;
286  }
287 }
288 
289 } // namespace openPMD
290 
291 #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
std::map< std::tuple< int, std::string, uint8_t, std::string, int, Datatype, typename decltype(Series::iterations)::key_type >, std::pair< Duration, Duration > > durations
Time needed for writing and reading per compression strategy and level.
Definition: MPIBenchmarkReport.hpp:63
Datatype
Concrete datatype of an object available at runtime.
Definition: Datatype.hpp:45
The report for a single benchmark produced by <openPMD/benchmark/mpi/MPIBenchmark>.
Definition: MPIBenchmarkReport.hpp:44
Public definitions of openPMD-api.
Definition: Date.cpp:28
std::pair< Duration, Duration > getReport(int rank, std::string compression, uint8_t level, std::string extension, int threadSize, Datatype dt, typename decltype(Series::iterations)::key_type iterations)
Retrieve the time measured for a certain compression strategy.
Definition: MPIBenchmarkReport.hpp:266