OGS
transformMeshToNodePartitionedMesh.cpp
Go to the documentation of this file.
1
11
12#include <mpi.h>
13
14#include <numeric>
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>
21#include <vector>
22
23#include "BaseLib/Logging.h"
25#include "MeshLib/Mesh.h"
26#include "MeshLib/Node.h"
29
30namespace
31{
33 MeshLib::NodePartitionedMesh const& bulk_mesh,
34 std::vector<std::size_t> const& local_bulk_node_ids_for_subdomain,
35 std::size_t const subdomain_node_id)
36{
37 return (!bulk_mesh.isGhostNode(
38 local_bulk_node_ids_for_subdomain[subdomain_node_id]));
39};
40} // namespace
41
42namespace MeshLib
43{
44std::pair<int, int> getMPIRankAndSize(MPI_Comm const& mpi_comm)
45{
46 int mpi_comm_size;
47 MPI_Comm_size(mpi_comm, &mpi_comm_size);
48 int mpi_comm_rank;
49 MPI_Comm_rank(mpi_comm, &mpi_comm_rank);
50 return {mpi_comm_rank, mpi_comm_size};
51}
52
53std::pair<std::vector<Node*>, std::vector<Element*>> copyNodesAndElements(
54 std::vector<Element*> const& input_elements)
55{
56 auto elements = cloneElements(input_elements);
57
58 // original node ids to newly created nodes.
59 std::unordered_map<std::size_t, Node*> id_node_hash_map;
60 id_node_hash_map.reserve(
61 elements.size()); // There will be at least one node per element.
62
63 for (auto& e : elements)
64 {
65 // For each node find a cloned node in map or create if there is none.
66 unsigned const n_nodes = e->getNumberOfNodes();
67 for (unsigned i = 0; i < n_nodes; ++i)
68 {
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())
72 {
73 auto new_node_in_map = id_node_hash_map[n->getID()] =
74 new Node(*n);
75 e->setNode(i, new_node_in_map);
76 }
77 else
78 {
79 e->setNode(i, it->second);
80 }
81 }
82 }
83
84 std::map<std::size_t, Node*> nodes_map;
85 for (auto const& n : id_node_hash_map)
86 {
87 nodes_map[n.first] = n.second;
88 }
89
90 // Copy the unique nodes pointers.
91 auto element_nodes =
92 nodes_map | ranges::views::values | ranges::to<std::vector>;
93
94 return std::make_pair(element_nodes, elements);
95}
96
97unsigned long computeNumberOfRegularNodes(NodePartitionedMesh const* bulk_mesh,
98 Mesh const* subdomain_mesh)
99{
100 auto const subdomain_nodes = subdomain_mesh->getNodes();
101 auto const local_bulk_node_ids_for_subdomain =
102 *bulkNodeIDs(*subdomain_mesh);
103 auto const subdomain_node_ids =
104 subdomain_nodes | MeshLib::views::ids | ranges::to<std::vector>;
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));
109
110 DBUG("[{}] number of regular nodes: {}", subdomain_mesh->getName(),
111 number_of_regular_nodes);
112 return number_of_regular_nodes;
113}
114
115std::vector<std::size_t>
117 NodePartitionedMesh const* bulk_mesh, Mesh const* subdomain_mesh)
118{
119 // create/fill global nodes information of the sub-domain partition
120 auto const subdomain_nodes = subdomain_mesh->getNodes();
121 auto const local_bulk_node_ids_for_subdomain =
122 *bulkNodeIDs(*subdomain_mesh);
123
124 // global node ids for the sub-domain:
125 // | regular base node ids | ghost base node ids | regular higher
126 // order node ids | ghost higher order node ids |
127 // | partition offset + local ids | regular node ids from other
128 // partitions |
129
130 // In order to compute the global node ids for the subdomain mesh nodes the
131 // sum of the number of regular nodes of the previous partitions is required
132 // first.
133 // Count the local regular base nodes to compute offset for other subsequent
134 // partitions
135 std::size_t const number_of_regular_nodes =
136 computeNumberOfRegularNodes(bulk_mesh, subdomain_mesh);
137
138 // in the following information exchange with other ranks is required
139 MPI_Comm mpi_comm = MPI_COMM_WORLD;
140 auto const [mpi_comm_rank, mpi_comm_size] = getMPIRankAndSize(mpi_comm);
141
142 // send own number of regular nodes to all others
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,
146 mpi_comm);
147 // compute the 'offset' in the global_node_ids
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));
153
154 // add the offset to the partitioned-owned subdomain
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(),
160 partition_offset);
161 // set the global id for the regular base nodes
162 for (auto const id : subdomain_nodes | MeshLib::views::ids)
163 {
164 if (isRegularNode(*bulk_mesh, local_bulk_node_ids_for_subdomain, id))
165 {
166 subdomain_global_node_ids.emplace_back(partition_offset + id);
167 }
168 }
169 return subdomain_global_node_ids;
170}
171
173 NodePartitionedMesh const* bulk_mesh,
174 Mesh const* subdomain_mesh,
175 std::vector<std::size_t> const& global_regular_base_node_ids)
176{
177 // in the following information exchange with other ranks is required
178 MPI_Comm const mpi_comm = MPI_COMM_WORLD;
179 auto const [mpi_comm_rank, mpi_comm_size] = getMPIRankAndSize(mpi_comm);
180
181 // count regular nodes that is the offset in the mapping
182 auto const local_bulk_node_ids_for_subdomain =
183 *bulkNodeIDs(*subdomain_mesh);
184
185 // compute mapping subdomain ghost node id <-> bulk ghost node id
186 std::vector<std::size_t> subdomain_node_id_to_bulk_node_id;
187 for (auto const id : subdomain_mesh->getNodes() | MeshLib::views::ids)
188 {
189 if (!isRegularNode(*bulk_mesh, local_bulk_node_ids_for_subdomain, id))
190 {
191 auto const bulk_node_id = local_bulk_node_ids_for_subdomain[id];
192 // this puts the global bulk node id at pos:
193 // subdomain_node->getID() - number_of_regular_nodes
194 subdomain_node_id_to_bulk_node_id.push_back(
195 bulk_mesh->getGlobalNodeID(bulk_node_id));
196 }
197 }
198
199 // send ids of bulk ghost nodes belonging to the subdomain mesh to all
200 // other ranks and at the same time receive from all other ranks
201 // first send the sizes to all other to are able to allocate buffer
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);
206
207 DBUG("[{}] global_number_of_subdomain_node_id_to_bulk_node_id: '{}' ",
208 subdomain_mesh->getName(),
209 global_number_of_subdomain_node_id_to_bulk_node_id);
210
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,
213 mpi_comm);
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(), /* sendbuf */
221 size, /* sendcount */
222 MPI_UNSIGNED_LONG, /* sendtype */
223 ghost_node_ids_of_all_ranks.data(), /* recvbuf (out) */
224 numbers_of_ids_at_ranks.data(), /* recvcounts */
225 offsets.data(), /* displs */
226 MPI_UNSIGNED_LONG, /* recvtype */
227 mpi_comm);
228
229 // construct a map for fast search of local bulk node ids
230 std::map<std::size_t, std::size_t> global_to_local_bulk_node_ids;
231 for (auto const id : bulk_mesh->getNodes() | MeshLib::views::ids)
232 {
233 if (!bulk_mesh->isGhostNode(id))
234 {
235 global_to_local_bulk_node_ids[bulk_mesh->getGlobalNodeID(id)] = id;
236 }
237 }
238
239 // construct a map for fast search of local sub-domain node ids
240 std::map<std::size_t, std::size_t> local_subdomain_node_ids;
241 for (auto const id : subdomain_mesh->getNodes() | MeshLib::views::ids)
242 {
243 local_subdomain_node_ids[local_bulk_node_ids_for_subdomain[id]] = id;
244 }
245
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());
249 // search in all ranks within the bulk ids for the corresponding id
250 for (int rank = 0; rank < mpi_comm_size; ++rank)
251 {
252 if (rank == mpi_comm_rank)
253 {
254 continue;
255 }
256 for (int i = offsets[rank]; i < offsets[rank + 1]; ++i)
257 {
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())
261 {
262 continue;
263 }
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];
268 DBUG(
269 "[{}] found global subdomain node id: '{}' for global bulk "
270 "node id {} ",
271 subdomain_mesh->getName(),
272 local_subdomain_node_ids_of_all_ranks[i],
273 ghost_node_ids_of_all_ranks[i]);
274 }
275 }
276
277 // send the computed ids back
278 std::vector<std::size_t> computed_global_ids_for_subdomain_ghost_nodes(
279 global_number_of_subdomain_node_id_to_bulk_node_id);
280 MPI_Allreduce(
281 local_subdomain_node_ids_of_all_ranks.data(), /* sendbuf */
282 computed_global_ids_for_subdomain_ghost_nodes
283 .data(), /* recvbuf (out) */
284 global_number_of_subdomain_node_id_to_bulk_node_id, /* sendcount */
285 MPI_UNSIGNED_LONG, /* sendtype */
286 MPI_MAX, /* operation */
287 mpi_comm);
288
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;
295}
296
298 Mesh const* subdomain_mesh)
299{
300 auto const number_of_regular_base_nodes =
301 subdomain_mesh->computeNumberOfBaseNodes();
302
303 MPI_Comm const mpi_comm = MPI_COMM_WORLD;
304 auto const [mpi_comm_rank, mpi_comm_size] = getMPIRankAndSize(mpi_comm);
305
306 std::vector<std::size_t> gathered_number_of_regular_base_nodes(
307 mpi_comm_size);
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);
311
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));
317
318 return numbers_of_regular_base_nodes_at_rank;
319}
320
321// similar to the above only with regular higher order nodes
323 Mesh const* subdomain_mesh)
324{
325 // in the following information exchange with other ranks is required
326 MPI_Comm const mpi_comm = MPI_COMM_WORLD;
327 auto [mpi_comm_rank, mpi_comm_size] = getMPIRankAndSize(mpi_comm);
328
329 auto const number_of_regular_base_nodes =
330 subdomain_mesh->computeNumberOfBaseNodes();
331
332 std::vector<std::size_t> gathered_number_of_regular_higher_order_nodes(
333 mpi_comm_size);
334 auto const number_of_regular_higher_order_nodes =
335 subdomain_mesh->getNumberOfNodes() - number_of_regular_base_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);
339
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);
342 std::partial_sum(
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));
346
347 return numbers_of_regular_higher_order_nodes_at_rank;
348}
349
351 NodePartitionedMesh const* bulk_mesh, Mesh const* subdomain_mesh)
352{
353 auto global_regular_base_node_ids =
355 subdomain_mesh);
356 auto const local_ids_for_subdomain_ghost_nodes =
358 bulk_mesh, subdomain_mesh, global_regular_base_node_ids);
359
360 // At the moment higher order bulk meshes aren't handled correctly.
361 // If it should become necessary in the future, the necessary things must be
362 // implemented at this point.
363 /*
364 auto const regular_higher_order_node_global_node_ids =
365 computeRegularHigherOrderNodeGlobalNodeIDsofSubDomainPartition(
366 bulk_mesh, subdomain_mesh);
367 auto const ghost_higher_order_node_global_node_ids =
368 computeGhostHigherOrderNodeGlobalNodeIDsofSubDomainPartition(
369 bulk_mesh, subdomain_mesh);
370 */
371
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());
377
378 return global_node_ids;
379}
380
381unsigned long getNumberOfGlobalNodes(Mesh const* subdomain_mesh)
382{
383 // sum all nodes over all partitions in number_of_global_nodes
384 unsigned long number_of_local_nodes = subdomain_mesh->getNodes().size();
385 unsigned long number_of_global_nodes = 0;
386
387 MPI_Comm mpi_comm = MPI_COMM_WORLD;
388
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;
394}
395
396std::unique_ptr<NodePartitionedMesh> transformMeshToNodePartitionedMesh(
397 NodePartitionedMesh const* const bulk_mesh,
398 Mesh const* const subdomain_mesh)
399{
400 DBUG("Creating NodePartitionedMesh from '{}'", subdomain_mesh->getName());
401
402 auto const subdomain_global_node_ids =
403 computeGlobalNodeIDsOfSubDomainPartition(bulk_mesh, subdomain_mesh);
404
405 // according to comment in struct PartitionedMeshInfo this value
406 // is unused
407 unsigned long number_of_global_base_nodes = 0;
408
409 unsigned long const number_of_global_nodes =
410 getNumberOfGlobalNodes(subdomain_mesh);
411 auto numbers_of_regular_base_nodes_at_rank =
413 auto numbers_of_regular_higher_order_nodes_at_rank =
415
416 auto const number_of_regular_nodes =
417 computeNumberOfRegularNodes(bulk_mesh, subdomain_mesh);
418 DBUG("[{}] number_of_regular_nodes: {}", subdomain_mesh->getName(),
419 number_of_regular_nodes);
420 auto const [nodes, elements] =
421 copyNodesAndElements(subdomain_mesh->getElements());
422
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));
429}
430
431} // namespace MeshLib
Definition of Duplicate functions.
Definition of the Element class.
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Definition: Logging.h:30
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
Definition: Point3dWithID.h:62
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
Definition: Mesh.h:102
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition: Mesh.h:105
std::size_t computeNumberOfBaseNodes() const
Get the number of base nodes.
Definition: Mesh.cpp:228
Properties & getProperties()
Definition: Mesh.h:130
const std::string getName() const
Get name of the mesh.
Definition: Mesh.h:99
std::size_t getNumberOfNodes() const
Get the number of nodes.
Definition: Mesh.h:96
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.
Definition: Mesh.h:217
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)
Definition: Mesh.cpp:282
std::unique_ptr< NodePartitionedMesh > transformMeshToNodePartitionedMesh(NodePartitionedMesh const *const bulk_mesh, Mesh const *const subdomain_mesh)
bool isRegularNode(MeshLib::NodePartitionedMesh const &bulk_mesh, std::vector< std::size_t > const &local_bulk_node_ids_for_subdomain, std::size_t const subdomain_node_id)