55 #include <Kokkos_ArithTraits.hpp>
57 #if defined(HAVE_MUELU_TPETRA) && defined(HAVE_MPI)
58 #include <Xpetra_TpetraImport.hpp>
59 #include <Tpetra_Import.hpp>
60 #include <Tpetra_Distributor.hpp>
70 namespace PerfDetails {
71 template <
class Scalar,
class Node>
74 using impl_scalar_type =
typename Kokkos::ArithTraits<Scalar>::val_type;
76 using exec_space =
typename Node::execution_space;
77 using memory_space =
typename Node::memory_space;
78 using range_policy = Kokkos::RangePolicy<exec_space>;
80 Kokkos::View<impl_scalar_type *, memory_space> a(
"a", VECTOR_SIZE);
81 Kokkos::View<impl_scalar_type *, memory_space> b(
"b", VECTOR_SIZE);
82 Kokkos::View<impl_scalar_type *, memory_space> c(
"c", VECTOR_SIZE);
83 double total_test_time = 0.0;
88 "stream/fill", range_policy(0, VECTOR_SIZE), KOKKOS_LAMBDA(
const size_t i) {
89 a(i) = ONE * (double)i;
94 using clock = std::chrono::high_resolution_clock;
98 for (
int i = 0; i < KERNEL_REPEATS; i++) {
100 Kokkos::parallel_for(
101 "stream/add", range_policy(0, VECTOR_SIZE), KOKKOS_LAMBDA(
const size_t j) {
105 exec_space().fence();
107 double my_test_time = std::chrono::duration<double>(stop - start).count();
108 total_test_time += my_test_time;
111 return total_test_time / KERNEL_REPEATS;
114 template <
class Scalar,
class Node>
117 using impl_scalar_type =
typename Kokkos::ArithTraits<Scalar>::val_type;
119 using exec_space =
typename Node::execution_space;
120 using memory_space =
typename Node::memory_space;
121 using range_policy = Kokkos::RangePolicy<exec_space>;
123 Kokkos::View<impl_scalar_type *, memory_space> a(
"a", VECTOR_SIZE);
124 Kokkos::View<impl_scalar_type *, memory_space> b(
"b", VECTOR_SIZE);
125 double total_test_time = 0.0;
129 Kokkos::parallel_for(
130 "stream/fill", range_policy(0, VECTOR_SIZE), KOKKOS_LAMBDA(
const size_t i) {
133 exec_space().fence();
135 using clock = std::chrono::high_resolution_clock;
138 for (
int i = 0; i < KERNEL_REPEATS; i++) {
139 start = clock::now();
140 Kokkos::parallel_for(
141 "stream/copy", range_policy(0, VECTOR_SIZE), KOKKOS_LAMBDA(
const size_t j) {
145 exec_space().fence();
147 double my_test_time = std::chrono::duration<double>(stop - start).count();
148 total_test_time += my_test_time;
151 return total_test_time / KERNEL_REPEATS;
154 double table_lookup(
const std::vector<int> &x,
const std::vector<double> &y,
int value) {
159 int N = (int)x.size();
161 for (; hi < N; hi++) {
170 }
else if (hi == N) {
174 int run = x[hi] - x[hi - 1];
175 double rise = y[hi] - y[hi - 1];
176 double slope = rise / run;
177 int diff = value - x[hi - 1];
179 return y[hi - 1] + slope * diff;
183 int run = x[hi] - x[hi - 1];
184 double rise = y[hi] - y[hi - 1];
185 double slope = rise / run;
186 int diff = value - x[hi - 1];
188 return y[hi - 1] + slope * diff;
193 const double GB = 1024.0 * 1024.0 * 1024.0;
195 double time_per_call = time / num_calls;
196 return memory_per_call_bytes /
GB / time_per_call;
199 template <
class exec_space,
class memory_space>
205 if (nproc < 2)
return;
207 const int buff_size = (int)pow(2, MAX_SIZE);
209 sizes.resize(MAX_SIZE + 1);
210 times.resize(MAX_SIZE + 1);
213 Kokkos::View<char *, memory_space> r_buf(
"recv", buff_size), s_buf(
"send", buff_size);
214 Kokkos::deep_copy(s_buf, 1);
219 int buddy = odd ? rank - 1 : rank + 1;
221 for (
int i = 0; i < MAX_SIZE + 1; i++) {
222 int msg_size = (int)pow(2, i);
225 double t0 = MPI_Wtime();
226 for (
int j = 0; j < KERNEL_REPEATS; j++) {
229 comm.
send(msg_size, (
char *)s_buf.data(), buddy);
230 comm.
receive(buddy, msg_size, (
char *)r_buf.data());
232 comm.
receive(buddy, msg_size, (
char *)r_buf.data());
233 comm.
send(msg_size, (
char *)s_buf.data(), buddy);
238 double time_per_call = (MPI_Wtime() - t0) / (2.0 * KERNEL_REPEATS);
240 times[i] = time_per_call;
247 template <
class exec_space,
class memory_space,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
249 int nproc =
import->getSourceMap()->getComm()->getSize();
250 if (nproc < 2)
return;
251 #if defined(HAVE_MUELU_TPETRA) && defined(HAVE_MPI)
257 MPI_Comm communicator = *mcomm->getRawMpiComm();
259 if (Ximport.
is_null() || mcomm.is_null())
return;
260 auto Timport = Ximport->getTpetra_Import();
261 auto distor = Timport->getDistributor();
266 int num_recvs = (int)distor.getNumReceives();
267 int num_sends = (int)distor.getNumSends();
269 const int buff_size_per_msg = (int)pow(2, MAX_SIZE);
270 sizes.resize(MAX_SIZE + 1);
271 times.resize(MAX_SIZE + 1);
274 Kokkos::View<char *, memory_space> f_recv_buf(
"forward_recv", buff_size_per_msg * num_recvs), f_send_buf(
"forward_send", buff_size_per_msg * num_sends);
275 Kokkos::View<char *, memory_space> r_recv_buf(
"reverse_recv", buff_size_per_msg * num_sends), r_send_buf(
"reverse_send", buff_size_per_msg * num_recvs);
276 Kokkos::deep_copy(f_send_buf, 1);
277 Kokkos::deep_copy(r_send_buf, 1);
279 std::vector<MPI_Request> requests(num_sends + num_recvs);
280 std::vector<MPI_Status> status(num_sends + num_recvs);
282 for (
int i = 0; i < MAX_SIZE + 1; i++) {
283 int msg_size = (int)pow(2, i);
285 MPI_Barrier(communicator);
287 double t0 = MPI_Wtime();
288 for (
int j = 0; j < KERNEL_REPEATS; j++) {
291 for (
int r = 0; r < num_recvs; r++) {
292 const int tag = 1000 + j;
293 MPI_Irecv(f_recv_buf.data() + msg_size * r, msg_size, MPI_CHAR, procsFrom[r], tag, communicator, &requests[ct]);
296 for (
int s = 0; s < num_sends; s++) {
297 const int tag = 1000 + j;
298 MPI_Isend(f_send_buf.data() + msg_size * s, msg_size, MPI_CHAR, procsTo[s], tag, communicator, &requests[ct]);
302 MPI_Waitall(ct, requests.data(), status.data());
306 for (
int r = 0; r < num_sends; r++) {
307 const int tag = 2000 + j;
308 MPI_Irecv(r_recv_buf.data() + msg_size * r, msg_size, MPI_CHAR, procsTo[r], tag, communicator, &requests[ct]);
311 for (
int s = 0; s < num_recvs; s++) {
312 const int tag = 2000 + j;
313 MPI_Isend(r_send_buf.data() + msg_size * s, msg_size, MPI_CHAR, procsFrom[s], tag, communicator, &requests[ct]);
317 MPI_Waitall(ct, requests.data(), status.data());
320 double time_per_call = (MPI_Wtime() - t0) / (2.0 * KERNEL_REPEATS);
322 times[i] = time_per_call;
330 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
332 : launch_and_wait_latency_(-1.0) {}
338 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
341 launch_latency_make_table(KERNEL_REPEATS);
342 double latency = launch_latency_lookup();
344 if (LOG_MAX_SIZE < 2)
347 stream_sizes_.resize(LOG_MAX_SIZE + 1);
348 stream_copy_times_.resize(LOG_MAX_SIZE + 1);
349 stream_add_times_.resize(LOG_MAX_SIZE + 1);
350 latency_corrected_stream_copy_times_.resize(LOG_MAX_SIZE + 1);
351 latency_corrected_stream_add_times_.resize(LOG_MAX_SIZE + 1);
353 for (
int i = 0; i < LOG_MAX_SIZE + 1; i++) {
354 int size = (int)pow(2, i);
355 double c_time = PerfDetails::stream_vector_copy<Scalar, Node>(KERNEL_REPEATS, size);
356 double a_time = PerfDetails::stream_vector_add<Scalar, Node>(KERNEL_REPEATS, size);
358 stream_sizes_[i] = size;
361 stream_copy_times_[i] = c_time / 2.0;
362 stream_add_times_[i] = a_time / 3.0;
366 latency_corrected_stream_copy_times_[i] = (c_time - latency <= 0.0) ? c_time / 2.0 : ((c_time - latency) / 2.0);
367 latency_corrected_stream_add_times_[i] = (a_time - latency <= 0.0) ? a_time / 3.0 : ((a_time - latency) / 3.0);
371 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
377 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
383 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
386 return std::min(stream_vector_copy_lookup(SIZE_IN_BYTES), stream_vector_add_lookup(SIZE_IN_BYTES));
389 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
395 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
401 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
404 return std::min(latency_corrected_stream_vector_copy_lookup(SIZE_IN_BYTES), latency_corrected_stream_vector_add_lookup(SIZE_IN_BYTES));
407 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
409 print_stream_vector_table_impl(out,
false, prefix);
412 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
414 print_stream_vector_table_impl(out,
true, prefix);
417 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
420 std::ios old_format(NULL);
421 old_format.copyfmt(out);
424 << setw(20) <<
"Length in Scalars" << setw(1) <<
" "
425 << setw(20) <<
"COPY (us)" << setw(1) <<
" "
426 << setw(20) <<
"ADD (us)" << setw(1) <<
" "
427 << setw(20) <<
"COPY (GB/s)" << setw(1) <<
" "
428 << setw(20) <<
"ADD (GB/s)" << std::endl;
431 << setw(20) <<
"-----------------" << setw(1) <<
" "
432 << setw(20) <<
"---------" << setw(1) <<
" "
433 << setw(20) <<
"--------" << setw(1) <<
" "
434 << setw(20) <<
"-----------" << setw(1) <<
" "
435 << setw(20) <<
"----------" << std::endl;
437 for (
int i = 0; i < (int)stream_sizes_.size(); i++) {
438 int size = stream_sizes_[i];
439 double c_time = use_latency_correction ? latency_corrected_stream_copy_times_[i] : stream_copy_times_[i];
440 double a_time = use_latency_correction ? latency_corrected_stream_add_times_[i] : stream_add_times_[i];
446 << setw(20) << size << setw(1) <<
" "
447 << setw(20) << fixed << setprecision(4) << (c_time * 1e6) << setw(1) <<
" "
448 << setw(20) << fixed << setprecision(4) << (a_time * 1e6) << setw(1) <<
" "
449 << setw(20) << fixed << setprecision(4) << c_bw << setw(1) <<
" "
450 << setw(20) << fixed << setprecision(4) << a_bw << std::endl;
453 out.copyfmt(old_format);
460 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
462 PerfDetails::pingpong_basic<Kokkos::HostSpace::execution_space, Kokkos::HostSpace::memory_space>(KERNEL_REPEATS, LOG_MAX_SIZE, *comm, pingpong_sizes_, pingpong_host_times_);
464 PerfDetails::pingpong_basic<typename Node::execution_space, typename Node::memory_space>(KERNEL_REPEATS, LOG_MAX_SIZE, *comm, pingpong_sizes_, pingpong_device_times_);
467 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
473 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
479 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
481 if (pingpong_sizes_.size() == 0)
return;
484 std::ios old_format(NULL);
485 old_format.copyfmt(out);
488 << setw(20) <<
"Message Size" << setw(1) <<
" "
489 << setw(20) <<
"Host (us)" << setw(1) <<
" "
490 << setw(20) <<
"Device (us)" << std::endl;
493 << setw(20) <<
"------------" << setw(1) <<
" "
494 << setw(20) <<
"---------" << setw(1) <<
" "
495 << setw(20) <<
"-----------" << std::endl;
497 for (
int i = 0; i < (int)pingpong_sizes_.size(); i++) {
498 int size = pingpong_sizes_[i];
499 double h_time = pingpong_host_times_[i];
500 double d_time = pingpong_device_times_[i];
503 << setw(20) << size << setw(1) <<
" "
504 << setw(20) << fixed << setprecision(4) << (h_time * 1e6) << setw(1) <<
" "
505 << setw(20) << fixed << setprecision(4) << (d_time * 1e6) << setw(1) << std::endl;
508 out.copyfmt(old_format);
514 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
516 PerfDetails::halopong_basic<Kokkos::HostSpace::execution_space, Kokkos::HostSpace::memory_space>(KERNEL_REPEATS, LOG_MAX_SIZE,
import, halopong_sizes_, halopong_host_times_);
518 PerfDetails::halopong_basic<typename Node::execution_space, typename Node::memory_space>(KERNEL_REPEATS, LOG_MAX_SIZE,
import, halopong_sizes_, halopong_device_times_);
521 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
527 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
533 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
535 if (halopong_sizes_.size() == 0)
return;
538 std::ios old_format(NULL);
539 old_format.copyfmt(out);
542 << setw(20) <<
"Message Size" << setw(1) <<
" "
543 << setw(20) <<
"Host (us)" << setw(1) <<
" "
544 << setw(20) <<
"Device (us)" << std::endl;
547 << setw(20) <<
"------------" << setw(1) <<
" "
548 << setw(20) <<
"---------" << setw(1) <<
" "
549 << setw(20) <<
"-----------" << std::endl;
551 for (
int i = 0; i < (int)halopong_sizes_.size(); i++) {
552 int size = halopong_sizes_[i];
553 double h_time = halopong_host_times_[i];
554 double d_time = halopong_device_times_[i];
557 << setw(20) << size << setw(1) <<
" "
558 << setw(20) << fixed << setprecision(4) << (h_time * 1e6) << setw(1) <<
" "
559 << setw(20) << fixed << setprecision(4) << (d_time * 1e6) << setw(1) << std::endl;
562 out.copyfmt(old_format);
569 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
571 using exec_space =
typename Node::execution_space;
572 using range_policy = Kokkos::RangePolicy<exec_space>;
573 using clock = std::chrono::high_resolution_clock;
575 double total_test_time = 0;
577 for (
int i = 0; i < KERNEL_REPEATS; i++) {
578 start = clock::now();
579 Kokkos::parallel_for(
580 "empty kernel", range_policy(0, 1), KOKKOS_LAMBDA(
const size_t j) {
583 exec_space().fence();
585 double my_test_time = std::chrono::duration<double>(stop - start).count();
586 total_test_time += my_test_time;
589 launch_and_wait_latency_ = total_test_time / KERNEL_REPEATS;
592 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
595 return launch_and_wait_latency_;
598 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
601 std::ios old_format(NULL);
602 old_format.copyfmt(out);
605 << setw(20) <<
"Launch+Wait Latency (us)" << setw(1) <<
" "
606 << setw(20) << fixed << setprecision(4) << (launch_and_wait_latency_ * 1e6) << std::endl;
608 out.copyfmt(old_format);
void halopong_basic(int KERNEL_REPEATS, int MAX_SIZE, const RCP< const Xpetra::Import< LocalOrdinal, GlobalOrdinal, Node > > &import, std::vector< int > &sizes, std::vector< double > ×)
virtual int getSize() const =0
virtual int getRank() const =0
double convert_time_to_bandwidth_gbs(double time, int num_calls, double memory_per_call_bytes)
double stream_vector_copy(int KERNEL_REPEATS, int VECTOR_SIZE)
double pingpong_device_lookup(int SIZE_IN_BYTES)
virtual int receive(const int sourceRank, const Ordinal bytes, char recvBuffer[]) const =0
void stream_vector_make_table(int KERNEL_REPEATS, int LOG_MAX_SIZE=20)
void print_latency_corrected_stream_vector_table(std::ostream &out, const std::string &prefix="")
double latency_corrected_stream_vector_add_lookup(int SIZE_IN_BYTES)
double latency_corrected_stream_vector_lookup(int SIZE_IN_BYTES)
void pingpong_basic(int KERNEL_REPEATS, int MAX_SIZE, const Teuchos::Comm< int > &comm, std::vector< int > &sizes, std::vector< double > ×)
double table_lookup(const std::vector< int > &x, const std::vector< double > &y, int value)
virtual void barrier() const =0
double stream_vector_lookup(int SIZE_IN_BYTES)
MueLu::DefaultScalar Scalar
double stream_vector_copy_lookup(int SIZE_IN_BYTES)
double halopong_host_lookup(int SIZE_IN_BYTES_PER_MESSAGE)
double halopong_device_lookup(int SIZE_IN_BYTES_PER_MESSAGE)
void print_stream_vector_table_impl(std::ostream &out, bool use_latency_correction, const std::string &prefix)
double latency_corrected_stream_vector_copy_lookup(int SIZE_IN_BYTES)
void print_launch_latency_table(std::ostream &out, const std::string &prefix="")
void print_halopong_table(std::ostream &out, const std::string &prefix="")
void print_pingpong_table(std::ostream &out, const std::string &prefix="")
double pingpong_host_lookup(int SIZE_IN_BYTES)
void halopong_make_table(int KERNEL_REPEATS, int LOG_MAX_SIZE, const RCP< const Xpetra::Import< LocalOrdinal, GlobalOrdinal, Node > > &import)
virtual void send(const Ordinal bytes, const char sendBuffer[], const int destRank) const =0
void launch_latency_make_table(int KERNEL_REPEATS)
double stream_vector_add(int KERNEL_REPEATS, int VECTOR_SIZE)
double launch_latency_lookup()
double stream_vector_add_lookup(int SIZE_IN_BYTES)
void print_stream_vector_table(std::ostream &out, const std::string &prefix="")
void pingpong_make_table(int KERNEL_REPEATS, int LOG_MAX_SIZE, const RCP< const Teuchos::Comm< int > > &comm)