OGS 6.2.1-97-g73d1aeda3
NodeWiseMeshPartitioner.cpp
Go to the documentation of this file.
1 
16 
17 #include <limits>
18 #include <numeric>
19 #include <unordered_map>
20 
21 #include <logog/include/logog.hpp>
22 
23 #include "BaseLib/Error.h"
24 #include "BaseLib/Stream.h"
25 
27 
28 namespace ApplicationUtils
29 {
30 struct NodeStruct
31 {
33  double const x_,
34  double const y_,
35  double const z_)
36  : id(id_), x(x_), y(y_), z(z_)
37  {
38  }
39 
41  double x;
42  double y;
43  double z;
44 };
45 
47  MeshLib::MeshItemType const item_type) const
48 {
49  if (item_type == MeshLib::MeshItemType::Node)
50  {
51  return nodes.size();
52  }
53 
54  if (item_type == MeshLib::MeshItemType::Cell)
55  {
56  return regular_elements.size() + ghost_elements.size();
57  }
58  OGS_FATAL("Mesh items other than nodes and cells are not supported.");
59 }
60 
62  std::ostream& os, std::vector<std::size_t> const& global_node_ids) const
63 {
64  std::vector<NodeStruct> nodes_buffer;
65  nodes_buffer.reserve(nodes.size());
66 
67  for (const auto* node : nodes)
68  {
69  double const* coords = node->getCoords();
70  nodes_buffer.emplace_back(global_node_ids[node->getID()], coords[0],
71  coords[1], coords[2]);
72  }
73  return os.write(reinterpret_cast<const char*>(nodes_buffer.data()),
74  sizeof(NodeStruct) * nodes_buffer.size());
75 }
76 
82  std::vector<const MeshLib::Element*> const& elements)
83 {
84  return 3 * elements.size() +
85  std::accumulate(begin(elements), end(elements), 0,
86  [](auto const nnodes, auto const* e) {
87  return nnodes + e->getNumberOfNodes();
88  });
89 }
90 
91 std::ostream& Partition::writeConfigBinary(std::ostream& os) const
92 {
93  long const data[] = {
94  static_cast<long>(nodes.size()),
95  static_cast<long>(number_of_base_nodes),
96  static_cast<long>(regular_elements.size()),
97  static_cast<long>(ghost_elements.size()),
98  static_cast<long>(number_of_non_ghost_base_nodes),
99  static_cast<long>(number_of_non_ghost_nodes),
100  static_cast<long>(number_of_mesh_base_nodes),
101  static_cast<long>(number_of_mesh_all_nodes),
102  static_cast<long>(
103  getNumberOfIntegerVariablesOfElements(regular_elements)),
104  static_cast<long>(
105  getNumberOfIntegerVariablesOfElements(ghost_elements)),
106  };
107 
108  return os.write(reinterpret_cast<const char*>(data), sizeof(data));
109 }
110 
111 std::size_t nodeIdBulkMesh(
112  MeshLib::Node const& node,
113  std::vector<std::size_t> const* node_id_mapping = nullptr)
114 {
115  return node_id_mapping ? (*node_id_mapping)[node.getID()] : node.getID();
116 }
117 
118 std::size_t partitionLookup(
119  MeshLib::Node const& node,
120  std::vector<std::size_t> const& partition_ids,
121  std::vector<std::size_t> const* node_id_mapping = nullptr)
122 {
123  auto node_id = [&node_id_mapping](MeshLib::Node const& n) {
124  return nodeIdBulkMesh(n, node_id_mapping);
125  };
126 
127  return partition_ids[node_id(node)];
128 }
129 
137 std::pair<std::vector<MeshLib::Node*>, std::vector<MeshLib::Node*>>
139  std::size_t const part_id,
140  const bool is_mixed_high_order_linear_elems,
141  std::size_t const n_base_nodes,
142  std::vector<MeshLib::Node*> const& nodes,
143  std::vector<std::size_t> const& partition_ids,
144  std::vector<std::size_t> const* node_id_mapping = nullptr)
145 {
146  auto node_id = [&node_id_mapping](MeshLib::Node const& n) {
147  return nodeIdBulkMesh(n, node_id_mapping);
148  };
149 
150  // Find nodes belonging to a given partition id.
151  std::vector<MeshLib::Node*> partition_nodes;
152  copy_if(begin(nodes), end(nodes), std::back_inserter(partition_nodes),
153  [&](auto const& n) {
154  return partitionLookup(*n, partition_ids, node_id_mapping) ==
155  part_id;
156  });
157 
158  // Space for resulting vectors.
159  std::vector<MeshLib::Node*> base_nodes;
160  base_nodes.reserve(partition_nodes.size() /
161  2); // if linear mesh, then one reallocation, no realloc
162  // for higher order elements meshes.
163  std::vector<MeshLib::Node*> higher_order_nodes;
164  higher_order_nodes.reserve(
165  partition_nodes.size() /
166  2); // if linear mesh, then wasted space, good estimate for quadratic
167  // order mesh, and realloc needed for higher order element meshes.
168 
169  // Split the nodes into base nodes and extra nodes.
170  partition_copy(begin(partition_nodes), end(partition_nodes),
171  std::back_inserter(base_nodes),
172  std::back_inserter(higher_order_nodes),
173  [&](MeshLib::Node* const n) {
174  return !is_mixed_high_order_linear_elems ||
175  node_id(*n) > n_base_nodes;
176  });
177 
178  return {base_nodes, higher_order_nodes};
179 }
180 
181 std::ptrdiff_t numberOfRegularNodes(
182  MeshLib::Element const& e, std::size_t const part_id,
183  std::vector<std::size_t> const& partition_ids,
184  std::vector<std::size_t> const* node_id_mapping = nullptr)
185 {
186  return std::count_if(e.getNodes(), e.getNodes() + e.getNumberOfNodes(),
187  [&](MeshLib::Node* const n) {
188  return partitionLookup(*n, partition_ids,
189  node_id_mapping) == part_id;
190  });
191 }
192 
197 std::tuple<std::vector<MeshLib::Element const*>,
198  std::vector<MeshLib::Element const*>>
200  std::size_t const part_id,
201  std::vector<MeshLib::Element*> const& elements,
202  std::vector<std::size_t> const& partition_ids,
203  std::vector<std::size_t> const* node_id_mapping = nullptr)
204 {
205  std::vector<MeshLib::Element const*> regular_elements;
206  std::vector<MeshLib::Element const*> ghost_elements;
207 
208  for (auto elem : elements)
209  {
210  auto const regular_nodes = numberOfRegularNodes(
211  *elem, part_id, partition_ids, node_id_mapping);
212 
213  if (regular_nodes == 0)
214  {
215  continue;
216  }
217 
218  if (regular_nodes ==
219  static_cast<std::ptrdiff_t>(elem->getNumberOfNodes()))
220  {
221  regular_elements.push_back(elem);
222  }
223  else
224  {
225  ghost_elements.push_back(elem);
226  }
227  }
228  return std::tuple<std::vector<MeshLib::Element const*>,
229  std::vector<MeshLib::Element const*>>{regular_elements,
230  ghost_elements};
231 }
232 
237 std::tuple<std::vector<MeshLib::Node*>, std::vector<MeshLib::Node*>>
239  std::size_t const part_id,
240  const bool is_mixed_high_order_linear_elems,
241  std::size_t const number_of_base_nodes,
242  std::vector<MeshLib::Node*> const& nodes,
243  std::vector<MeshLib::Element const*> const& ghost_elements,
244  std::vector<std::size_t> const& partition_ids,
245  std::vector<std::size_t> const* node_id_mapping = nullptr)
246 {
247  auto node_id = [&node_id_mapping](MeshLib::Node const& n) {
248  return nodeIdBulkMesh(n, node_id_mapping);
249  };
250 
251  std::vector<MeshLib::Node*> base_nodes;
252  std::vector<MeshLib::Node*> ghost_nodes;
253 
254  std::vector<bool> nodes_reserved(nodes.size(), false);
255  for (const auto* ghost_elem : ghost_elements)
256  {
257  for (unsigned i = 0; i < ghost_elem->getNumberOfNodes(); i++)
258  {
259  auto const& n = ghost_elem->getNode(i);
260  if (nodes_reserved[n->getID()])
261  {
262  continue;
263  }
264 
265  if (partitionLookup(*n, partition_ids, node_id_mapping) != part_id)
266  {
267  if (!is_mixed_high_order_linear_elems ||
268  node_id(*n) > number_of_base_nodes)
269  {
270  base_nodes.push_back(nodes[n->getID()]);
271  }
272  else
273  {
274  ghost_nodes.push_back(nodes[n->getID()]);
275  }
276  nodes_reserved[n->getID()] = true;
277  }
278  }
279  }
280  return std::tuple<std::vector<MeshLib::Node*>, std::vector<MeshLib::Node*>>{
281  base_nodes, ghost_nodes};
282 }
283 
285  std::size_t const part_id, const bool is_mixed_high_order_linear_elems)
286 {
287  auto& partition = _partitions[part_id];
288  std::vector<MeshLib::Node*> higher_order_regular_nodes;
289  std::tie(partition.nodes, higher_order_regular_nodes) =
290  findNonGhostNodesInPartition(part_id, is_mixed_high_order_linear_elems,
291  _mesh->getNumberOfBaseNodes(),
292  _mesh->getNodes(), _nodes_partition_ids);
293 
294  partition.number_of_non_ghost_base_nodes = partition.nodes.size();
295  partition.number_of_non_ghost_nodes =
296  partition.number_of_non_ghost_base_nodes +
297  higher_order_regular_nodes.size();
298 
299  std::tie(partition.regular_elements, partition.ghost_elements) =
300  findElementsInPartition(part_id, _mesh->getElements(),
301  _nodes_partition_ids);
302  std::vector<MeshLib::Node*> base_ghost_nodes;
303  std::vector<MeshLib::Node*> higher_order_ghost_nodes;
304  std::tie(base_ghost_nodes, higher_order_ghost_nodes) =
305  findGhostNodesInPartition(part_id, is_mixed_high_order_linear_elems,
306  _mesh->getNumberOfBaseNodes(),
307  _mesh->getNodes(), partition.ghost_elements,
308  _nodes_partition_ids);
309 
310  std::copy(begin(base_ghost_nodes), end(base_ghost_nodes),
311  std::back_inserter(partition.nodes));
312 
313  partition.number_of_base_nodes = partition.nodes.size();
314 
315  if (is_mixed_high_order_linear_elems)
316  {
317  std::copy(begin(higher_order_regular_nodes),
318  end(higher_order_regular_nodes),
319  std::back_inserter(partition.nodes));
320  std::copy(begin(higher_order_ghost_nodes),
321  end(higher_order_ghost_nodes),
322  std::back_inserter(partition.nodes));
323  }
324 
325  // Set the node numbers of base and all mesh nodes.
326  partition.number_of_mesh_base_nodes = _mesh->getNumberOfBaseNodes();
327  partition.number_of_mesh_all_nodes = _mesh->getNumberOfNodes();
328 }
329 
332 template <typename T>
334  Partition const& p,
335  std::size_t const offset,
336  MeshLib::PropertyVector<T> const& pv,
337  MeshLib::PropertyVector<T>& partitioned_pv)
338 {
339  auto const& nodes = p.nodes;
340  auto const nnodes = nodes.size();
341  for (std::size_t i = 0; i < nnodes; ++i)
342  {
343  const auto global_id = nodes[i]->getID();
344  partitioned_pv[offset + i] = pv[global_id];
345  }
346  return nnodes;
347 }
348 
352 template <typename T>
354  Partition const& p,
355  std::size_t const offset,
356  MeshLib::PropertyVector<T> const& pv,
357  MeshLib::PropertyVector<T>& partitioned_pv)
358 {
359  std::size_t const n_regular(p.regular_elements.size());
360  for (std::size_t i = 0; i < n_regular; ++i)
361  {
362  const auto id = p.regular_elements[i]->getID();
363  partitioned_pv[offset + i] = pv[id];
364  }
365 
366  std::size_t const n_ghost(p.ghost_elements.size());
367  for (std::size_t i = 0; i < n_ghost; ++i)
368  {
369  const auto id = p.ghost_elements[i]->getID();
370  partitioned_pv[offset + n_regular + i] = pv[id];
371  }
372  return n_regular + n_ghost;
373 }
374 
375 template <typename T>
376 bool copyPropertyVector(MeshLib::Properties const& original_properties,
377  MeshLib::Properties& partitioned_properties,
378  std::vector<Partition> const& partitions,
379  std::string const& name,
380  std::size_t const total_number_of_tuples)
381 {
382  if (!original_properties.existsPropertyVector<T>(name))
383  {
384  return false;
385  }
386 
387  auto const& pv = original_properties.getPropertyVector<T>(name);
388  auto partitioned_pv = partitioned_properties.createNewPropertyVector<T>(
389  name, pv->getMeshItemType(), pv->getNumberOfComponents());
390  partitioned_pv->resize(total_number_of_tuples *
391  pv->getNumberOfComponents());
392 
393  auto copy_property_vector_values = [&](Partition const& p,
394  std::size_t offset) {
395  if (pv->getMeshItemType() == MeshLib::MeshItemType::Node)
396  {
397  return copyNodePropertyVectorValues(p, offset, *pv,
398  *partitioned_pv);
399  }
400  if (pv->getMeshItemType() == MeshLib::MeshItemType::Cell)
401  {
402  return copyCellPropertyVectorValues(p, offset, *pv,
403  *partitioned_pv);
404  }
405  OGS_FATAL(
406  "Copying of property vector values for mesh item type %s is "
407  "not implemented.",
408  pv->getMeshItemType());
409  };
410 
411  std::size_t position_offset(0);
412  for (auto p : partitions)
413  {
414  position_offset += copy_property_vector_values(p, position_offset);
415  }
416  return true;
417 }
418 
424 template <typename Function>
425 void applyToPropertyVectors(std::vector<std::string> const& property_names,
426  Function f)
427 {
428  for (auto const& name : property_names)
429  {
430  // Open question, why is the 'unsigned long' case not compiling giving
431  // an error "expected '(' for function-style cast or type construction"
432  // with clang-7, and "error C4576: a parenthesized type followed by an
433  // initializer list is a non-standard explicit type conversion syntax"
434  // with MSVC-15.
435  bool success =
436  f(double{}, name) || f(float{}, name) || f(int{}, name) ||
437  f(long{}, name) || f(unsigned{}, name) ||
438  f(static_cast<unsigned long>(0), name) || f(std::size_t{}, name);
439  if (!success)
440  {
441  OGS_FATAL("Could not apply function to PropertyVector '%s'.",
442  name.c_str());
443  }
444  }
445 }
446 
447 void processProperties(MeshLib::Properties const& properties,
448  MeshLib::MeshItemType const mesh_item_type,
449  std::vector<Partition> const& partitions,
450  MeshLib::Properties& partitioned_properties)
451 {
452  std::size_t const total_number_of_tuples =
453  std::accumulate(begin(partitions), end(partitions), 0,
454  [&](std::size_t const sum, Partition const& p) {
455  return sum + p.numberOfMeshItems(mesh_item_type);
456  });
457 
458  DBUG(
459  "total number of tuples define on mesh item type '%d' after "
460  "partitioning: %d ",
461  mesh_item_type, total_number_of_tuples);
462 
463  // 1 create new PV
464  // 2 resize the PV with total_number_of_tuples
465  // 3 copy the values according to the partition info
466  applyToPropertyVectors(properties.getPropertyVectorNames(mesh_item_type),
467  [&](auto type, std::string const& name) {
468  return copyPropertyVector<decltype(type)>(
469  properties, partitioned_properties,
470  partitions, name, total_number_of_tuples);
471  });
472 }
473 
475  MeshLib::Properties const& properties,
476  std::vector<Partition> const& partitions)
477 {
478  MeshLib::Properties partitioned_properties;
479  processProperties(properties, MeshLib::MeshItemType::Node, partitions,
480  partitioned_properties);
481  processProperties(properties, MeshLib::MeshItemType::Cell, partitions,
482  partitioned_properties);
483  return partitioned_properties;
484 }
485 
487  const bool is_mixed_high_order_linear_elems)
488 {
489  for (std::size_t part_id = 0; part_id < _partitions.size(); part_id++)
490  {
491  INFO("Processing partition: %d", part_id);
492  processPartition(part_id, is_mixed_high_order_linear_elems);
493  }
494 
495  renumberNodeIndices(is_mixed_high_order_linear_elems);
496 
497  _partitioned_properties =
498  partitionProperties(_mesh->getProperties(), _partitions);
499 }
500 
502  MeshLib::PropertyVector<std::size_t>* const bulk_node_ids_pv,
503  std::vector<Partition> const& local_partitions) const
504 {
505  if (bulk_node_ids_pv == nullptr)
506  {
507  return;
508  }
509 
510  auto& bulk_node_ids = *bulk_node_ids_pv;
511 
512  std::size_t offset = 0; // offset in property vector for current partition
513 
514  assert(_partitions.size() == local_partitions.size());
515  int const n_partitions = static_cast<int>(_partitions.size());
516  for (int partition_id = 0; partition_id < n_partitions; ++partition_id)
517  {
518  auto const& bulk_partition = _partitions[partition_id];
519  auto const& local_partition = local_partitions[partition_id];
520 
521  // Create global-to-local node id mapping for the bulk partition.
522  auto const& bulk_nodes = bulk_partition.nodes;
523  auto const n_bulk_nodes = bulk_nodes.size();
524  std::map<std::size_t, std::size_t> global_to_local;
525  for (std::size_t local_node_id = 0; local_node_id < n_bulk_nodes;
526  ++local_node_id)
527  {
528  global_to_local[bulk_nodes[local_node_id]->getID()] = local_node_id;
529  }
530 
531  auto const& local_nodes = local_partition.nodes;
532  auto const n_local_nodes = local_nodes.size();
533  for (std::size_t local_node_id = 0; local_node_id < n_local_nodes;
534  ++local_node_id)
535  {
536  bulk_node_ids[offset + local_node_id] =
537  global_to_local[bulk_node_ids[offset + local_node_id]];
538  }
539  offset += n_local_nodes;
540  }
541 }
542 
544  MeshLib::PropertyVector<std::size_t>* const bulk_element_ids_pv,
545  std::vector<Partition> const& local_partitions) const
546 {
547  if (bulk_element_ids_pv == nullptr)
548  {
549  return;
550  }
551 
552  auto& bulk_element_ids = *bulk_element_ids_pv;
553 
554  std::size_t offset = 0; // offset in property vector for current partition
555 
556  assert(_partitions.size() == local_partitions.size());
557  int const n_partitions = static_cast<int>(_partitions.size());
558  for (int partition_id = 0; partition_id < n_partitions; ++partition_id)
559  {
560  auto const& bulk_partition = _partitions[partition_id];
561  auto const& local_partition = local_partitions[partition_id];
562 
563  // Create global-to-local element id mapping for the bulk partition.
564  std::map<std::size_t, std::size_t> global_to_local;
565  auto map_elements =
566  [&global_to_local](
567  std::vector<MeshLib::Element const*> const& elements,
568  std::size_t const offset) {
569  auto const n_elements = elements.size();
570  for (std::size_t e = 0; e < n_elements; ++e)
571  {
572  global_to_local[elements[e]->getID()] = offset + e;
573  }
574  };
575 
576  map_elements(bulk_partition.regular_elements, 0);
577  map_elements(bulk_partition.ghost_elements,
578  bulk_partition.regular_elements.size());
579 
580  // Renumber the local bulk_element_ids map.
581  auto renumber_elements =
582  [&bulk_element_ids, &global_to_local](
583  std::vector<MeshLib::Element const*> const& elements,
584  std::size_t const offset) {
585  auto const n_elements = elements.size();
586  for (std::size_t e = 0; e < n_elements; ++e)
587  {
588  bulk_element_ids[offset + e] =
589  global_to_local[bulk_element_ids[offset + e]];
590  }
591  return n_elements;
592  };
593 
594  offset += renumber_elements(local_partition.regular_elements, offset);
595  offset += renumber_elements(local_partition.ghost_elements, offset);
596  }
597 }
598 
600  MeshLib::Mesh const& mesh,
601  bool const is_mixed_high_order_linear_elems) const
602 {
603  auto const& bulk_node_ids =
604  mesh.getProperties().getPropertyVector<std::size_t>(
605  "bulk_node_ids", MeshLib::MeshItemType::Node, 1);
606 
607  std::vector<Partition> partitions(_partitions.size());
608  for (std::size_t part_id = 0; part_id < _partitions.size(); part_id++)
609  {
610  auto& partition = partitions[part_id];
611  INFO("Processing partition: %d", part_id);
612  // Set the node numbers of base and all mesh nodes.
613  partition.number_of_mesh_base_nodes = mesh.getNumberOfBaseNodes();
614  partition.number_of_mesh_all_nodes = mesh.getNumberOfNodes();
615 
616  std::vector<MeshLib::Node*> higher_order_regular_nodes;
617  std::tie(partition.nodes, higher_order_regular_nodes) =
619  part_id, is_mixed_high_order_linear_elems,
620  mesh.getNumberOfBaseNodes(), mesh.getNodes(),
621  _nodes_partition_ids, bulk_node_ids);
622 
623  partition.number_of_non_ghost_base_nodes = partition.nodes.size();
624  partition.number_of_non_ghost_nodes =
625  partition.number_of_non_ghost_base_nodes +
626  higher_order_regular_nodes.size();
627 
628  std::tie(partition.regular_elements, partition.ghost_elements) =
629  findElementsInPartition(part_id, mesh.getElements(),
630  _nodes_partition_ids, bulk_node_ids);
631 
632  std::vector<MeshLib::Node*> base_ghost_nodes;
633  std::vector<MeshLib::Node*> higher_order_ghost_nodes;
634  std::tie(base_ghost_nodes, higher_order_ghost_nodes) =
635  findGhostNodesInPartition(part_id, is_mixed_high_order_linear_elems,
636  mesh.getNumberOfBaseNodes(),
637  mesh.getNodes(), partition.ghost_elements,
638  _nodes_partition_ids, bulk_node_ids);
639 
640  std::copy(begin(base_ghost_nodes), end(base_ghost_nodes),
641  std::back_inserter(partition.nodes));
642 
643  partition.number_of_base_nodes = partition.nodes.size();
644 
645  if (is_mixed_high_order_linear_elems)
646  {
647  std::copy(begin(higher_order_regular_nodes),
648  end(higher_order_regular_nodes),
649  std::back_inserter(partition.nodes));
650  std::copy(begin(higher_order_ghost_nodes),
651  end(higher_order_ghost_nodes),
652  std::back_inserter(partition.nodes));
653  }
654  }
655 
656  return partitions;
657 }
658 
660  const bool is_mixed_high_order_linear_elems)
661 {
662  std::size_t node_global_id_offset = 0;
663  // Renumber the global indices.
664  // -- Base nodes
665  for (auto& partition : _partitions)
666  {
667  for (std::size_t i = 0; i < partition.number_of_non_ghost_base_nodes;
668  i++)
669  {
670  _nodes_global_ids[partition.nodes[i]->getID()] =
671  node_global_id_offset;
672  node_global_id_offset++;
673  }
674  }
675 
676  if (!is_mixed_high_order_linear_elems)
677  {
678  return;
679  }
680 
681  // -- Nodes for high order elements.
682  for (auto& partition : _partitions)
683  {
684  const std::size_t end_id = partition.number_of_base_nodes +
685  partition.number_of_non_ghost_nodes -
686  partition.number_of_non_ghost_base_nodes;
687  for (std::size_t i = partition.number_of_base_nodes; i < end_id; i++)
688  {
689  _nodes_global_ids[partition.nodes[i]->getID()] =
690  node_global_id_offset;
691  node_global_id_offset++;
692  }
693  }
694 }
695 
696 template <typename T>
697 void writePropertyVectorValuesBinary(std::ostream& os,
698  MeshLib::PropertyVector<T> const& pv)
699 {
700  os.write(reinterpret_cast<const char*>(pv.data()), pv.size() * sizeof(T));
701 }
702 
703 template <typename T>
705  MeshLib::Properties const& partitioned_properties, std::string const& name,
706  std::ostream& out_val, std::ostream& out_meta)
707 {
708  if (!partitioned_properties.existsPropertyVector<T>(name))
709  {
710  return false;
711  }
712 
714  pvmd.property_name = name;
715  auto* pv = partitioned_properties.getPropertyVector<T>(name);
717  pvmd.number_of_components = pv->getNumberOfComponents();
718  pvmd.number_of_tuples = pv->getNumberOfTuples();
719  writePropertyVectorValuesBinary(out_val, *pv);
721  return true;
722 }
723 
724 void writePropertiesBinary(const std::string& file_name_base,
725  MeshLib::Properties const& partitioned_properties,
726  std::vector<Partition> const& partitions,
727  MeshLib::MeshItemType const mesh_item_type)
728 {
729  auto const& property_names =
730  partitioned_properties.getPropertyVectorNames(mesh_item_type);
731  if (property_names.empty())
732  {
733  return;
734  }
735 
736  auto const file_name_infix = toString(mesh_item_type);
737 
738  auto const file_name_cfg = file_name_base + "_partitioned_" +
739  file_name_infix + "_properties_cfg" +
740  std::to_string(partitions.size()) + ".bin";
741  std::ofstream out(file_name_cfg, std::ios::binary);
742  if (!out)
743  {
744  OGS_FATAL("Could not open file '%s' for output.",
745  file_name_cfg.c_str());
746  }
747 
748  auto const file_name_val = file_name_base + "_partitioned_" +
749  file_name_infix + "_properties_val" +
750  std::to_string(partitions.size()) + ".bin";
751  std::ofstream out_val(file_name_val, std::ios::binary);
752  if (!out_val)
753  {
754  OGS_FATAL("Could not open file '%s' for output.",
755  file_name_val.c_str());
756  }
757 
758  std::size_t const number_of_properties(property_names.size());
760 
761  applyToPropertyVectors(property_names,
762  [&](auto type, std::string const& name) {
763  return writePropertyVectorBinary<decltype(type)>(
764  partitioned_properties, name, out_val, out);
765  });
766 
767  unsigned long offset = 0;
768  for (const auto& partition : partitions)
769  {
771  offset, static_cast<unsigned long>(
772  partition.numberOfMeshItems(mesh_item_type))};
773  DBUG(
774  "Write meta data for node-based PropertyVector: global offset %d, "
775  "number of tuples %d",
776  pvpmd.offset, pvpmd.number_of_tuples);
778  offset += pvpmd.number_of_tuples;
779  }
780 }
781 
783 {
787 
788  std::ostream& writeConfigBinary(std::ostream& os) const;
789 };
790 
791 std::ostream& ConfigOffsets::writeConfigBinary(std::ostream& os) const
792 {
793  os.write(reinterpret_cast<const char*>(this), sizeof(ConfigOffsets));
794 
795  static long reserved = 0; // Value reserved in the binary format, not used
796  // in the partitioning process.
797  return os.write(reinterpret_cast<const char*>(&reserved), sizeof(long));
798 }
799 
801 {
802  long node;
805 };
806 
808 {
809  return {static_cast<long>(partition.nodes.size()),
810  static_cast<long>(partition.regular_elements.size() +
812  partition.regular_elements)),
813  static_cast<long>(partition.ghost_elements.size() +
815  partition.ghost_elements))};
816 }
817 
819  PartitionOffsets const& offsets)
820 {
821  return {
822  static_cast<long>(oldConfig.node_rank_offset +
823  offsets.node * sizeof(NodeStruct)),
824  // Offset the ending entry of the element integer variales of
825  // the non-ghost elements of this partition in the vector of elem_info.
826  static_cast<long>(oldConfig.element_rank_offset +
827  offsets.regular_elements * sizeof(long)),
828 
829  // Offset the ending entry of the element integer variales of
830  // the ghost elements of this partition in the vector of elem_info.
831  static_cast<long>(oldConfig.ghost_element_rank_offset +
832  offsets.ghost_elements * sizeof(long))};
833 }
834 
840 std::tuple<std::vector<long>, std::vector<long>> writeConfigDataBinary(
841  const std::string& file_name_base, std::vector<Partition> const& partitions)
842 {
843  auto const file_name_cfg = file_name_base + "_partitioned_msh_cfg" +
844  std::to_string(partitions.size()) + ".bin";
845  std::ofstream of_bin_cfg(file_name_cfg, std::ios::binary);
846  if (!of_bin_cfg)
847  {
848  OGS_FATAL("Could not open file '%s' for output.",
849  file_name_cfg.c_str());
850  }
851 
852  std::vector<long> num_elem_integers;
853  num_elem_integers.reserve(partitions.size());
854  std::vector<long> num_g_elem_integers;
855  num_g_elem_integers.reserve(partitions.size());
856 
857  ConfigOffsets config_offsets = {0, 0, 0}; // 0 for first partition.
858  for (const auto& partition : partitions)
859  {
860  partition.writeConfigBinary(of_bin_cfg);
861 
862  config_offsets.writeConfigBinary(of_bin_cfg);
863  auto const& new_offsets = computePartitionElementOffsets(partition);
864  config_offsets = incrementConfigOffsets(config_offsets, new_offsets);
865 
866  num_elem_integers.push_back(new_offsets.regular_elements);
867  num_g_elem_integers.push_back(new_offsets.ghost_elements);
868  }
869 
870  return std::make_tuple(num_elem_integers, num_g_elem_integers);
871 }
872 
881  const MeshLib::Element& elem,
882  const std::unordered_map<std::size_t, long>& local_node_ids,
883  std::vector<long>& elem_info,
884  long& counter)
885 {
886  unsigned mat_id = 0; // TODO: Material ID to be set from the mesh data
887  const long nn = elem.getNumberOfNodes();
888  elem_info[counter++] = mat_id;
889  elem_info[counter++] = static_cast<long>(elem.getCellType());
890  elem_info[counter++] = nn;
891 
892  for (long i = 0; i < nn; i++)
893  {
894  auto const& n = *elem.getNode(i);
895  elem_info[counter++] = local_node_ids.at(n.getID());
896  }
897 }
898 
900 std::unordered_map<std::size_t, long> enumerateLocalNodeIds(
901  std::vector<MeshLib::Node*> const& nodes)
902 {
903  std::unordered_map<std::size_t, long> local_ids;
904  local_ids.reserve(nodes.size());
905 
906  long local_node_id = 0;
907  for (const auto* node : nodes)
908  {
909  local_ids[node->getID()] = local_node_id++;
910  }
911  return local_ids;
912 }
913 
920 void writeElementsBinary(std::string const& file_name_base,
921  std::vector<Partition> const& partitions,
922  std::vector<long> const& num_elem_integers,
923  std::vector<long> const& num_g_elem_integers)
924 {
925  const std::string npartitions_str = std::to_string(partitions.size());
926 
927  auto const file_name_ele =
928  file_name_base + "_partitioned_msh_ele" + npartitions_str + ".bin";
929  std::ofstream element_info_os(file_name_ele, std::ios::binary);
930  if (!element_info_os)
931  {
932  OGS_FATAL("Could not open file '%s' for output.",
933  file_name_ele.c_str());
934  }
935 
936  auto const file_name_ele_g =
937  file_name_base + "_partitioned_msh_ele_g" + npartitions_str + ".bin";
938  std::ofstream ghost_element_info_os(file_name_ele_g, std::ios::binary);
939  if (!ghost_element_info_os)
940  {
941  OGS_FATAL("Could not open file '%s' for output.",
942  file_name_ele_g.c_str());
943  }
944 
945  for (std::size_t i = 0; i < partitions.size(); i++)
946  {
947  const auto& partition = partitions[i];
948  auto const local_node_ids = enumerateLocalNodeIds(partition.nodes);
949 
950  // A vector contians all element integer variables of
951  // the non-ghost elements of this partition
952  std::vector<long> ele_info(num_elem_integers[i]);
953 
954  // Non-ghost elements.
955  long counter = partition.regular_elements.size();
956 
957  for (std::size_t j = 0; j < partition.regular_elements.size(); j++)
958  {
959  const auto* elem = partition.regular_elements[j];
960  ele_info[j] = counter;
961  getElementIntegerVariables(*elem, local_node_ids, ele_info,
962  counter);
963  }
964  // Write vector data of non-ghost elements
965  element_info_os.write(reinterpret_cast<const char*>(ele_info.data()),
966  ele_info.size() * sizeof(long));
967 
968  // Ghost elements
969  ele_info.resize(num_g_elem_integers[i]);
970 
971  counter = partition.ghost_elements.size();
972 
973  for (std::size_t j = 0; j < partition.ghost_elements.size(); j++)
974  {
975  const auto* elem = partition.ghost_elements[j];
976  ele_info[j] = counter;
977  getElementIntegerVariables(*elem, local_node_ids, ele_info,
978  counter);
979  }
980  // Write vector data of ghost elements
981  ghost_element_info_os.write(
982  reinterpret_cast<const char*>(ele_info.data()),
983  ele_info.size() * sizeof(long));
984  }
985 }
986 
991 void writeNodesBinary(const std::string& file_name_base,
992  std::vector<Partition> const& partitions,
993  std::vector<std::size_t> const& global_node_ids)
994 {
995  auto const file_name = file_name_base + "_partitioned_msh_nod" +
996  std::to_string(partitions.size()) + ".bin";
997  std::ofstream os(file_name, std::ios::binary);
998  if (!os)
999  {
1000  OGS_FATAL("Could not open file '%s' for output.", file_name.c_str());
1001  }
1002 
1003  for (const auto& partition : partitions)
1004  {
1005  partition.writeNodesBinary(os, global_node_ids);
1006  }
1007 }
1008 
1009 void NodeWiseMeshPartitioner::writeBinary(const std::string& file_name_base)
1010 {
1011  writePropertiesBinary(file_name_base, _partitioned_properties, _partitions,
1013  writePropertiesBinary(file_name_base, _partitioned_properties, _partitions,
1015 
1016  const auto elem_integers =
1017  writeConfigDataBinary(file_name_base, _partitions);
1018 
1019  const std::vector<IntegerType>& num_elem_integers =
1020  std::get<0>(elem_integers);
1021  const std::vector<IntegerType>& num_g_elem_integers =
1022  std::get<1>(elem_integers);
1023  writeElementsBinary(file_name_base, _partitions, num_elem_integers,
1024  num_g_elem_integers);
1025 
1026  writeNodesBinary(file_name_base, _partitions, _nodes_global_ids);
1027 }
1028 
1030  std::string const& output_filename_base,
1031  std::vector<Partition> const& partitions,
1032  MeshLib::Properties const& partitioned_properties) const
1033 {
1034  writeNodesBinary(output_filename_base, partitions, _nodes_global_ids);
1035 
1036  const auto elem_integers =
1037  writeConfigDataBinary(output_filename_base, partitions);
1038 
1039  const std::vector<IntegerType>& num_elem_integers =
1040  std::get<0>(elem_integers);
1041  const std::vector<IntegerType>& num_g_elem_integers =
1042  std::get<1>(elem_integers);
1043  writeElementsBinary(output_filename_base, partitions, num_elem_integers,
1044  num_g_elem_integers);
1045 
1046  writePropertiesBinary(output_filename_base, partitioned_properties,
1047  partitions, MeshLib::MeshItemType::Node);
1048  writePropertiesBinary(output_filename_base, partitioned_properties,
1049  partitions, MeshLib::MeshItemType::Cell);
1050 }
1051 
1053  const std::string& file_name_base)
1054 {
1055  const std::string fname = file_name_base + "_partitioned_cfg" +
1056  std::to_string(_npartitions) + ".msh";
1057  std::fstream os_subd_head(fname, std::ios::out | std::ios::trunc);
1058  const std::string mesh_info =
1059  "Subdomain mesh ("
1060  "Number of nodes; Number of base nodes;"
1061  " Number of regular elements; Number of ghost elements;"
1062  " Number of non-ghost base nodes; Number of non-ghost nodes"
1063  " Number of base nodes of the global mesh;"
1064  " Number of nodes of the global mesh;"
1065  " Number of integer variables to define non-ghost elements;"
1066  " Number of integer variables to define ghost elements.)";
1067  os_subd_head << mesh_info << "\n";
1068  os_subd_head << _npartitions << "\n";
1069 
1070  for (const auto& partition : _partitions)
1071  {
1072  os_subd_head << partition.nodes.size();
1073  os_subd_head << " " << partition.number_of_base_nodes;
1074  os_subd_head << " " << partition.regular_elements.size();
1075  os_subd_head << " " << partition.ghost_elements.size();
1076  os_subd_head << " " << partition.number_of_non_ghost_base_nodes;
1077  os_subd_head << " " << partition.number_of_non_ghost_nodes;
1078  os_subd_head << " " << _mesh->getNumberOfBaseNodes();
1079  os_subd_head << " " << _mesh->getNumberOfNodes();
1080  os_subd_head << " "
1082  partition.regular_elements);
1083  os_subd_head << " "
1085  partition.ghost_elements)
1086  << " 0\n";
1087  }
1088 }
1089 
1091  const std::string& file_name_base)
1092 {
1093  const std::string fname = file_name_base + "_partitioned_elems_" +
1094  std::to_string(_npartitions) + ".msh";
1095  std::fstream os_subd(fname, std::ios::out | std::ios::trunc);
1096  for (const auto& partition : _partitions)
1097  {
1098  // Set the local node indices of the current partition.
1099  IntegerType node_local_id_offset = 0;
1100  std::vector<IntegerType> nodes_local_ids(_mesh->getNumberOfNodes(), -1);
1101  for (const auto* node : partition.nodes)
1102  {
1103  nodes_local_ids[node->getID()] = node_local_id_offset;
1104  ++node_local_id_offset;
1105  }
1106 
1107  for (const auto* elem : partition.regular_elements)
1108  {
1109  writeLocalElementNodeIndices(os_subd, *elem, nodes_local_ids);
1110  }
1111  for (const auto* elem : partition.ghost_elements)
1112  {
1113  writeLocalElementNodeIndices(os_subd, *elem, nodes_local_ids);
1114  }
1115  os_subd << std::endl;
1116  }
1117 }
1118 
1119 void NodeWiseMeshPartitioner::writeNodesASCII(const std::string& file_name_base)
1120 {
1121  const std::string fname = file_name_base + "_partitioned_nodes_" +
1122  std::to_string(_npartitions) + ".msh";
1123  std::fstream os_subd_node(fname, std::ios::out | std::ios::trunc);
1124  os_subd_node.precision(std::numeric_limits<double>::digits10);
1125  os_subd_node.setf(std::ios::scientific);
1126 
1127  for (const auto& partition : _partitions)
1128  {
1129  for (const auto* node : partition.nodes)
1130  {
1131  double const* coords = node->getCoords();
1132  os_subd_node << _nodes_global_ids[node->getID()] << " " << coords[0]
1133  << " " << coords[1] << " " << coords[2] << "\n";
1134  }
1135  os_subd_node << std::endl;
1136  }
1137 }
1138 
1139 void NodeWiseMeshPartitioner::writeASCII(const std::string& file_name_base)
1140 {
1141  writeConfigDataASCII(file_name_base);
1142  writeElementsASCII(file_name_base);
1143  writeNodesASCII(file_name_base);
1144 }
1145 
1147  std::ostream& os,
1148  const MeshLib::Element& elem,
1149  const std::vector<IntegerType>& local_node_ids)
1150 {
1151  unsigned mat_id = 0; // TODO: Material ID to be set from the mesh data
1152  os << mat_id << " " << static_cast<unsigned>(elem.getCellType()) << " "
1153  << elem.getNumberOfNodes() << " ";
1154  for (unsigned i = 0; i < elem.getNumberOfNodes(); i++)
1155  {
1156  os << " " << local_node_ids[elem.getNodeIndex(i)];
1157  }
1158  os << "\n";
1159 }
1160 
1161 } // namespace ApplicationUtils
std::pair< std::vector< MeshLib::Node * >, std::vector< MeshLib::Node * > > findNonGhostNodesInPartition(std::size_t const part_id, const bool is_mixed_high_order_linear_elems, std::size_t const n_base_nodes, std::vector< MeshLib::Node *> const &nodes, std::vector< std::size_t > const &partition_ids, std::vector< std::size_t > const *node_id_mapping=nullptr)
std::size_t getNumberOfNodes() const
Get the number of nodes.
Definition: Mesh.h:99
void writePropertiesBinary(const std::string &file_name_base, MeshLib::Properties const &partitioned_properties, std::vector< Partition > const &partitions, MeshLib::MeshItemType const mesh_item_type)
bool writePropertyVectorBinary(MeshLib::Properties const &partitioned_properties, std::string const &name, std::ostream &out_val, std::ostream &out_meta)
NodeStruct(NodeWiseMeshPartitioner::IntegerType const id_, double const x_, double const y_, double const z_)
std::size_t copyNodePropertyVectorValues(Partition const &p, std::size_t const offset, MeshLib::PropertyVector< T > const &pv, MeshLib::PropertyVector< T > &partitioned_pv)
void writePropertyVectorMetaDataBinary(std::ostream &os, PropertyVectorMetaData const &pvmd)
std::size_t getNumberOfBaseNodes() const
Get the number of base nodes.
Definition: Mesh.h:126
PropertyVector< T > * createNewPropertyVector(std::string const &name, MeshItemType mesh_item_type, std::size_t n_components=1)
Definition: Properties.h:15
std::ostream & writeNodesBinary(std::ostream &os, std::vector< std::size_t > const &global_node_ids) const
Declare a class to perform node wise mesh partitioning.
std::vector< std::string > getPropertyVectorNames() const
Definition: Properties.cpp:37
Implementation of the VtuInterface class.
void writePropertyVectorValuesBinary(std::ostream &os, MeshLib::PropertyVector< T > const &pv)
std::size_t partitionLookup(MeshLib::Node const &node, std::vector< std::size_t > const &partition_ids, std::vector< std::size_t > const *node_id_mapping=nullptr)
void writeElementsBinary(std::string const &file_name_base, std::vector< Partition > const &partitions, std::vector< long > const &num_elem_integers, std::vector< long > const &num_g_elem_integers)
void processPartition(std::size_t const part_id, const bool is_mixed_high_order_linear_elems)
void writeNodesASCII(const std::string &file_name_base)
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
Definition: Mesh.h:105
NodeWiseMeshPartitioner::IntegerType getNumberOfIntegerVariablesOfElements(std::vector< const MeshLib::Element *> const &elements)
void writePropertyVectorPartitionMetaData(std::ostream &os, PropertyVectorPartitionMetaData const &pvpmd)
std::ostream & writeConfigBinary(std::ostream &os) const
std::function< bool(const double t, MappedConstVector< N > const &y, MappedVector< N > &ydot)> Function
MeshLib::Properties & getProperties()
Definition: Mesh.h:134
std::tuple< std::vector< MeshLib::Node * >, std::vector< MeshLib::Node * > > findGhostNodesInPartition(std::size_t const part_id, const bool is_mixed_high_order_linear_elems, std::size_t const number_of_base_nodes, std::vector< MeshLib::Node *> const &nodes, std::vector< MeshLib::Element const *> const &ghost_elements, std::vector< std::size_t > const &partition_ids, std::vector< std::size_t > const *node_id_mapping=nullptr)
std::unordered_map< std::size_t, long > enumerateLocalNodeIds(std::vector< MeshLib::Node *> const &nodes)
Generates a mapping of given node ids to a new local (renumbered) node ids.
void applyToPropertyVectors(std::vector< std::string > const &property_names, Function f)
std::vector< const MeshLib::Element * > ghost_elements
std::size_t getID() const
Definition: Point3dWithID.h:62
bool copyPropertyVector(MeshLib::Properties const &original_properties, MeshLib::Properties &partitioned_properties, std::vector< Partition > const &partitions, std::string const &name, std::size_t const total_number_of_tuples)
void renumberBulkNodeIdsProperty(MeshLib::PropertyVector< std::size_t > *const bulk_node_ids, std::vector< Partition > const &local_partitions) const
void writeLocalElementNodeIndices(std::ostream &os, const MeshLib::Element &elem, const std::vector< IntegerType > &local_node_ids)
Write local indices of element nodes to a ASCII file.
Property manager on mesh items. Class Properties manages scalar, vector or matrix properties...
Definition: Properties.h:37
std::size_t nodeIdBulkMesh(MeshLib::Node const &node, std::vector< std::size_t > const *node_id_mapping=nullptr)
MeshLib::Properties partitionProperties(MeshLib::Properties const &properties, std::vector< Partition > const &partitions)
Creates partitioned mesh properties for nodes and cells.
void renumberBulkElementIdsProperty(MeshLib::PropertyVector< std::size_t > *const bulk_element_ids_pv, std::vector< Partition > const &local_partitions) const
std::ptrdiff_t numberOfRegularNodes(MeshLib::Element const &e, std::size_t const part_id, std::vector< std::size_t > const &partition_ids, std::vector< std::size_t > const *node_id_mapping=nullptr)
void getElementIntegerVariables(const MeshLib::Element &elem, const std::unordered_map< std::size_t, long > &local_node_ids, std::vector< long > &elem_info, long &counter)
void writeASCII(const std::string &file_name_base)
MeshItemType
Definition: Location.h:21
std::vector< Partition > partitionOtherMesh(MeshLib::Mesh const &mesh, bool const is_mixed_high_order_linear_elems) const
void writeBinary(const std::string &file_name_base)
void writeElementsASCII(const std::string &file_name_base)
std::size_t getNodeIndex(unsigned i) const
Definition: Element.cpp:182
void partitionByMETIS(const bool is_mixed_high_order_linear_elems)
std::ostream & writeConfigBinary(std::ostream &os) const
virtual CellType getCellType() const =0
std::tuple< std::vector< long >, std::vector< long > > writeConfigDataBinary(const std::string &file_name_base, std::vector< Partition > const &partitions)
std::tuple< std::vector< MeshLib::Element const * >, std::vector< MeshLib::Element const * > > findElementsInPartition(std::size_t const part_id, std::vector< MeshLib::Element *> const &elements, std::vector< std::size_t > const &partition_ids, std::vector< std::size_t > const *node_id_mapping=nullptr)
PropertyVector< T > const * getPropertyVector(std::string const &name) const
Definition: Properties.h:119
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition: Mesh.h:108
bool existsPropertyVector(std::string const &name) const
Definition: Properties.h:79
void renumberNodeIndices(const bool is_mixed_high_order_linear_elems)
void writeNodesBinary(const std::string &file_name_base, std::vector< Partition > const &partitions, std::vector< std::size_t > const &global_node_ids)
std::vector< MeshLib::Node * > nodes
nodes.
void processProperties(MeshLib::Properties const &properties, MeshLib::MeshItemType const mesh_item_type, std::vector< Partition > const &partitions, MeshLib::Properties &partitioned_properties)
#define OGS_FATAL(fmt,...)
Definition: Error.h:63
virtual unsigned getNumberOfNodes() const =0
Returns the number of all nodes including both linear and nonlinear nodes.
std::size_t numberOfMeshItems(MeshLib::MeshItemType const item_type) const
ConfigOffsets incrementConfigOffsets(ConfigOffsets const &oldConfig, PartitionOffsets const &offsets)
void writeValueBinary(std::ostream &out, T const &val)
write value as binary into the given output stream
Definition: FileTools.h:39
void copy(MatrixOrVector const &x, MatrixOrVector &y)
Copies x to y.
Definition: LinAlg.h:36
Node *const * getNodes() const
Get array of element nodes.
Definition: Element.h:78
const char * toString(mgis::behaviour::Behaviour::Kinematic kin)
Converts MGIS kinematic to a string representation.
Definition: MFront.cpp:95
void writeOtherMesh(std::string const &output_filename_base, std::vector< Partition > const &partitions, MeshLib::Properties const &partitioned_properties) const
std::size_t size() const
Method returns the number of tuples times the number of tuple components.
std::size_t copyCellPropertyVectorValues(Partition const &p, std::size_t const offset, MeshLib::PropertyVector< T > const &pv, MeshLib::PropertyVector< T > &partitioned_pv)
std::vector< const MeshLib::Element * > regular_elements
Non ghost elements.
PartitionOffsets computePartitionElementOffsets(Partition const &partition)
NodeWiseMeshPartitioner::IntegerType id
const Node * getNode(unsigned i) const
Definition: Element.cpp:158
void writeConfigDataASCII(const std::string &file_name_base)