From f1c1fcdbc51de5ac1e28ca7add7056071ca56959 Mon Sep 17 00:00:00 2001 From: Eldar Date: Thu, 12 Jul 2018 08:19:58 -0500 Subject: [PATCH] ReduceScatter for compute p2p comm pattern --- source/base/mpi.cc | 48 ++++++++++++++++++++++++++ tests/mpi/point_to_point_pattern_01.cc | 11 +++--- 2 files changed, 55 insertions(+), 4 deletions(-) diff --git a/source/base/mpi.cc b/source/base/mpi.cc index 707a19fa8c..81d3b3aaf5 100644 --- a/source/base/mpi.cc +++ b/source/base/mpi.cc @@ -197,7 +197,54 @@ namespace Utilities "There is no point in communicating with ourselves.")); } +# if DEAL_II_MPI_VERSION_GTE(2, 2) + // Get a binary vector from the `destinations` + std::vector dest_vector(n_procs); + for (const auto &el : destinations) + dest_vector[el] = 1; + + // Find how many processes will send to this one + // by reducing with sum and then scattering the + // results over all processes + unsigned int n_recv_from; + const int ierr = MPI_Reduce_scatter_block( + &dest_vector[0], &n_recv_from, 1, MPI_UNSIGNED, MPI_SUM, mpi_comm); + AssertThrowMPI(ierr); + + // Send myid to every process in `destinations` vector... + std::vector send_requests(destinations.size()); + for (const auto &el : destinations) + MPI_Isend(&myid, + 1, + MPI_UNSIGNED, + el, + 32766, + mpi_comm, + &send_requests[&el - &destinations[0]]); + + // if no one to receive from, return an empty vector + if (n_recv_from == 0) + return std::vector(); + + // ...otherwise receive `n_recv_from` times from the processes + // who communicate with this one. Store the obtained id's + // in the resulting vector + std::vector origins(n_recv_from); + for (auto &el : origins) + MPI_Recv(&el, + 1, + MPI_UNSIGNED, + MPI_ANY_SOURCE, + 32766, + mpi_comm, + MPI_STATUS_IGNORE); + + MPI_Waitall(destinations.size(), + send_requests.data(), + MPI_STATUSES_IGNORE); + return origins; +# else // let all processors communicate the maximal number of destinations // they have const unsigned int max_n_destinations = @@ -242,6 +289,7 @@ namespace Utilities break; return origins; +# endif } diff --git a/tests/mpi/point_to_point_pattern_01.cc b/tests/mpi/point_to_point_pattern_01.cc index bfc40dc43a..caad9858ad 100644 --- a/tests/mpi/point_to_point_pattern_01.cc +++ b/tests/mpi/point_to_point_pattern_01.cc @@ -1,6 +1,6 @@ // --------------------------------------------------------------------- // -// Copyright (C) 2009 - 2017 by the deal.II authors +// Copyright (C) 2009 - 2018 by the deal.II authors // // This file is part of the deal.II library. // @@ -15,13 +15,15 @@ -// check if mpi is working +// Check if point-to-point communication +// works correctly with MPI_Reduce_scatter_block +// implementation #include -#include "../tests.h" +#include -//#include +#include "../tests.h" void @@ -74,6 +76,7 @@ test_mpi() std::vector origins = Utilities::MPI::compute_point_to_point_communication_pattern(MPI_COMM_WORLD, destinations); + std::sort(origins.begin(), origins.end()); if (myid == 0) { -- 2.39.5