DOLFINx 0.10.0.0
DOLFINx C++ interface
Loading...
Searching...
No Matches
MPI.h
1// Copyright (C) 2007-2023 Magnus Vikstrøm, Garth N. Wells and Paul T. Kühner
2//
3// This file is part of DOLFINx (https://www.fenicsproject.org)
4//
5// SPDX-License-Identifier: LGPL-3.0-or-later
6
7#pragma once
8
9#include "Timer.h"
10#include "log.h"
11#include "types.h"
12#include <algorithm>
13#include <array>
14#include <cassert>
15#include <complex>
16#include <cstdint>
17#include <dolfinx/graph/AdjacencyList.h>
18#include <numeric>
19#include <set>
20#include <span>
21#include <tuple>
22#include <type_traits>
23#include <utility>
24#include <vector>
25
26#define MPICH_IGNORE_CXX_SEEK 1
27#include <mpi.h>
28
30namespace dolfinx::MPI
31{
33enum class tag : int
34{
35 consensus_pcx = 1200,
36 consensus_pex = 1201,
37 consensus_nbx = 1202,
38};
39
42class Comm
43{
44public:
46 explicit Comm(MPI_Comm comm, bool duplicate = true);
47
49 Comm(const Comm& comm) noexcept;
50
52 Comm(Comm&& comm) noexcept;
53
54 // Disable copy assignment operator
55 Comm& operator=(const Comm& comm) = delete;
56
58 Comm& operator=(Comm&& comm) noexcept;
59
61 ~Comm();
62
64 MPI_Comm comm() const noexcept;
65
66private:
67 // MPI communicator
68 MPI_Comm _comm;
69};
70
72int rank(MPI_Comm comm);
73
76int size(MPI_Comm comm);
77
82void check_error(MPI_Comm comm, int code);
83
90constexpr std::array<std::int64_t, 2> local_range(int rank, std::int64_t N,
91 int size)
92{
93 assert(rank >= 0);
94 assert(N >= 0);
95 assert(size > 0);
96
97 // Compute number of items per rank and remainder
98 const std::int64_t n = N / size;
99 const std::int64_t r = N % size;
100
101 // Compute local range
102 if (rank < r)
103 return {rank * (n + 1), rank * (n + 1) + n + 1};
104 else
105 return {rank * n + r, rank * n + r + n};
106}
107
114constexpr int index_owner(int size, std::size_t index, std::size_t N)
115{
116 assert(index < N);
117
118 // Compute number of items per rank and remainder
119 const std::size_t n = N / size;
120 const std::size_t r = N % size;
121
122 if (index < r * (n + 1))
123 {
124 // First r ranks own n + 1 indices
125 return index / (n + 1);
126 }
127 else
128 {
129 // Remaining ranks own n indices
130 return r + (index - r * (n + 1)) / n;
131 }
132}
133
157std::vector<int> compute_graph_edges_pcx(MPI_Comm comm,
158 std::span<const int> edges);
159
185std::vector<int>
186compute_graph_edges_nbx(MPI_Comm comm, std::span<const int> edges,
187 int tag = static_cast<int>(tag::consensus_nbx));
188
208template <typename U>
209std::pair<std::vector<std::int32_t>,
210 std::vector<typename std::remove_reference_t<typename U::value_type>>>
211distribute_to_postoffice(MPI_Comm comm, const U& x,
212 std::array<std::int64_t, 2> shape,
213 std::int64_t rank_offset);
214
236template <typename U>
237std::vector<typename std::remove_reference_t<typename U::value_type>>
238distribute_from_postoffice(MPI_Comm comm, std::span<const std::int64_t> indices,
239 const U& x, std::array<std::int64_t, 2> shape,
240 std::int64_t rank_offset);
241
262template <typename U>
263std::vector<typename std::remove_reference_t<typename U::value_type>>
264distribute_data(MPI_Comm comm0, std::span<const std::int64_t> indices,
265 MPI_Comm comm1, const U& x, int shape1);
266
267template <typename T>
268struct dependent_false : std::false_type
269{
270};
271
273
275template <typename T>
277
280template <typename T>
282
285#define MAP_TO_MPI_TYPE(cpp_t, mpi_t) \
286 template <> \
287 struct mpi_type_mapping<cpp_t> \
288 { \
289 static inline MPI_Datatype type = mpi_t; \
290 };
291
295MAP_TO_MPI_TYPE(float, MPI_FLOAT)
296MAP_TO_MPI_TYPE(double, MPI_DOUBLE)
297MAP_TO_MPI_TYPE(std::complex<float>, MPI_C_FLOAT_COMPLEX)
298MAP_TO_MPI_TYPE(std::complex<double>, MPI_C_DOUBLE_COMPLEX)
299MAP_TO_MPI_TYPE(std::int8_t, MPI_INT8_T)
300MAP_TO_MPI_TYPE(std::int16_t, MPI_INT16_T)
301MAP_TO_MPI_TYPE(std::int32_t, MPI_INT32_T)
302MAP_TO_MPI_TYPE(std::int64_t, MPI_INT64_T)
303MAP_TO_MPI_TYPE(std::uint8_t, MPI_UINT8_T)
304MAP_TO_MPI_TYPE(std::uint16_t, MPI_UINT16_T)
305MAP_TO_MPI_TYPE(std::uint32_t, MPI_UINT32_T)
306MAP_TO_MPI_TYPE(std::uint64_t, MPI_UINT64_T)
309
310//---------------------------------------------------------------------------
311template <typename U>
312std::pair<std::vector<std::int32_t>,
313 std::vector<typename std::remove_reference_t<typename U::value_type>>>
314distribute_to_postoffice(MPI_Comm comm, const U& x,
315 std::array<std::int64_t, 2> shape,
316 std::int64_t rank_offset)
317{
318 assert(rank_offset >= 0 or x.empty());
319 using T = typename std::remove_reference_t<typename U::value_type>;
320
321 const int size = dolfinx::MPI::size(comm);
322 const int rank = dolfinx::MPI::rank(comm);
323 assert(x.size() % shape[1] == 0);
324 const std::int32_t shape0_local = x.size() / shape[1];
325
326 spdlog::debug("Sending data to post offices (distribute_to_postoffice)");
327
328 // Post office ranks will receive data from this rank
329 std::vector<int> row_to_dest(shape0_local);
330 for (std::int32_t i = 0; i < shape0_local; ++i)
331 {
332 int dest = MPI::index_owner(size, i + rank_offset, shape[0]);
333 row_to_dest[i] = dest;
334 }
335
336 // Build list of (dest, positions) for each row that doesn't belong to
337 // this rank, then sort
338 std::vector<std::array<std::int32_t, 2>> dest_to_index;
339 dest_to_index.reserve(shape0_local);
340 for (std::int32_t i = 0; i < shape0_local; ++i)
341 {
342 std::size_t idx = i + rank_offset;
343 if (int dest = MPI::index_owner(size, idx, shape[0]); dest != rank)
344 dest_to_index.push_back({dest, i});
345 }
346 std::ranges::sort(dest_to_index);
347
348 // Build list of neighbour src ranks and count number of items (rows
349 // of x) to receive from each src post office (by neighbourhood rank)
350 std::vector<int> dest;
351 std::vector<std::int32_t> num_items_per_dest,
352 pos_to_neigh_rank(shape0_local, -1);
353 {
354 auto it = dest_to_index.begin();
355 while (it != dest_to_index.end())
356 {
357 const int neigh_rank = dest.size();
358
359 // Store global rank
360 dest.push_back((*it)[0]);
361
362 // Find iterator to next global rank
363 auto it1
364 = std::find_if(it, dest_to_index.end(),
365 [r = dest.back()](auto& idx) { return idx[0] != r; });
366
367 // Store number of items for current rank
368 num_items_per_dest.push_back(std::distance(it, it1));
369
370 // Map from local x index to local destination rank
371 for (auto e = it; e != it1; ++e)
372 pos_to_neigh_rank[(*e)[1]] = neigh_rank;
373
374 // Advance iterator
375 it = it1;
376 }
377 }
378
379 // Determine source ranks
380 const std::vector<int> src = MPI::compute_graph_edges_nbx(comm, dest);
381 spdlog::info(
382 "Number of neighbourhood source ranks in distribute_to_postoffice: {}",
383 static_cast<int>(src.size()));
384
385 // Create neighbourhood communicator for sending data to post offices
386 MPI_Comm neigh_comm;
387 int err = MPI_Dist_graph_create_adjacent(
388 comm, src.size(), src.data(), MPI_UNWEIGHTED, dest.size(), dest.data(),
389 MPI_UNWEIGHTED, MPI_INFO_NULL, false, &neigh_comm);
390 dolfinx::MPI::check_error(comm, err);
391
392 // Compute send displacements
393 std::vector<std::int32_t> send_disp = {0};
394 std::partial_sum(num_items_per_dest.begin(), num_items_per_dest.end(),
395 std::back_inserter(send_disp));
396
397 // Pack send buffers
398 std::vector<T> send_buffer_data(shape[1] * send_disp.back());
399 std::vector<std::int64_t> send_buffer_index(send_disp.back());
400 {
401 std::vector<std::int32_t> send_offsets = send_disp;
402 for (std::int32_t i = 0; i < shape0_local; ++i)
403 {
404 if (int neigh_dest = pos_to_neigh_rank[i]; neigh_dest != -1)
405 {
406 std::size_t pos = send_offsets[neigh_dest];
407 send_buffer_index[pos] = i + rank_offset;
408 std::copy_n(std::next(x.begin(), i * shape[1]), shape[1],
409 std::next(send_buffer_data.begin(), shape[1] * pos));
410 ++send_offsets[neigh_dest];
411 }
412 }
413 }
414
415 // Send number of items to post offices (destination) that I will be
416 // sending
417 std::vector<int> num_items_recv(src.size());
418 num_items_per_dest.reserve(1);
419 num_items_recv.reserve(1);
420 err = MPI_Neighbor_alltoall(num_items_per_dest.data(), 1, MPI_INT,
421 num_items_recv.data(), 1, MPI_INT, neigh_comm);
422 dolfinx::MPI::check_error(comm, err);
423
424 // Prepare receive displacement and buffers
425 std::vector<std::int32_t> recv_disp(num_items_recv.size() + 1, 0);
426 std::partial_sum(num_items_recv.begin(), num_items_recv.end(),
427 std::next(recv_disp.begin()));
428
429 // Send/receive global indices
430 std::vector<std::int64_t> recv_buffer_index(recv_disp.back());
431 err = MPI_Neighbor_alltoallv(
432 send_buffer_index.data(), num_items_per_dest.data(), send_disp.data(),
433 MPI_INT64_T, recv_buffer_index.data(), num_items_recv.data(),
434 recv_disp.data(), MPI_INT64_T, neigh_comm);
435 dolfinx::MPI::check_error(comm, err);
436
437 // Send/receive data (x)
438 MPI_Datatype compound_type;
439 MPI_Type_contiguous(shape[1], dolfinx::MPI::mpi_t<T>, &compound_type);
440 MPI_Type_commit(&compound_type);
441 std::vector<T> recv_buffer_data(shape[1] * recv_disp.back());
442 err = MPI_Neighbor_alltoallv(
443 send_buffer_data.data(), num_items_per_dest.data(), send_disp.data(),
444 compound_type, recv_buffer_data.data(), num_items_recv.data(),
445 recv_disp.data(), compound_type, neigh_comm);
446 dolfinx::MPI::check_error(comm, err);
447 err = MPI_Type_free(&compound_type);
448 dolfinx::MPI::check_error(comm, err);
449 err = MPI_Comm_free(&neigh_comm);
450 dolfinx::MPI::check_error(comm, err);
451
452 spdlog::debug("Completed send data to post offices.");
453
454 // Convert to local indices
455 const std::int64_t r0 = MPI::local_range(rank, shape[0], size)[0];
456 std::vector<std::int32_t> index_local(recv_buffer_index.size());
457 std::ranges::transform(recv_buffer_index, index_local.begin(),
458 [r0](auto idx) { return idx - r0; });
459
460 return {index_local, recv_buffer_data};
461}
462//---------------------------------------------------------------------------
463template <typename U>
464std::vector<typename std::remove_reference_t<typename U::value_type>>
465distribute_from_postoffice(MPI_Comm comm, std::span<const std::int64_t> indices,
466 const U& x, std::array<std::int64_t, 2> shape,
467 std::int64_t rank_offset)
468{
469 assert(rank_offset >= 0 or x.empty());
470 using T = typename std::remove_reference_t<typename U::value_type>;
471
472 common::Timer timer("Distribute row-wise data (scalable)");
473 assert(shape[1] > 0);
474
475 const int size = dolfinx::MPI::size(comm);
476 const int rank = dolfinx::MPI::rank(comm);
477 assert(x.size() % shape[1] == 0);
478 const std::int64_t shape0_local = x.size() / shape[1];
479
480 // 0. Send x data to/from post offices
481
482 // Send receive x data to post office (only for rows that need to be
483 // communicated)
484 auto [post_indices, post_x] = dolfinx::MPI::distribute_to_postoffice(
485 comm, x, {shape[0], shape[1]}, rank_offset);
486 assert(post_indices.size() == post_x.size() / shape[1]);
487
488 // 1. Send request to post office ranks for data
489
490 // Build list of (src, global index, global, index position) for each
491 // entry in 'indices' that doesn't belong to this rank, then sort
492 std::vector<std::tuple<int, std::int64_t, std::int32_t>> src_to_index;
493 for (std::size_t i = 0; i < indices.size(); ++i)
494 {
495 std::size_t idx = indices[i];
496 if (int src = dolfinx::MPI::index_owner(size, idx, shape[0]); src != rank)
497 src_to_index.push_back({src, idx, i});
498 }
499 std::ranges::sort(src_to_index);
500
501 // Build list is neighbour src ranks and count number of items (rows
502 // of x) to receive from each src post office (by neighbourhood rank)
503 std::vector<std::int32_t> num_items_per_src;
504 std::vector<int> src;
505 {
506 auto it = src_to_index.begin();
507 while (it != src_to_index.end())
508 {
509 src.push_back(std::get<0>(*it));
510 auto it1
511 = std::find_if(it, src_to_index.end(), [r = src.back()](auto& idx)
512 { return std::get<0>(idx) != r; });
513 num_items_per_src.push_back(std::distance(it, it1));
514 it = it1;
515 }
516 }
517
518 // Determine 'delivery' destination ranks (ranks that want data from
519 // me)
520 const std::vector<int> dest
522 spdlog::info(
523 "Neighbourhood destination ranks from post office in "
524 "distribute_data (rank, num dests, num dests/mpi_size): {}, {}, {}",
525 rank, static_cast<int>(dest.size()),
526 static_cast<double>(dest.size()) / size);
527
528 // Create neighbourhood communicator for sending data to post offices
529 // (src), and receiving data form my send my post office
530 MPI_Comm neigh_comm0;
531 int err = MPI_Dist_graph_create_adjacent(
532 comm, dest.size(), dest.data(), MPI_UNWEIGHTED, src.size(), src.data(),
533 MPI_UNWEIGHTED, MPI_INFO_NULL, false, &neigh_comm0);
534 dolfinx::MPI::check_error(comm, err);
535
536 // Communicate number of requests to each source
537 std::vector<int> num_items_recv(dest.size());
538 num_items_per_src.reserve(1);
539 num_items_recv.reserve(1);
540 err = MPI_Neighbor_alltoall(num_items_per_src.data(), 1, MPI_INT,
541 num_items_recv.data(), 1, MPI_INT, neigh_comm0);
542 dolfinx::MPI::check_error(comm, err);
543
544 // Prepare send/receive displacements
545 std::vector<std::int32_t> send_disp = {0};
546 std::partial_sum(num_items_per_src.begin(), num_items_per_src.end(),
547 std::back_inserter(send_disp));
548 std::vector<std::int32_t> recv_disp = {0};
549 std::partial_sum(num_items_recv.begin(), num_items_recv.end(),
550 std::back_inserter(recv_disp));
551
552 // Pack my requested indices (global) in send buffer ready to send to
553 // post offices
554 assert(send_disp.back() == (int)src_to_index.size());
555 std::vector<std::int64_t> send_buffer_index(src_to_index.size());
556 std::ranges::transform(src_to_index, send_buffer_index.begin(),
557 [](auto x) { return std::get<1>(x); });
558
559 // Prepare the receive buffer
560 std::vector<std::int64_t> recv_buffer_index(recv_disp.back());
561 err = MPI_Neighbor_alltoallv(
562 send_buffer_index.data(), num_items_per_src.data(), send_disp.data(),
563 MPI_INT64_T, recv_buffer_index.data(), num_items_recv.data(),
564 recv_disp.data(), MPI_INT64_T, neigh_comm0);
565 dolfinx::MPI::check_error(comm, err);
566
567 err = MPI_Comm_free(&neigh_comm0);
568 dolfinx::MPI::check_error(comm, err);
569
570 // 2. Send data (rows of x) from post office back to requesting ranks
571 // (transpose of the preceding communication pattern operation)
572
573 // Build map from local index to post_indices position. Set to -1 for
574 // data that was already on this rank and was therefore was not
575 // sent/received via a postoffice.
576 const std::array<std::int64_t, 2> postoffice_range
577 = dolfinx::MPI::local_range(rank, shape[0], size);
578 std::vector<std::int32_t> post_indices_map(
579 postoffice_range[1] - postoffice_range[0], -1);
580 for (std::size_t i = 0; i < post_indices.size(); ++i)
581 {
582 assert(post_indices[i] < (int)post_indices_map.size());
583 post_indices_map[post_indices[i]] = i;
584 }
585
586 // Build send buffer
587 std::vector<T> send_buffer_data(shape[1] * recv_disp.back());
588 for (std::size_t p = 0; p < recv_disp.size() - 1; ++p)
589 {
590 int offset = recv_disp[p];
591 for (std::int32_t i = recv_disp[p]; i < recv_disp[p + 1]; ++i)
592 {
593 std::int64_t index = recv_buffer_index[i];
594 if (index >= rank_offset and index < (rank_offset + shape0_local))
595 {
596 // I already had this index before any communication
597 std::int32_t local_index = index - rank_offset;
598 std::copy_n(std::next(x.begin(), shape[1] * local_index), shape[1],
599 std::next(send_buffer_data.begin(), shape[1] * offset));
600 }
601 else
602 {
603 // Take from my 'post bag'
604 auto local_index = index - postoffice_range[0];
605 std::int32_t pos = post_indices_map[local_index];
606 assert(pos != -1);
607 std::copy_n(std::next(post_x.begin(), shape[1] * pos), shape[1],
608 std::next(send_buffer_data.begin(), shape[1] * offset));
609 }
610
611 ++offset;
612 }
613 }
614
615 err = MPI_Dist_graph_create_adjacent(
616 comm, src.size(), src.data(), MPI_UNWEIGHTED, dest.size(), dest.data(),
617 MPI_UNWEIGHTED, MPI_INFO_NULL, false, &neigh_comm0);
618 dolfinx::MPI::check_error(comm, err);
619
620 MPI_Datatype compound_type0;
621 MPI_Type_contiguous(shape[1], dolfinx::MPI::mpi_t<T>, &compound_type0);
622 MPI_Type_commit(&compound_type0);
623
624 std::vector<T> recv_buffer_data(shape[1] * send_disp.back());
625 err = MPI_Neighbor_alltoallv(
626 send_buffer_data.data(), num_items_recv.data(), recv_disp.data(),
627 compound_type0, recv_buffer_data.data(), num_items_per_src.data(),
628 send_disp.data(), compound_type0, neigh_comm0);
629 dolfinx::MPI::check_error(comm, err);
630
631 err = MPI_Type_free(&compound_type0);
632 dolfinx::MPI::check_error(comm, err);
633 err = MPI_Comm_free(&neigh_comm0);
634 dolfinx::MPI::check_error(comm, err);
635
636 std::vector<std::int32_t> index_pos_to_buffer(indices.size(), -1);
637 for (std::size_t i = 0; i < src_to_index.size(); ++i)
638 index_pos_to_buffer[std::get<2>(src_to_index[i])] = i;
639
640 // Extra data to return
641 std::vector<T> x_new(shape[1] * indices.size());
642 for (std::size_t i = 0; i < indices.size(); ++i)
643 {
644 const std::int64_t index = indices[i];
645 if (index >= rank_offset and index < (rank_offset + shape0_local))
646 {
647 // Had data from the start in x
648 auto local_index = index - rank_offset;
649 std::copy_n(std::next(x.begin(), shape[1] * local_index), shape[1],
650 std::next(x_new.begin(), shape[1] * i));
651 }
652 else
653 {
654 if (int src = dolfinx::MPI::index_owner(size, index, shape[0]);
655 src == rank)
656 {
657 // In my post office bag
658 auto local_index = index - postoffice_range[0];
659 std::int32_t pos = post_indices_map[local_index];
660 assert(pos != -1);
661 std::copy_n(std::next(post_x.begin(), shape[1] * pos), shape[1],
662 std::next(x_new.begin(), shape[1] * i));
663 }
664 else
665 {
666 // In my received post
667 std::int32_t pos = index_pos_to_buffer[i];
668 assert(pos != -1);
669 std::copy_n(std::next(recv_buffer_data.begin(), shape[1] * pos),
670 shape[1], std::next(x_new.begin(), shape[1] * i));
671 }
672 }
673 }
674
675 return x_new;
676}
677//---------------------------------------------------------------------------
678template <typename U>
679std::vector<typename std::remove_reference_t<typename U::value_type>>
680distribute_data(MPI_Comm comm0, std::span<const std::int64_t> indices,
681 MPI_Comm comm1, const U& x, int shape1)
682{
683 assert(shape1 > 0);
684 assert(x.size() % shape1 == 0);
685 const std::int64_t shape0_local = x.size() / shape1;
686
687 int err;
688 std::int64_t shape0 = 0;
689 err = MPI_Allreduce(&shape0_local, &shape0, 1, MPI_INT64_T, MPI_SUM, comm0);
690 dolfinx::MPI::check_error(comm0, err);
691
692 std::int64_t rank_offset = -1;
693 if (comm1 != MPI_COMM_NULL)
694 {
695 rank_offset = 0;
696 err = MPI_Exscan(&shape0_local, &rank_offset, 1,
697 dolfinx::MPI::mpi_t<std::int64_t>, MPI_SUM, comm1);
698 dolfinx::MPI::check_error(comm1, err);
699 }
700 else
701 {
702 rank_offset = -1;
703 if (!x.empty())
704 throw std::runtime_error("Non-empty data on null MPI communicator");
705 }
706
707 return distribute_from_postoffice(comm0, indices, x, {shape0, shape1},
708 rank_offset);
709}
710//---------------------------------------------------------------------------
711
712} // namespace dolfinx::MPI
Comm(MPI_Comm comm, bool duplicate=true)
Duplicate communicator and wrap duplicate.
Definition MPI.cpp:12
~Comm()
Destructor (frees wrapped communicator)
Definition MPI.cpp:36
MPI_Comm comm() const noexcept
Return the underlying MPI_Comm object.
Definition MPI.cpp:62
Timer for measuring and logging elapsed time durations.
Definition Timer.h:41
MPI support functionality.
Definition MPI.h:31
MPI_Datatype mpi_t
Retrieves the MPI data type associated to the provided type.
Definition MPI.h:281
std::vector< typename std::remove_reference_t< typename U::value_type > > distribute_data(MPI_Comm comm0, std::span< const std::int64_t > indices, MPI_Comm comm1, const U &x, int shape1)
Distribute rows of a rectangular data array to ranks where they are required (scalable version).
Definition MPI.h:680
std::pair< std::vector< std::int32_t >, std::vector< typename std::remove_reference_t< typename U::value_type > > > distribute_to_postoffice(MPI_Comm comm, const U &x, std::array< std::int64_t, 2 > shape, std::int64_t rank_offset)
Distribute row data to 'post office' ranks.
Definition MPI.h:314
std::vector< int > compute_graph_edges_nbx(MPI_Comm comm, std::span< const int > edges, int tag=static_cast< int >(tag::consensus_nbx))
Determine incoming graph edges using the NBX consensus algorithm.
Definition MPI.cpp:162
constexpr int index_owner(int size, std::size_t index, std::size_t N)
Return which rank owns index in global range [0, N - 1] (inverse of MPI::local_range).
Definition MPI.h:114
std::vector< typename std::remove_reference_t< typename U::value_type > > distribute_from_postoffice(MPI_Comm comm, std::span< const std::int64_t > indices, const U &x, std::array< std::int64_t, 2 > shape, std::int64_t rank_offset)
Distribute rows of a rectangular data array from post office ranks to ranks where they are required.
Definition MPI.h:465
std::vector< int > compute_graph_edges_pcx(MPI_Comm comm, std::span< const int > edges)
Determine incoming graph edges using the PCX consensus algorithm.
Definition MPI.cpp:95
void check_error(MPI_Comm comm, int code)
Check MPI error code. If the error code is not equal to MPI_SUCCESS, then std::abort is called.
Definition MPI.cpp:80
int size(MPI_Comm comm)
Definition MPI.cpp:72
int rank(MPI_Comm comm)
Return process rank for the communicator.
Definition MPI.cpp:64
constexpr std::array< std::int64_t, 2 > local_range(int rank, std::int64_t N, int size)
Return local range for the calling process, partitioning the global [0, N - 1] range across all ranks...
Definition MPI.h:90
tag
MPI communication tags.
Definition MPI.h:34
Definition MPI.h:269
MPI Type.
Definition MPI.h:276