From 2eef9e012c3577ee96d580791c086dfc3fcac48a Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Tue, 1 Aug 2023 09:32:28 +0200 Subject: [PATCH] Add NoncontiguousPartitioner::import_from_ghosted_array() --- .../base/mpi_noncontiguous_partitioner.h | 57 ++++++ .../mpi_noncontiguous_partitioner.templates.h | 181 ++++++++++++++++++ .../mpi_noncontiguous_partitioner.inst.in | 6 + .../base/mpi_noncontiguous_partitioner_04.cc | 91 +++++++++ ...ncontiguous_partitioner_04.mpirun=2.output | 7 + 5 files changed, 342 insertions(+) create mode 100644 tests/base/mpi_noncontiguous_partitioner_04.cc create mode 100644 tests/base/mpi_noncontiguous_partitioner_04.mpirun=2.output diff --git a/include/deal.II/base/mpi_noncontiguous_partitioner.h b/include/deal.II/base/mpi_noncontiguous_partitioner.h index 5a8c2657a7..2d579031ae 100644 --- a/include/deal.II/base/mpi_noncontiguous_partitioner.h +++ b/include/deal.II/base/mpi_noncontiguous_partitioner.h @@ -166,6 +166,63 @@ namespace Utilities const ArrayView & ghost_array, std::vector & requests) const; + /** + * Similar to the above functions but for importing vector entries + * from @p ghost_array to @p locally_owned_storage. + * + * @note In contrast to the functions in + * Utilities::MPI::Partitioner, this function expects that + * locally_owned_storage is empty. + */ + template + void + import_from_ghosted_array( + const VectorOperation::values vector_operation, + const ArrayView & ghost_array, + const ArrayView & locally_owned_storage) const; + + /** + * Similar to the above function with the difference that + * users can provide temporaty arrays. This function calls + * import_from_ghosted_array_start() and + * import_from_ghosted_array_finish() in sequence. + */ + template + void + import_from_ghosted_array(const VectorOperation::values vector_operation, + const unsigned int communication_channel, + const ArrayView & ghost_array, + const ArrayView & temporary_storage, + const ArrayView & locally_owned_storage, + std::vector &requests) const; + + /** + * Start update for importig values: Data is packed, non-blocking send + * and receives are started. + */ + template + void + import_from_ghosted_array_start( + const VectorOperation::values vector_operation, + const unsigned int communication_channel, + const ArrayView & ghost_array, + const ArrayView & temporary_storage, + std::vector & requests) const; + + /** + * Finish update for importing values. The method waits until all data has + * been sent and received. Once data from any process is received it is + * processed and placed at the right position of the vector + * @p locally_owned_storage. + */ + template + void + import_from_ghosted_array_finish( + const VectorOperation::values vector_operation, + const ArrayView &temporary_storage, + const ArrayView & locally_owned_storage, + std::vector & requests) const; + /** * Returns the number of processes this process sends data to and the * number of processes this process receives data from. diff --git a/include/deal.II/base/mpi_noncontiguous_partitioner.templates.h b/include/deal.II/base/mpi_noncontiguous_partitioner.templates.h index 9e83c2819c..d6e2484dc9 100644 --- a/include/deal.II/base/mpi_noncontiguous_partitioner.templates.h +++ b/include/deal.II/base/mpi_noncontiguous_partitioner.templates.h @@ -188,6 +188,187 @@ namespace Utilities #endif } + + + template + void + NoncontiguousPartitioner::import_from_ghosted_array( + const VectorOperation::values vector_operation, + const ArrayView & src, + const ArrayView & dst) const + { + // allocate internal memory since needed + if (requests.size() != send_ranks.size() + recv_ranks.size()) + requests.resize(send_ranks.size() + recv_ranks.size()); + + if (this->buffers.size() != send_ptr.back() * sizeof(Number)) + this->buffers.resize(this->temporary_storage_size() * sizeof(Number)); + + // perform actual exchange + this->template import_from_ghosted_array( + vector_operation, + 0, + src, + ArrayView(reinterpret_cast(this->buffers.data()), + send_ptr.back()), + dst, + this->requests); + } + + + + template + void + NoncontiguousPartitioner::import_from_ghosted_array( + const VectorOperation::values vector_operation, + const unsigned int communication_channel, + const ArrayView & ghost_array, + const ArrayView & temporary_storage, + const ArrayView & locally_owned_array, + std::vector & requests) const + { + this->template import_from_ghosted_array_start( + vector_operation, + communication_channel, + ghost_array, + temporary_storage, + requests); + this->template import_from_ghosted_array_finish( + vector_operation, temporary_storage, locally_owned_array, requests); + } + + + + template + void + NoncontiguousPartitioner::import_from_ghosted_array_start( + const VectorOperation::values vector_operation, + const unsigned int communication_channel, + const ArrayView & src, + const ArrayView & buffers, + std::vector & requests) const + { +#ifndef DEAL_II_WITH_MPI + (void)vector_operation; + (void)communication_channel; + (void)src; + (void)buffers; + (void)requests; + Assert(false, ExcNeedsMPI()); +#else + (void)vector_operation; // nothing to do here + + AssertDimension(requests.size(), recv_ranks.size() + send_ranks.size()); + + const auto tag = + communication_channel + + internal::Tags::noncontiguous_partitioner_update_ghost_values_start; + + AssertIndexRange( + tag, + internal::Tags::noncontiguous_partitioner_update_ghost_values_end + 1); + + // post recv + AssertIndexRange(send_ranks.size(), send_ptr.size()); + for (types::global_dof_index i = 0; i < send_ranks.size(); ++i) + { + const int ierr = + MPI_Irecv(buffers.data() + send_ptr[i], + send_ptr[i + 1] - send_ptr[i], + Utilities::MPI::mpi_type_id_for_type, + send_ranks[i], + tag, + communicator, + &requests[i + send_ranks.size()]); + AssertThrowMPI(ierr); + } + + // pack data and send away + for (types::global_dof_index i = 0; i < recv_ranks.size(); ++i) + { + AssertIndexRange(i + 1, recv_ptr.size()); + for (types::global_dof_index j = recv_ptr[i], c = 0; + j < recv_ptr[i + 1]; + j++) + buffers[recv_ptr[i] + c++] = src[recv_indices[j]]; + + // send data + Assert((recv_ptr[i] < buffers.size()) || + (recv_ptr[i] == buffers.size() && + recv_ptr[i + 1] == recv_ptr[i]), + ExcMessage("The input buffer doesn't contain enough entries")); + const int ierr = + MPI_Isend(buffers.data() + recv_ptr[i], + recv_ptr[i + 1] - recv_ptr[i], + Utilities::MPI::mpi_type_id_for_type, + recv_ranks[i], + tag, + communicator, + &requests[i]); + AssertThrowMPI(ierr); + } +#endif + } + + + + template + void + NoncontiguousPartitioner::import_from_ghosted_array_finish( + const VectorOperation::values vector_operation, + const ArrayView &buffers, + const ArrayView & dst, + std::vector & requests) const + { +#ifndef DEAL_II_WITH_MPI + (void)vector_operation; + (void)buffers; + (void)dst; + (void)requests; + Assert(false, ExcNeedsMPI()); +#else + Assert(vector_operation == VectorOperation::add, ExcNotImplemented()); + + AssertDimension(requests.size(), recv_ranks.size() + send_ranks.size()); + + Assert(std::accumulate(dst.begin(), dst.end(), 0) == 0, + ExcMessage("The destination vector has to be empty.")); + + // wait that all data packages have been received + // note: this for-loop cold be merged with the next for-loop, + // however, for this send_indices would be needed to stored + // rank by rank + for (types::global_dof_index proc = 0; proc < send_ranks.size(); ++proc) + { + int i; + MPI_Status status; + const auto ierr = MPI_Waitany(recv_ranks.size(), + requests.data() + send_ranks.size(), + &i, + &status); + AssertThrowMPI(ierr); + } + + // write data into destination vector + for (types::global_dof_index i = 0, k = 0; i < send_ranks.size(); ++i) + { + // collect data to be send + for (types::global_dof_index j = send_ptr[i]; j < send_ptr[i + 1]; + j++) + { + AssertIndexRange(k, send_indices.size()); + dst[send_indices[k]] += buffers[j]; + ++k; + } + } + + // wait that all data packages have been sent + const int ierr = + MPI_Waitall(send_ranks.size(), requests.data(), MPI_STATUSES_IGNORE); + AssertThrowMPI(ierr); +#endif + } + } // namespace MPI } // namespace Utilities diff --git a/source/base/mpi_noncontiguous_partitioner.inst.in b/source/base/mpi_noncontiguous_partitioner.inst.in index 08b73a779f..e54bd6f5e4 100644 --- a/source/base/mpi_noncontiguous_partitioner.inst.in +++ b/source/base/mpi_noncontiguous_partitioner.inst.in @@ -25,6 +25,12 @@ for (S : REAL_SCALARS) NoncontiguousPartitioner::export_to_ghosted_array( const ArrayView &src, const ArrayView & dst) const; + + template void + NoncontiguousPartitioner::import_from_ghosted_array( + const VectorOperation::values vector_operation, + const ArrayView & src, + const ArrayView & dst) const; \} \} } diff --git a/tests/base/mpi_noncontiguous_partitioner_04.cc b/tests/base/mpi_noncontiguous_partitioner_04.cc new file mode 100644 index 0000000000..e6ba6c3f67 --- /dev/null +++ b/tests/base/mpi_noncontiguous_partitioner_04.cc @@ -0,0 +1,91 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2023 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. +// +// --------------------------------------------------------------------- + + +// Test Utilities::MPI::NoncontiguousPartitioner::import_from_ghosted_array(). + +#include +#include + +#include "../tests.h" + + +void +test(const MPI_Comm comm) +{ + std::vector index_set_has; + std::vector index_set_want; + + if (Utilities::MPI::this_mpi_process(comm) == 0) + { + index_set_has.push_back(0); + index_set_has.push_back(1); + index_set_has.push_back(2); + + index_set_want.push_back(0); + index_set_want.push_back(1); + index_set_want.push_back(2); + index_set_want.push_back(3); + index_set_want.push_back(5); + } + else + { + index_set_has.push_back(3); + index_set_has.push_back(4); + index_set_has.push_back(5); + + index_set_want.push_back(0); + index_set_want.push_back(2); + index_set_want.push_back(3); + index_set_want.push_back(4); + index_set_want.push_back(5); + } + + Utilities::MPI::NoncontiguousPartitioner vector(index_set_has, + index_set_want, + comm); + + AlignedVector src(index_set_want.size()); + AlignedVector dst(index_set_has.size()); + + for (unsigned int i = 0; i < src.size(); ++i) + src[i] = i + Utilities::MPI::this_mpi_process(comm) * 100 + 1; + + vector.import_from_ghosted_array(VectorOperation::add, + ArrayView(src.data(), src.size()), + ArrayView(dst.data(), dst.size())); + + for (size_t i = 0; i < src.size(); ++i) + deallog << static_cast(src[i]) << ' '; + deallog << std::endl; + for (size_t i = 0; i < dst.size(); ++i) + deallog << static_cast(dst[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; + + { + deallog.push("all"); + test(comm); + deallog.pop(); + } +} diff --git a/tests/base/mpi_noncontiguous_partitioner_04.mpirun=2.output b/tests/base/mpi_noncontiguous_partitioner_04.mpirun=2.output new file mode 100644 index 0000000000..44b1f20143 --- /dev/null +++ b/tests/base/mpi_noncontiguous_partitioner_04.mpirun=2.output @@ -0,0 +1,7 @@ + +DEAL:0:all::1 2 3 4 5 +DEAL:0:all::102 2 105 + +DEAL:1:all::101 102 103 104 105 +DEAL:1:all::107 104 110 + -- 2.39.5