15#include <range/v3/numeric.hpp>
16#include <range/v3/range/conversion.hpp>
17#include <range/v3/view/enumerate.hpp>
18#include <range/v3/view/indirect.hpp>
19#include <range/v3/view/map.hpp>
20#include <unordered_map>
34 std::vector<std::size_t>
const& local_bulk_node_ids_for_subdomain,
35 std::size_t
const subdomain_node_id)
38 local_bulk_node_ids_for_subdomain[subdomain_node_id]));
47 MPI_Comm_size(mpi_comm, &mpi_comm_size);
49 MPI_Comm_rank(mpi_comm, &mpi_comm_rank);
50 return {mpi_comm_rank, mpi_comm_size};
54 std::vector<Element*>
const& input_elements)
59 std::unordered_map<std::size_t, Node*> id_node_hash_map;
60 id_node_hash_map.reserve(
63 for (
auto& e : elements)
66 unsigned const n_nodes = e->getNumberOfNodes();
67 for (
unsigned i = 0; i < n_nodes; ++i)
69 Node const* n = e->getNode(i);
70 auto const it = id_node_hash_map.find(n->
getID());
71 if (it == id_node_hash_map.end())
73 auto new_node_in_map = id_node_hash_map[n->
getID()] =
75 e->setNode(i, new_node_in_map);
79 e->setNode(i, it->second);
84 std::map<std::size_t, Node*> nodes_map;
85 for (
auto const& n : id_node_hash_map)
87 nodes_map[n.first] = n.second;
92 nodes_map | ranges::views::values | ranges::to<std::vector>;
94 return std::make_pair(element_nodes, elements);
98 Mesh const* subdomain_mesh)
100 auto const subdomain_nodes = subdomain_mesh->
getNodes();
101 auto const local_bulk_node_ids_for_subdomain =
103 auto const subdomain_node_ids =
105 unsigned long const number_of_regular_nodes =
106 std::count_if(subdomain_node_ids.begin(), subdomain_node_ids.end(),
107 std::bind_front(isRegularNode, *bulk_mesh,
108 local_bulk_node_ids_for_subdomain));
110 DBUG(
"[{}] number of regular nodes: {}", subdomain_mesh->
getName(),
111 number_of_regular_nodes);
112 return number_of_regular_nodes;
115std::vector<std::size_t>
120 auto const subdomain_nodes = subdomain_mesh->
getNodes();
121 auto const local_bulk_node_ids_for_subdomain =
135 std::size_t
const number_of_regular_nodes =
139 MPI_Comm mpi_comm = MPI_COMM_WORLD;
143 std::vector<std::size_t> gathered_number_of_regular_nodes(mpi_comm_size);
144 MPI_Allgather(&number_of_regular_nodes, 1, MPI_UNSIGNED_LONG,
145 gathered_number_of_regular_nodes.data(), 1, MPI_UNSIGNED_LONG,
148 std::vector<std::size_t> numbers_of_regular_nodes_at_rank;
149 numbers_of_regular_nodes_at_rank.push_back(0);
150 std::partial_sum(begin(gathered_number_of_regular_nodes),
151 end(gathered_number_of_regular_nodes),
152 back_inserter(numbers_of_regular_nodes_at_rank));
155 std::vector<std::size_t> subdomain_global_node_ids;
156 subdomain_global_node_ids.reserve(subdomain_nodes.size());
157 auto const partition_offset =
158 numbers_of_regular_nodes_at_rank[mpi_comm_rank];
159 DBUG(
"[{}] partition offset: {}", subdomain_mesh->
getName(),
164 if (isRegularNode(*bulk_mesh, local_bulk_node_ids_for_subdomain,
id))
166 subdomain_global_node_ids.emplace_back(partition_offset +
id);
169 return subdomain_global_node_ids;
174 Mesh const* subdomain_mesh,
175 std::vector<std::size_t>
const& global_regular_base_node_ids)
178 MPI_Comm
const mpi_comm = MPI_COMM_WORLD;
182 auto const local_bulk_node_ids_for_subdomain =
186 std::vector<std::size_t> subdomain_node_id_to_bulk_node_id;
189 if (!isRegularNode(*bulk_mesh, local_bulk_node_ids_for_subdomain,
id))
191 auto const bulk_node_id = local_bulk_node_ids_for_subdomain[id];
194 subdomain_node_id_to_bulk_node_id.push_back(
202 auto const size = subdomain_node_id_to_bulk_node_id.size();
203 std::size_t global_number_of_subdomain_node_id_to_bulk_node_id = 0;
204 MPI_Allreduce(&size, &global_number_of_subdomain_node_id_to_bulk_node_id, 1,
205 MPI_UNSIGNED_LONG, MPI_SUM, mpi_comm);
207 DBUG(
"[{}] global_number_of_subdomain_node_id_to_bulk_node_id: '{}' ",
209 global_number_of_subdomain_node_id_to_bulk_node_id);
211 std::vector<int> numbers_of_ids_at_ranks(mpi_comm_size);
212 MPI_Allgather(&size, 1, MPI_INT, numbers_of_ids_at_ranks.data(), 1, MPI_INT,
214 std::vector<int> offsets;
215 offsets.push_back(0);
216 std::partial_sum(begin(numbers_of_ids_at_ranks),
217 end(numbers_of_ids_at_ranks), back_inserter(offsets));
218 std::vector<std::size_t> ghost_node_ids_of_all_ranks(
219 global_number_of_subdomain_node_id_to_bulk_node_id);
220 MPI_Allgatherv(subdomain_node_id_to_bulk_node_id.data(),
223 ghost_node_ids_of_all_ranks.data(),
224 numbers_of_ids_at_ranks.data(),
230 std::map<std::size_t, std::size_t> global_to_local_bulk_node_ids;
240 std::map<std::size_t, std::size_t> local_subdomain_node_ids;
243 local_subdomain_node_ids[local_bulk_node_ids_for_subdomain[id]] = id;
246 std::vector<std::size_t> local_subdomain_node_ids_of_all_ranks(
247 global_number_of_subdomain_node_id_to_bulk_node_id,
248 std::numeric_limits<std::size_t>::max());
250 for (
int rank = 0; rank < mpi_comm_size; ++rank)
252 if (rank == mpi_comm_rank)
256 for (
int i = offsets[rank]; i < offsets[rank + 1]; ++i)
258 auto it = global_to_local_bulk_node_ids.find(
259 ghost_node_ids_of_all_ranks[i]);
260 if (it == global_to_local_bulk_node_ids.end())
264 auto const local_bulk_node_id = it->second;
265 local_subdomain_node_ids_of_all_ranks[i] =
266 global_regular_base_node_ids
267 [local_subdomain_node_ids.find(local_bulk_node_id)->second];
269 "[{}] found global subdomain node id: '{}' for global bulk "
272 local_subdomain_node_ids_of_all_ranks[i],
273 ghost_node_ids_of_all_ranks[i]);
278 std::vector<std::size_t> computed_global_ids_for_subdomain_ghost_nodes(
279 global_number_of_subdomain_node_id_to_bulk_node_id);
281 local_subdomain_node_ids_of_all_ranks.data(),
282 computed_global_ids_for_subdomain_ghost_nodes
284 global_number_of_subdomain_node_id_to_bulk_node_id,
289 std::vector<std::size_t> global_ids_for_subdomain_ghost_nodes(
290 computed_global_ids_for_subdomain_ghost_nodes.begin() +
291 offsets[mpi_comm_rank],
292 computed_global_ids_for_subdomain_ghost_nodes.begin() +
293 offsets[mpi_comm_rank + 1]);
294 return global_ids_for_subdomain_ghost_nodes;
298 Mesh const* subdomain_mesh)
300 auto const number_of_regular_base_nodes =
303 MPI_Comm
const mpi_comm = MPI_COMM_WORLD;
306 std::vector<std::size_t> gathered_number_of_regular_base_nodes(
308 MPI_Allgather(&number_of_regular_base_nodes, 1, MPI_UNSIGNED_LONG,
309 gathered_number_of_regular_base_nodes.data(), 1,
310 MPI_UNSIGNED_LONG, mpi_comm);
312 std::vector<std::size_t> numbers_of_regular_base_nodes_at_rank;
313 numbers_of_regular_base_nodes_at_rank.push_back(0);
314 std::partial_sum(begin(gathered_number_of_regular_base_nodes),
315 end(gathered_number_of_regular_base_nodes),
316 back_inserter(numbers_of_regular_base_nodes_at_rank));
318 return numbers_of_regular_base_nodes_at_rank;
323 Mesh const* subdomain_mesh)
326 MPI_Comm
const mpi_comm = MPI_COMM_WORLD;
329 auto const number_of_regular_base_nodes =
332 std::vector<std::size_t> gathered_number_of_regular_higher_order_nodes(
334 auto const number_of_regular_higher_order_nodes =
336 MPI_Allgather(&number_of_regular_higher_order_nodes, 1, MPI_UNSIGNED_LONG,
337 gathered_number_of_regular_higher_order_nodes.data(), 1,
338 MPI_UNSIGNED_LONG, mpi_comm);
340 std::vector<std::size_t> numbers_of_regular_higher_order_nodes_at_rank;
341 numbers_of_regular_higher_order_nodes_at_rank.push_back(0);
343 begin(gathered_number_of_regular_higher_order_nodes),
344 end(gathered_number_of_regular_higher_order_nodes),
345 back_inserter(numbers_of_regular_higher_order_nodes_at_rank));
347 return numbers_of_regular_higher_order_nodes_at_rank;
353 auto global_regular_base_node_ids =
356 auto const local_ids_for_subdomain_ghost_nodes =
358 bulk_mesh, subdomain_mesh, global_regular_base_node_ids);
372 std::vector<std::size_t> global_node_ids(
373 std::move(global_regular_base_node_ids));
374 global_node_ids.insert(global_node_ids.end(),
375 local_ids_for_subdomain_ghost_nodes.begin(),
376 local_ids_for_subdomain_ghost_nodes.end());
378 return global_node_ids;
384 unsigned long number_of_local_nodes = subdomain_mesh->
getNodes().size();
385 unsigned long number_of_global_nodes = 0;
387 MPI_Comm mpi_comm = MPI_COMM_WORLD;
389 MPI_Allreduce(&number_of_local_nodes, &number_of_global_nodes, 1,
390 MPI_UNSIGNED_LONG, MPI_SUM, mpi_comm);
391 DBUG(
"[{}] number_of_global_nodes: {}'", subdomain_mesh->
getName(),
392 number_of_global_nodes);
393 return number_of_global_nodes;
398 Mesh const*
const subdomain_mesh)
400 DBUG(
"Creating NodePartitionedMesh from '{}'", subdomain_mesh->
getName());
402 auto const subdomain_global_node_ids =
407 unsigned long number_of_global_base_nodes = 0;
409 unsigned long const number_of_global_nodes =
411 auto numbers_of_regular_base_nodes_at_rank =
413 auto numbers_of_regular_higher_order_nodes_at_rank =
416 auto const number_of_regular_nodes =
418 DBUG(
"[{}] number_of_regular_nodes: {}", subdomain_mesh->
getName(),
419 number_of_regular_nodes);
420 auto const [nodes, elements] =
423 return std::make_unique<NodePartitionedMesh>(
424 subdomain_mesh->
getName(), nodes, subdomain_global_node_ids, elements,
425 subdomain_mesh->
getProperties(), number_of_global_base_nodes,
426 number_of_global_nodes, number_of_regular_nodes,
427 std::move(numbers_of_regular_base_nodes_at_rank),
428 std::move(numbers_of_regular_higher_order_nodes_at_rank));
Definition of Duplicate functions.
Definition of the Element class.
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Definition of the Mesh class.
Definition of mesh class for partitioned mesh (by node) for parallel computing within the framework o...
Definition of the Node class.
std::size_t getID() const
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
std::size_t computeNumberOfBaseNodes() const
Get the number of base nodes.
Properties & getProperties()
const std::string getName() const
Get name of the mesh.
std::size_t getNumberOfNodes() const
Get the number of nodes.
bool isGhostNode(const std::size_t node_id) const
Check whether a node with ID of node_id is a ghost node.
std::size_t getGlobalNodeID(const std::size_t node_id) const
Get the global node ID of a node with its local ID.
constexpr ranges::views::view_closure ids
For an element of a range view return its id.
std::pair< std::vector< Node * >, std::vector< Element * > > copyNodesAndElements(std::vector< Element * > const &input_elements)
unsigned long getNumberOfGlobalNodes(Mesh const *subdomain_mesh)
std::vector< std::size_t > computeNumberOfRegularHigherOrderNodesAtRank(Mesh const *subdomain_mesh)
std::vector< std::size_t > computeGlobalNodeIDsOfSubDomainPartition(NodePartitionedMesh const *bulk_mesh, Mesh const *subdomain_mesh)
std::vector< std::size_t > computeNumberOfRegularBaseNodesAtRank(Mesh const *subdomain_mesh)
std::vector< std::size_t > computeRegularBaseNodeGlobalNodeIDsOfSubDomainPartition(NodePartitionedMesh const *bulk_mesh, Mesh const *subdomain_mesh)
std::pair< int, int > getMPIRankAndSize(MPI_Comm const &mpi_comm)
unsigned long computeNumberOfRegularNodes(NodePartitionedMesh const *bulk_mesh, Mesh const *subdomain_mesh)
std::vector< Element * > cloneElements(std::vector< Element * > const &elements)
Clones a vector of elements using the Element::clone() function.
std::vector< std::size_t > computeGhostBaseNodeGlobalNodeIDsOfSubDomainPartition(NodePartitionedMesh const *bulk_mesh, Mesh const *subdomain_mesh, std::vector< std::size_t > const &global_regular_base_node_ids)
PropertyVector< std::size_t > const * bulkNodeIDs(Mesh const &mesh)
std::unique_ptr< NodePartitionedMesh > transformMeshToNodePartitionedMesh(NodePartitionedMesh const *const bulk_mesh, Mesh const *const subdomain_mesh)