From: Wolfgang Bangerth Date: Sat, 2 Nov 2024 22:13:13 +0000 (-0600) Subject: Use Peter's new function in a number of places. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=749c8d7085122b68479612b666dc5df4a27cbdc0;p=dealii.git Use Peter's new function in a number of places. --- diff --git a/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h b/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h index 87b04c4f31..9d77b486df 100644 --- a/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h +++ b/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h @@ -851,30 +851,17 @@ namespace internal is_dst_remote_potentially_relevant.nth_index_in_set(i)); } - // determine owner of remote cells - std::vector is_dst_remote_owners( - is_dst_remote.n_elements()); - - Utilities::MPI::internal::ComputeIndexOwner::ConsensusAlgorithmsPayload - process(is_dst_locally_owned, - is_dst_remote, - communicator, - is_dst_remote_owners, - true); + // determine owners of remote cells + const auto [is_dst_remote_owners, targets_with_indexset] = + Utilities::MPI::compute_index_owner_and_requesters(is_dst_locally_owned, + is_dst_remote, + communicator); - Utilities::MPI::ConsensusAlgorithms::Selector< - std::vector< - std::pair>, - std::vector> - consensus_algorithm; - consensus_algorithm.run(process, communicator); this->is_dst_locally_owned = is_dst_locally_owned; this->is_dst_remote = is_dst_remote; this->is_src_locally_owned = is_src_locally_owned; - const auto targets_with_indexset = process.get_requesters(); - #ifndef DEAL_II_WITH_MPI Assert(targets_with_indexset.empty(), ExcInternalError()); #else diff --git a/source/base/mpi_noncontiguous_partitioner.cc b/source/base/mpi_noncontiguous_partitioner.cc index 178f224964..28ccc8e8a8 100644 --- a/source/base/mpi_noncontiguous_partitioner.cc +++ b/source/base/mpi_noncontiguous_partitioner.cc @@ -104,23 +104,10 @@ namespace Utilities requests.clear(); // set up communication pattern - std::vector owning_ranks_of_ghosts( - indexset_want.n_elements()); - - // set up dictionary - Utilities::MPI::internal::ComputeIndexOwner::ConsensusAlgorithmsPayload - process(indexset_has, - indexset_want, - communicator, - owning_ranks_of_ghosts, - true); - - Utilities::MPI::ConsensusAlgorithms::Selector< - std::vector< - std::pair>, - std::vector> - consensus_algorithm; - consensus_algorithm.run(process, communicator); + const auto [owning_ranks_of_ghosts, targets_with_indexset] = + Utilities::MPI::compute_index_owner_and_requesters(indexset_has, + indexset_want, + communicator); // set up map of processes from where this rank will receive values { @@ -146,8 +133,6 @@ namespace Utilities } { - const auto targets_with_indexset = process.get_requesters(); - send_ptr.push_back(recv_ptr.back()); for (const auto &target_with_indexset : targets_with_indexset) { diff --git a/source/base/partitioner.cc b/source/base/partitioner.cc index 8e537fcde3..dcf4257236 100644 --- a/source/base/partitioner.cc +++ b/source/base/partitioner.cc @@ -262,26 +262,10 @@ namespace Utilities local_range_data.second = my_shift + old_locally_owned_size; } - std::vector owning_ranks_of_ghosts( - ghost_indices_data.n_elements()); - - // set up dictionary - internal::ComputeIndexOwner::ConsensusAlgorithmsPayload process( - locally_owned_range_data, - ghost_indices_data, - communicator, - owning_ranks_of_ghosts, - /* track origins of ghosts*/ true); - - // read dictionary by communicating with the process who owns the index - // in the static partition (i.e. in the dictionary). This process - // returns the actual owner of the index. - ConsensusAlgorithms::Selector< - std::vector< - std::pair>, - std::vector> - consensus_algorithm; - consensus_algorithm.run(process, communicator); + const auto [owning_ranks_of_ghosts, import_data] = + Utilities::MPI::compute_index_owner_and_requesters( + locally_owned_range_data, ghost_indices_data, communicator); + { ghost_targets_data = {}; @@ -303,9 +287,6 @@ namespace Utilities } } - // find how much the individual processes that want import from me - std::map import_data = process.get_requesters(); - // count import requests and set up the compressed indices n_import_indices_data = 0; import_targets_data = {}; diff --git a/source/matrix_free/vector_data_exchange.cc b/source/matrix_free/vector_data_exchange.cc index cf7057ec36..72216f5943 100644 --- a/source/matrix_free/vector_data_exchange.cc +++ b/source/matrix_free/vector_data_exchange.cc @@ -428,22 +428,10 @@ namespace internal Utilities::MPI::mpi_processes_within_communicator(comm, comm_sm); // determine owners of ghost indices and determine requesters - std::vector owning_ranks_of_ghosts( - is_locally_ghost.n_elements()); - - Utilities::MPI::internal::ComputeIndexOwner::ConsensusAlgorithmsPayload - process(is_locally_owned, - is_locally_ghost, - comm, - owning_ranks_of_ghosts, - /*track_index_requests = */ true); - - Utilities::MPI::ConsensusAlgorithms::Selector< - std::vector< - std::pair>, - std::vector> - consensus_algorithm; - consensus_algorithm.run(process, comm); + const auto [owning_ranks_of_ghosts, rank_to_global_indices] = + Utilities::MPI::compute_index_owner_and_requesters(is_locally_owned, + is_locally_ghost, + comm); // decompress ghost_indices_within_larger_ghost_set for simpler // data access during setup @@ -525,8 +513,6 @@ namespace internal // process requesters { - const auto rank_to_global_indices = process.get_requesters(); - for (const auto &rank_and_global_indices : rank_to_global_indices) { const auto sm_ranks_ptr =