18#include <deal.II/base/mpi_noncontiguous_partitioner.templates.h>
40 const std::vector<types::global_dof_index> &indices_has,
41 const std::vector<types::global_dof_index> &indices_want,
49 std::pair<unsigned int, unsigned int>
106 std::vector<unsigned int> owning_ranks_of_ghosts(
111 process(indexset_has,
114 owning_ranks_of_ghosts,
118 std::pair<types::global_dof_index, types::global_dof_index>,
121 consensus_algorithm.run();
125 std::map<unsigned int, std::vector<types::global_dof_index>> recv_map;
127 for (
const auto &owner : owning_ranks_of_ghosts)
128 recv_map[owner] = std::vector<types::global_dof_index>();
132 recv_map[owning_ranks_of_ghosts[i]].push_back(i);
135 for (
const auto &target_with_indexset : recv_map)
137 recv_ranks.push_back(target_with_indexset.first);
139 for (
const auto cell_index : target_with_indexset.second)
150 for (
const auto &target_with_indexset : targets_with_indexset)
152 send_ranks.push_back(target_with_indexset.first);
154 for (
const auto cell_index : target_with_indexset.second)
166 const std::vector<types::global_dof_index> &indices_has,
167 const std::vector<types::global_dof_index> &indices_want,
172 std::vector<types::global_dof_index> indices_has_clean;
173 indices_has_clean.reserve(indices_has.size());
175 for (
const auto i : indices_has)
177 indices_has_clean.push_back(i);
179 std::vector<types::global_dof_index> indices_want_clean;
180 indices_want_clean.reserve(indices_want.size());
182 for (
const auto i : indices_want)
184 indices_want_clean.push_back(i);
188 indices_has_clean.empty() ?
190 (*std::max_element(indices_has_clean.begin(),
191 indices_has_clean.end()) +
195 indices_want_clean.empty() ?
197 (*std::max_element(indices_want_clean.begin(),
198 indices_want_clean.end()) +
207 index_set_has.
add_indices(indices_has_clean.begin(),
208 indices_has_clean.end());
211 index_set_want.
add_indices(indices_want_clean.begin(),
212 indices_want_clean.end());
220 std::vector<types::global_dof_index> temp_map_send(
228 i = temp_map_send[i];
232 std::vector<types::global_dof_index> temp_map_recv(
240 i = temp_map_recv[i];
246#include "mpi_noncontiguous_partitioner.inst"
size_type index_within_set(const size_type global_index) const
size_type n_elements() const
void add_indices(const ForwardIterator &begin, const ForwardIterator &end)
std::vector< MPI_Request > requests
std::vector< uint8_t > buffers
types::global_dof_index memory_consumption()
const MPI_Comm & get_mpi_communicator() const override
unsigned int temporary_storage_size() const
std::vector< types::global_dof_index > recv_ptr
std::vector< types::global_dof_index > send_ptr
std::vector< unsigned int > send_ranks
std::vector< types::global_dof_index > recv_indices
std::pair< unsigned int, unsigned int > n_targets() const
void reinit(const IndexSet &indexset_locally_owned, const IndexSet &indexset_ghost, const MPI_Comm &communicator) override
std::vector< types::global_dof_index > send_indices
std::vector< unsigned int > recv_ranks
NoncontiguousPartitioner()=default
std::map< unsigned int, IndexSet > get_requesters()
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
T max(const T &t, const MPI_Comm &mpi_communicator)
const types::global_dof_index invalid_dof_index
::VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)