From: Peter Munch Date: Tue, 17 Sep 2024 06:15:45 +0000 (+0200) Subject: NoncontiguousPartitioner: introduce n_components X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=2732e7632a6d55e74004792ee5e4ab79e63fed50;p=dealii.git NoncontiguousPartitioner: introduce n_components --- diff --git a/doc/news/changes/minor/20240925Munch b/doc/news/changes/minor/20240925Munch new file mode 100644 index 0000000000..3588a3a6d5 --- /dev/null +++ b/doc/news/changes/minor/20240925Munch @@ -0,0 +1,4 @@ +New: Utilities::MPI::NoncontiguousPartitioner::export_to_ghosted_array() +can now handle multiple components. +
+(Peter Munch, 2024/09/25) diff --git a/include/deal.II/base/mpi_noncontiguous_partitioner.h b/include/deal.II/base/mpi_noncontiguous_partitioner.h index a4bf8caa56..86cb314678 100644 --- a/include/deal.II/base/mpi_noncontiguous_partitioner.h +++ b/include/deal.II/base/mpi_noncontiguous_partitioner.h @@ -79,6 +79,14 @@ namespace Utilities * Fill the vector @p ghost_array according to the precomputed communication * pattern with values from @p locally_owned_array. * + * In the default case, only one object is communicated per entry + * (`n_components_templated == 1'). If you want to communicate more + * entries, you can increase the value of @p n_components_templated in the + * case that you know the size at compile time. If you want to set the + * size during runtime, you can set @p n_components. However, + * @p n_components_templated has to be set to `0` in this case. Either + * @p n_components_templated or @p n_components can be set. + * * @pre The vectors only have to provide a method begin(), which allows * to access their raw data. * @@ -91,11 +99,12 @@ namespace Utilities * functions separately and hereby overlap communication and * computation. */ - template + template void export_to_ghosted_array( const ArrayView &locally_owned_array, - const ArrayView &ghost_array) const; + const ArrayView &ghost_array, + const unsigned int n_components = 0) const; /** * Same as above but with an interface similar to @@ -111,14 +120,15 @@ namespace Utilities * @note Any value less than 10 is a valid value of * @p communication_channel. */ - template + template void export_to_ghosted_array( const unsigned int communication_channel, const ArrayView &locally_owned_array, const ArrayView &temporary_storage, const ArrayView &ghost_array, - std::vector &requests) const; + std::vector &requests, + const unsigned int n_components = 0) const; /** * Start update: Data is packed, non-blocking send and receives @@ -136,13 +146,14 @@ namespace Utilities * @note Any value less than 10 is a valid value of * @p communication_channel. */ - template + template void export_to_ghosted_array_start( const unsigned int communication_channel, const ArrayView &locally_owned_array, const ArrayView &temporary_storage, - std::vector &requests) const; + std::vector &requests, + const unsigned int n_components = 0) const; /** * Finish update. The method waits until all data has been sent and @@ -158,12 +169,13 @@ namespace Utilities * @pre The required size of the vectors are the same as in the functions * above. */ - template + template void export_to_ghosted_array_finish( const ArrayView &temporary_storage, const ArrayView &ghost_array, - std::vector &requests) const; + std::vector &requests, + const unsigned int n_components = 0) const; /** * Similar to the above functions but for importing vector entries diff --git a/include/deal.II/base/mpi_noncontiguous_partitioner.templates.h b/include/deal.II/base/mpi_noncontiguous_partitioner.templates.h index 2f04bec8b0..b9329c47b2 100644 --- a/include/deal.II/base/mpi_noncontiguous_partitioner.templates.h +++ b/include/deal.II/base/mpi_noncontiguous_partitioner.templates.h @@ -31,68 +31,89 @@ namespace Utilities { namespace MPI { - template + template void NoncontiguousPartitioner::export_to_ghosted_array( const ArrayView &src, - const ArrayView &dst) const + const ArrayView &dst, + const unsigned int n_components) const { + Assert((n_components_templated != 0) != (n_components != 0), + ExcNotImplemented()); + + const unsigned int n_components_to_be_used = + (n_components_templated != 0) ? n_components_templated : n_components; + // 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)); + if (this->buffers.size() != + send_ptr.back() * sizeof(Number) * n_components_to_be_used) + this->buffers.resize(this->temporary_storage_size() * sizeof(Number) * + n_components_to_be_used); // perform actual exchange - this->template export_to_ghosted_array( + this->template export_to_ghosted_array( 0, src, ArrayView(reinterpret_cast(this->buffers.data()), - send_ptr.back()), + send_ptr.back() * n_components_to_be_used), dst, - this->requests); + this->requests, + n_components); } - template + template void NoncontiguousPartitioner::export_to_ghosted_array( const unsigned int communication_channel, const ArrayView &locally_owned_array, const ArrayView &temporary_storage, const ArrayView &ghost_array, - std::vector &requests) const + std::vector &requests, + const unsigned int n_components) const { - this->template export_to_ghosted_array_start( + this->template export_to_ghosted_array_start( communication_channel, locally_owned_array, temporary_storage, - requests); - this->template export_to_ghosted_array_finish(temporary_storage, - ghost_array, - requests); + requests, + n_components); + this->template export_to_ghosted_array_finish( + temporary_storage, ghost_array, requests, n_components); } - template + template void NoncontiguousPartitioner::export_to_ghosted_array_start( const unsigned int communication_channel, const ArrayView &src, const ArrayView &buffers, - std::vector &requests) const + std::vector &requests, + const unsigned int n_components) const { #ifndef DEAL_II_WITH_MPI (void)communication_channel; (void)src; (void)buffers; (void)requests; + (void)n_components; Assert(false, ExcNeedsMPI()); #else AssertDimension(requests.size(), recv_ranks.size() + send_ranks.size()); + Assert((n_components_templated != 0) != (n_components != 0), + ExcNotImplemented()); + + const unsigned int n_components_to_be_used = + (n_components_templated != 0) ? n_components_templated : n_components; + const auto tag = communication_channel + internal::Tags::noncontiguous_partitioner_update_ghost_values_start; @@ -106,8 +127,8 @@ namespace Utilities for (types::global_dof_index i = 0; i < recv_ranks.size(); ++i) { const int ierr = - MPI_Irecv(buffers.data() + recv_ptr[i], - recv_ptr[i + 1] - recv_ptr[i], + MPI_Irecv(buffers.data() + recv_ptr[i] * n_components_to_be_used, + (recv_ptr[i + 1] - recv_ptr[i]) * n_components_to_be_used, Utilities::MPI::mpi_type_id_for_type, recv_ranks[i], tag, @@ -122,21 +143,24 @@ namespace Utilities { // collect data to be send for (types::global_dof_index j = send_ptr[i]; j < send_ptr[i + 1]; - j++) + ++j, ++k) { AssertIndexRange(k, send_indices.size()); - buffers[j] = src[send_indices[k]]; - ++k; + for (unsigned int comp = 0; comp < n_components_to_be_used; + ++comp) + buffers[j * n_components_to_be_used + comp] = + src[send_indices[k] * n_components_to_be_used + comp]; } // send data - Assert((send_ptr[i] < buffers.size()) || - (send_ptr[i] == buffers.size() && - send_ptr[i + 1] == send_ptr[i]), + Assert((send_ptr[i] * n_components_to_be_used < buffers.size()) || + (send_ptr[i] * n_components_to_be_used == buffers.size() && + send_ptr[i + 1] * n_components_to_be_used == + send_ptr[i] * n_components_to_be_used), ExcMessage("The input buffer doesn't contain enough entries")); const int ierr = - MPI_Isend(buffers.data() + send_ptr[i], - send_ptr[i + 1] - send_ptr[i], + MPI_Isend(buffers.data() + send_ptr[i] * n_components_to_be_used, + (send_ptr[i + 1] - send_ptr[i]) * n_components_to_be_used, Utilities::MPI::mpi_type_id_for_type, send_ranks[i], tag, @@ -149,19 +173,28 @@ namespace Utilities - template + template void NoncontiguousPartitioner::export_to_ghosted_array_finish( const ArrayView &buffers, const ArrayView &dst, - std::vector &requests) const + std::vector &requests, + const unsigned int n_components) const { #ifndef DEAL_II_WITH_MPI (void)buffers; (void)dst; (void)requests; + (void)n_components; Assert(false, ExcNeedsMPI()); #else + + Assert((n_components_templated != 0) != (n_components != 0), + ExcNotImplemented()); + + const unsigned int n_components_to_be_used = + (n_components_templated != 0) ? n_components_templated : n_components; + // receive all data packages and copy data from buffers for (types::global_dof_index proc = 0; proc < recv_ranks.size(); ++proc) { @@ -176,8 +209,10 @@ namespace Utilities AssertIndexRange(i + 1, recv_ptr.size()); for (types::global_dof_index j = recv_ptr[i], c = 0; j < recv_ptr[i + 1]; - j++) - dst[recv_indices[j]] = buffers[recv_ptr[i] + c++]; + ++j, ++c) + for (unsigned int comp = 0; comp < n_components_to_be_used; ++comp) + dst[recv_indices[j] * n_components_to_be_used + comp] = + buffers[(recv_ptr[i] + c) * n_components_to_be_used + comp]; } // wait that all data packages have been sent diff --git a/source/base/mpi_noncontiguous_partitioner.inst.in b/source/base/mpi_noncontiguous_partitioner.inst.in index 8a90e0b67b..70187d731f 100644 --- a/source/base/mpi_noncontiguous_partitioner.inst.in +++ b/source/base/mpi_noncontiguous_partitioner.inst.in @@ -23,7 +23,8 @@ for (S : REAL_SCALARS) template void NoncontiguousPartitioner::export_to_ghosted_array( const ArrayView &src, - const ArrayView &dst) const; + const ArrayView &dst, + const unsigned int) const; template void NoncontiguousPartitioner::import_from_ghosted_array( diff --git a/tests/base/mpi_noncontiguous_partitioner_05.cc b/tests/base/mpi_noncontiguous_partitioner_05.cc new file mode 100644 index 0000000000..5b9846c3a4 --- /dev/null +++ b/tests/base/mpi_noncontiguous_partitioner_05.cc @@ -0,0 +1,91 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2024 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + + +// Test Utilities::MPI::NoncontiguousPartitioner for non-contiguous index space +// and multiple components. + +#include +#include +#include + +#include "../tests.h" + + +void +test(const MPI_Comm comm) +{ + const unsigned int n_components = 2; + + IndexSet index_set_has(4); + IndexSet index_set_want(4); + + if (Utilities::MPI::this_mpi_process(comm) == 0) + { + index_set_has.add_index(1); + index_set_want.add_index(2); + } + else + { + index_set_has.add_index(2); + index_set_want.add_index(1); + index_set_want.add_index(2); + } + + Utilities::MPI::NoncontiguousPartitioner vector(index_set_has, + index_set_want, + comm); + + AlignedVector src(n_components * index_set_has.n_elements()); + AlignedVector dst(n_components * index_set_want.n_elements()); + + src[0] = Utilities::MPI::this_mpi_process(comm) * 100 + 1; + src[1] = Utilities::MPI::this_mpi_process(comm) * 100 + 2; + + vector.export_to_ghosted_array( + 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; + + dst.fill(0.0); + vector.export_to_ghosted_array( + ArrayView(src.data(), src.size()), + ArrayView(dst.data(), dst.size()), + n_components); + + 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_05.mpirun=2.output b/tests/base/mpi_noncontiguous_partitioner_05.mpirun=2.output new file mode 100644 index 0000000000..196985e76c --- /dev/null +++ b/tests/base/mpi_noncontiguous_partitioner_05.mpirun=2.output @@ -0,0 +1,9 @@ + +DEAL:0:all::1 2 +DEAL:0:all::101 102 +DEAL:0:all::101 102 + +DEAL:1:all::101 102 +DEAL:1:all::1 2 101 102 +DEAL:1:all::1 2 101 102 +