From 37948dfdfea47e9e5534f987701983be099aa794 Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Sat, 4 Apr 2020 20:18:24 +0200 Subject: [PATCH] Add Utilities::MPI::compute_union --- doc/news/changes/minor/20200404PeterMunch | 4 ++ include/deal.II/base/mpi.h | 11 ++++ include/deal.II/base/mpi.templates.h | 71 +++++++++++++++++++++++ include/deal.II/base/mpi_tags.h | 3 + source/base/mpi.cc | 4 ++ tests/mpi/union_01.cc | 51 ++++++++++++++++ tests/mpi/union_01.mpirun=2.output | 5 ++ 7 files changed, 149 insertions(+) create mode 100644 doc/news/changes/minor/20200404PeterMunch create mode 100644 tests/mpi/union_01.cc create mode 100644 tests/mpi/union_01.mpirun=2.output diff --git a/doc/news/changes/minor/20200404PeterMunch b/doc/news/changes/minor/20200404PeterMunch new file mode 100644 index 0000000000..20b83cd62a --- /dev/null +++ b/doc/news/changes/minor/20200404PeterMunch @@ -0,0 +1,4 @@ +New: Add new function Utilities::MPI::compute_union(), which computes the union of given input +vectors of all processes in a given MPI communicator. +
+(Peter Munch, 2020/04/04) diff --git a/include/deal.II/base/mpi.h b/include/deal.II/base/mpi.h index 6431e56ea1..8c118a978b 100644 --- a/include/deal.II/base/mpi.h +++ b/include/deal.II/base/mpi.h @@ -1024,6 +1024,17 @@ namespace Utilities const IndexSet &indices_to_look_up, const MPI_Comm &comm); + /** + * Compute the union of the input vectors @p vec of all processes in the + * MPI communicator @p comm. + * + * @note This is a collective operation. The result will available on all + * processes. + */ + template + std::vector + compute_union(const std::vector &vec, const MPI_Comm &comm); + #ifndef DOXYGEN // declaration for an internal function that lives in mpi.templates.h namespace internal diff --git a/include/deal.II/base/mpi.templates.h b/include/deal.II/base/mpi.templates.h index 4c0ea2c821..d99d6c4841 100644 --- a/include/deal.II/base/mpi.templates.h +++ b/include/deal.II/base/mpi.templates.h @@ -28,6 +28,7 @@ #include #include +#include #include DEAL_II_NAMESPACE_OPEN @@ -413,6 +414,76 @@ namespace Utilities internal::all_reduce(MPI_MIN, values, mpi_communicator, minima); } + + + template + std::vector + compute_union(const std::vector &vec, const MPI_Comm &comm) + { +#ifdef DEAL_II_WITH_MPI + // 1) collect vector entries and create union + std::vector result = vec; + + if (this_mpi_process(comm) == 0) + { + for (unsigned int i = 1; i < n_mpi_processes(comm); i++) + { + MPI_Status status; + const auto ierr_1 = MPI_Probe(MPI_ANY_SOURCE, + internal::Tags::compute_union, + comm, + &status); + AssertThrowMPI(ierr_1); + + int amount; + const auto ierr_2 = + MPI_Get_count(&status, + internal::mpi_type_id(vec.data()), + &amount); + AssertThrowMPI(ierr_2); + + const unsigned int old_size = result.size(); + result.resize(old_size + amount); + + const auto ierr_3 = MPI_Recv(result.data() + old_size, + amount, + MPI_INT, + status.MPI_SOURCE, + status.MPI_TAG, + comm, + MPI_STATUS_IGNORE); + AssertThrowMPI(ierr_3); + + std::sort(result.begin(), result.end()); + result.erase(std::unique(result.begin(), result.end()), + result.end()); + } + } + else + { + const auto ierr = MPI_Send(vec.data(), + vec.size(), + internal::mpi_type_id(vec.data()), + 0, + internal::Tags::compute_union, + comm); + AssertThrowMPI(ierr); + } + + // 2) broadcast result + int size = result.size(); + MPI_Bcast(&size, 1, MPI_INT, 0, comm); + result.resize(size); + MPI_Bcast( + result.data(), size, internal::mpi_type_id(vec.data()), 0, comm); + + return result; +#else + (void)comm; + return vec; +#endif + } + } // end of namespace MPI } // end of namespace Utilities diff --git a/include/deal.II/base/mpi_tags.h b/include/deal.II/base/mpi_tags.h index 548df86b96..4c5d2f5333 100644 --- a/include/deal.II/base/mpi_tags.h +++ b/include/deal.II/base/mpi_tags.h @@ -132,6 +132,9 @@ namespace Utilities /// NoncontiguousPartitioner::update_values noncontiguous_partitioner_update_values, + // Utilities::MPI::compute_union + compute_union, + }; } // namespace Tags } // namespace internal diff --git a/source/base/mpi.cc b/source/base/mpi.cc index d0a62bc9bb..d6784fe361 100644 --- a/source/base/mpi.cc +++ b/source/base/mpi.cc @@ -1126,6 +1126,10 @@ namespace Utilities locked = false; } + + template std::vector + compute_union(const std::vector &vec, const MPI_Comm &comm); + #include "mpi.inst" } // end of namespace MPI } // end of namespace Utilities diff --git a/tests/mpi/union_01.cc b/tests/mpi/union_01.cc new file mode 100644 index 0000000000..884e113499 --- /dev/null +++ b/tests/mpi/union_01.cc @@ -0,0 +1,51 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + +// Create a distributed triangulation with multigrid levels and copy it. + +#include + +#include "../tests.h" + +using namespace dealii; + +void +test(const MPI_Comm comm) +{ + std::vector vector; + + if (Utilities::MPI::this_mpi_process(comm) == 0) + vector = std::vector{2, 5, 7}; + else + vector = std::vector{5, 8, 9, 10}; + + const auto result = Utilities::MPI::compute_union(vector, comm); + + for (auto i : result) + deallog << i << " "; + deallog << std::endl; +} + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll all; + + const MPI_Comm comm = MPI_COMM_WORLD; + + test(comm); +} diff --git a/tests/mpi/union_01.mpirun=2.output b/tests/mpi/union_01.mpirun=2.output new file mode 100644 index 0000000000..ba6eac9eb4 --- /dev/null +++ b/tests/mpi/union_01.mpirun=2.output @@ -0,0 +1,5 @@ + +DEAL:0::2 5 7 8 9 10 + +DEAL:1::2 5 7 8 9 10 + -- 2.39.5