OGS 6.2.0-97-g4a610c866
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 (std::size_t elem_id = 0; elem_id < elements.size(); elem_id++)
209  {
210  const auto* elem = elements[elem_id];
211 
212  auto const regular_nodes = numberOfRegularNodes(
213  *elem, part_id, partition_ids, node_id_mapping);
214 
215  if (regular_nodes == 0)
216  {
217  continue;
218  }
219 
220  if (regular_nodes ==
221  static_cast<std::ptrdiff_t>(elem->getNumberOfNodes()))
222  {
223  regular_elements.push_back(elem);
224  }
225  else
226  {
227  ghost_elements.push_back(elem);
228  }
229  }
230  return std::tuple<std::vector<MeshLib::Element const*>,
231  std::vector<MeshLib::Element const*>>{regular_elements,
232  ghost_elements};
233 }
234 
239 std::tuple<std::vector<MeshLib::Node*>, std::vector<MeshLib::Node*>>
241  std::size_t const part_id,
242  const bool is_mixed_high_order_linear_elems,
243  std::size_t const number_of_base_nodes,
244  std::vector<MeshLib::Node*> const& nodes,
245  std::vector<MeshLib::Element const*> const& ghost_elements,
246  std::vector<std::size_t> const& partition_ids,
247  std::vector<std::size_t> const* node_id_mapping = nullptr)
248 {
249  auto node_id = [&node_id_mapping](MeshLib::Node const& n) {
250  return nodeIdBulkMesh(n, node_id_mapping);
251  };
252 
253  std::vector<MeshLib::Node*> base_nodes;
254  std::vector<MeshLib::Node*> ghost_nodes;
255 
256  std::vector<bool> nodes_reserved(nodes.size(), false);
257  for (const auto* ghost_elem : ghost_elements)
258  {
259  for (unsigned i = 0; i < ghost_elem->getNumberOfNodes(); i++)
260  {
261  auto const& n = ghost_elem->getNode(i);
262  if (nodes_reserved[n->getID()])
263  {
264  continue;
265  }
266 
267  if (partitionLookup(*n, partition_ids, node_id_mapping) != part_id)
268  {
269  if (!is_mixed_high_order_linear_elems ||
270  node_id(*n) > number_of_base_nodes)
271  {
272  base_nodes.push_back(nodes[n->getID()]);
273  }
274  else
275  {
276  ghost_nodes.push_back(nodes[n->getID()]);
277  }
278  nodes_reserved[n->getID()] = true;
279  }
280  }
281  }
282  return std::tuple<std::vector<MeshLib::Node*>, std::vector<MeshLib::Node*>>{
283  base_nodes, ghost_nodes};
284 }
285 
287  std::size_t const part_id, const bool is_mixed_high_order_linear_elems)
288 {
289  auto& partition = _partitions[part_id];
290  std::vector<MeshLib::Node*> higher_order_regular_nodes;
291  std::tie(partition.nodes, higher_order_regular_nodes) =
292  findNonGhostNodesInPartition(part_id, is_mixed_high_order_linear_elems,
293  _mesh->getNumberOfBaseNodes(),
294  _mesh->getNodes(), _nodes_partition_ids);
295 
296  partition.number_of_non_ghost_base_nodes = partition.nodes.size();
297  partition.number_of_non_ghost_nodes =
298  partition.number_of_non_ghost_base_nodes +
299  higher_order_regular_nodes.size();
300 
301  std::tie(partition.regular_elements, partition.ghost_elements) =
302  findElementsInPartition(part_id, _mesh->getElements(),
303  _nodes_partition_ids);
304  std::vector<MeshLib::Node*> base_ghost_nodes;
305  std::vector<MeshLib::Node*> higher_order_ghost_nodes;
306  std::tie(base_ghost_nodes, higher_order_ghost_nodes) =
307  findGhostNodesInPartition(part_id, is_mixed_high_order_linear_elems,
308  _mesh->getNumberOfBaseNodes(),
309  _mesh->getNodes(), partition.ghost_elements,
310  _nodes_partition_ids);
311 
312  std::copy(begin(base_ghost_nodes), end(base_ghost_nodes),
313  std::back_inserter(partition.nodes));
314 
315  partition.number_of_base_nodes = partition.nodes.size();
316 
317  if (is_mixed_high_order_linear_elems)
318  {
319  std::copy(begin(higher_order_regular_nodes),
320  end(higher_order_regular_nodes),
321  std::back_inserter(partition.nodes));
322  std::copy(begin(higher_order_ghost_nodes),
323  end(higher_order_ghost_nodes),
324  std::back_inserter(partition.nodes));
325  }
326 
327  // Set the node numbers of base and all mesh nodes.
328  partition.number_of_mesh_base_nodes = _mesh->getNumberOfBaseNodes();
329  partition.number_of_mesh_all_nodes = _mesh->getNumberOfNodes();
330 }
331 
334 template <typename T>
336  Partition const& p,
337  std::size_t const offset,
338  MeshLib::PropertyVector<T> const& pv,
339  MeshLib::PropertyVector<T>& partitioned_pv)
340 {
341  auto const& nodes = p.nodes;
342  auto const nnodes = nodes.size();
343  for (std::size_t i = 0; i < nnodes; ++i)
344  {
345  const auto global_id = nodes[i]->getID();
346  partitioned_pv[offset + i] = pv[global_id];
347  }
348  return nnodes;
349 }
350 
354 template <typename T>
356  Partition const& p,
357  std::size_t const offset,
358  MeshLib::PropertyVector<T> const& pv,
359  MeshLib::PropertyVector<T>& partitioned_pv)
360 {
361  std::size_t const n_regular(p.regular_elements.size());
362  for (std::size_t i = 0; i < n_regular; ++i)
363  {
364  const auto id = p.regular_elements[i]->getID();
365  partitioned_pv[offset + i] = pv[id];
366  }
367 
368  std::size_t const n_ghost(p.ghost_elements.size());
369  for (std::size_t i = 0; i < n_ghost; ++i)
370  {
371  const auto id = p.ghost_elements[i]->getID();
372  partitioned_pv[offset + n_regular + i] = pv[id];
373  }
374  return n_regular + n_ghost;
375 }
376 
377 template <typename T>
378 bool copyPropertyVector(MeshLib::Properties const& original_properties,
379  MeshLib::Properties& partitioned_properties,
380  std::vector<Partition> const& partitions,
381  std::string const& name,
382  std::size_t const total_number_of_tuples)
383 {
384  if (!original_properties.existsPropertyVector<T>(name))
385  {
386  return false;
387  }
388 
389  auto const& pv = original_properties.getPropertyVector<T>(name);
390  auto partitioned_pv = partitioned_properties.createNewPropertyVector<T>(
391  name, pv->getMeshItemType(), pv->getNumberOfComponents());
392  partitioned_pv->resize(total_number_of_tuples *
393  pv->getNumberOfComponents());
394 
395  auto copy_property_vector_values = [&](Partition const& p,
396  std::size_t offset) {
397  if (pv->getMeshItemType() == MeshLib::MeshItemType::Node)
398  {
399  return copyNodePropertyVectorValues(p, offset, *pv,
400  *partitioned_pv);
401  }
402  if (pv->getMeshItemType() == MeshLib::MeshItemType::Cell)
403  {
404  return copyCellPropertyVectorValues(p, offset, *pv,
405  *partitioned_pv);
406  }
407  OGS_FATAL(
408  "Copying of property vector values for mesh item type %s is "
409  "not implemented.",
410  pv->getMeshItemType());
411  };
412 
413  std::size_t position_offset(0);
414  for (auto p : partitions)
415  {
416  position_offset += copy_property_vector_values(p, position_offset);
417  }
418  return true;
419 }
420 
426 template <typename Function>
427 void applyToPropertyVectors(std::vector<std::string> const& property_names,
428  Function f)
429 {
430  for (auto const& name : property_names)
431  {
432  // Open question, why is the 'unsigned long' case not compiling giving
433  // an error "expected '(' for function-style cast or type construction"
434  // with clang-7, and "error C4576: a parenthesized type followed by an
435  // initializer list is a non-standard explicit type conversion syntax"
436  // with MSVC-15.
437  bool success =
438  f(double{}, name) || f(float{}, name) || f(int{}, name) ||
439  f(long{}, name) || f(unsigned{}, name) ||
440  f(static_cast<unsigned long>(0), name) || f(std::size_t{}, name);
441  if (!success)
442  {
443  OGS_FATAL("Could not apply function to PropertyVector '%s'.",
444  name.c_str());
445  }
446  }
447 }
448 
449 void processProperties(MeshLib::Properties const& properties,
450  MeshLib::MeshItemType const mesh_item_type,
451  std::vector<Partition> const& partitions,
452  MeshLib::Properties& partitioned_properties)
453 {
454  std::size_t const total_number_of_tuples =
455  std::accumulate(begin(partitions), end(partitions), 0,
456  [&](std::size_t const sum, Partition const& p) {
457  return sum + p.numberOfMeshItems(mesh_item_type);
458  });
459 
460  DBUG(
461  "total number of tuples define on mesh item type '%d' after "
462  "partitioning: %d ",
463  mesh_item_type, total_number_of_tuples);
464 
465  // 1 create new PV
466  // 2 resize the PV with total_number_of_tuples
467  // 3 copy the values according to the partition info
468  applyToPropertyVectors(properties.getPropertyVectorNames(mesh_item_type),
469  [&](auto type, std::string const& name) {
470  return copyPropertyVector<decltype(type)>(
471  properties, partitioned_properties,
472  partitions, name, total_number_of_tuples);
473  });
474 }
475 
477  MeshLib::Properties const& properties,
478  std::vector<Partition> const& partitions)
479 {
480  MeshLib::Properties partitioned_properties;
481  processProperties(properties, MeshLib::MeshItemType::Node, partitions,
482  partitioned_properties);
483  processProperties(properties, MeshLib::MeshItemType::Cell, partitions,
484  partitioned_properties);
485  return partitioned_properties;
486 }
487 
489  const bool is_mixed_high_order_linear_elems)
490 {
491  for (std::size_t part_id = 0; part_id < _partitions.size(); part_id++)
492  {
493  INFO("Processing partition: %d", part_id);
494  processPartition(part_id, is_mixed_high_order_linear_elems);
495  }
496 
497  renumberNodeIndices(is_mixed_high_order_linear_elems);
498 
499  _partitioned_properties =
500  partitionProperties(_mesh->getProperties(), _partitions);
501 }
502 
504  MeshLib::PropertyVector<std::size_t>* const bulk_node_ids_pv,
505  std::vector<Partition> const& local_partitions) const
506 {
507  if (bulk_node_ids_pv == nullptr)
508  {
509  return;
510  }
511 
512  auto& bulk_node_ids = *bulk_node_ids_pv;
513 
514  std::size_t offset = 0; // offset in property vector for current partition
515 
516  assert(_partitions.size() == local_partitions.size());
517  int const n_partitions = static_cast<int>(_partitions.size());
518  for (int partition_id = 0; partition_id < n_partitions; ++partition_id)
519  {
520  auto const& bulk_partition = _partitions[partition_id];
521  auto const& local_partition = local_partitions[partition_id];
522 
523  // Create global-to-local node id mapping for the bulk partition.
524  auto const& bulk_nodes = bulk_partition.nodes;
525  auto const n_bulk_nodes = bulk_nodes.size();
526  std::map<std::size_t, std::size_t> global_to_local;
527  for (std::size_t local_node_id = 0; local_node_id < n_bulk_nodes;
528  ++local_node_id)
529  {
530  global_to_local[bulk_nodes[local_node_id]->getID()] = local_node_id;
531  }
532 
533  auto const& local_nodes = local_partition.nodes;
534  auto const n_local_nodes = local_nodes.size();
535  for (std::size_t local_node_id = 0; local_node_id < n_local_nodes;
536  ++local_node_id)
537  {
538  bulk_node_ids[offset + local_node_id] =
539  global_to_local[bulk_node_ids[offset + local_node_id]];
540  }
541  offset += n_local_nodes;
542  }
543 }
544 
546  MeshLib::PropertyVector<std::size_t>* const bulk_element_ids_pv,
547  std::vector<Partition> const& local_partitions) const
548 {
549  if (bulk_element_ids_pv == nullptr)
550  {
551  return;
552  }
553 
554  auto& bulk_element_ids = *bulk_element_ids_pv;
555 
556  std::size_t offset = 0; // offset in property vector for current partition
557 
558  assert(_partitions.size() == local_partitions.size());
559  int const n_partitions = static_cast<int>(_partitions.size());
560  for (int partition_id = 0; partition_id < n_partitions; ++partition_id)
561  {
562  auto const& bulk_partition = _partitions[partition_id];
563  auto const& local_partition = local_partitions[partition_id];
564 
565  // Create global-to-local element id mapping for the bulk partition.
566  std::map<std::size_t, std::size_t> global_to_local;
567  auto map_elements =
568  [&global_to_local](
569  std::vector<MeshLib::Element const*> const& elements,
570  std::size_t const offset) {
571  auto const n_elements = elements.size();
572  for (std::size_t e = 0; e < n_elements; ++e)
573  {
574  global_to_local[elements[e]->getID()] = offset + e;
575  }
576  };
577 
578  map_elements(bulk_partition.regular_elements, 0);
579  map_elements(bulk_partition.ghost_elements,
580  bulk_partition.regular_elements.size());
581 
582  // Renumber the local bulk_element_ids map.
583  auto renumber_elements =
584  [&bulk_element_ids, &global_to_local](
585  std::vector<MeshLib::Element const*> const& elements,
586  std::size_t const offset) {
587  auto const n_elements = elements.size();
588  for (std::size_t e = 0; e < n_elements; ++e)
589  {
590  bulk_element_ids[offset + e] =
591  global_to_local[bulk_element_ids[offset + e]];
592  }
593  return n_elements;
594  };
595 
596  offset += renumber_elements(local_partition.regular_elements, offset);
597  offset += renumber_elements(local_partition.ghost_elements, offset);
598  }
599 }
600 
602  MeshLib::Mesh const& mesh,
603  bool const is_mixed_high_order_linear_elems) const
604 {
605  auto const& bulk_node_ids =
606  mesh.getProperties().getPropertyVector<std::size_t>(
607  "bulk_node_ids", MeshLib::MeshItemType::Node, 1);
608 
609  std::vector<Partition> partitions(_partitions.size());
610  for (std::size_t part_id = 0; part_id < _partitions.size(); part_id++)
611  {
612  auto& partition = partitions[part_id];
613  INFO("Processing partition: %d", part_id);
614  {
615  // Set the node numbers of base and all mesh nodes.
616  partition.number_of_mesh_base_nodes = mesh.getNumberOfBaseNodes();
617  partition.number_of_mesh_all_nodes = mesh.getNumberOfNodes();
618 
619  std::vector<MeshLib::Node*> higher_order_regular_nodes;
620  std::tie(partition.nodes, higher_order_regular_nodes) =
622  part_id, is_mixed_high_order_linear_elems,
623  mesh.getNumberOfBaseNodes(), mesh.getNodes(),
624  _nodes_partition_ids, bulk_node_ids);
625 
626  partition.number_of_non_ghost_base_nodes = partition.nodes.size();
627  partition.number_of_non_ghost_nodes =
628  partition.number_of_non_ghost_base_nodes +
629  higher_order_regular_nodes.size();
630 
631  std::tie(partition.regular_elements, partition.ghost_elements) =
632  findElementsInPartition(part_id, mesh.getElements(),
633  _nodes_partition_ids, bulk_node_ids);
634 
635  std::vector<MeshLib::Node*> base_ghost_nodes;
636  std::vector<MeshLib::Node*> higher_order_ghost_nodes;
637  std::tie(base_ghost_nodes, higher_order_ghost_nodes) =
639  part_id, is_mixed_high_order_linear_elems,
640  mesh.getNumberOfBaseNodes(), mesh.getNodes(),
641  partition.ghost_elements, _nodes_partition_ids,
642  bulk_node_ids);
643 
644  std::copy(begin(base_ghost_nodes), end(base_ghost_nodes),
645  std::back_inserter(partition.nodes));
646 
647  partition.number_of_base_nodes = partition.nodes.size();
648 
649  if (is_mixed_high_order_linear_elems)
650  {
651  std::copy(begin(higher_order_regular_nodes),
652  end(higher_order_regular_nodes),
653  std::back_inserter(partition.nodes));
654  std::copy(begin(higher_order_ghost_nodes),
655  end(higher_order_ghost_nodes),
656  std::back_inserter(partition.nodes));
657  }
658  }
659  }
660 
661  return partitions;
662 }
663 
665  const bool is_mixed_high_order_linear_elems)
666 {
667  std::size_t node_global_id_offset = 0;
668  // Renumber the global indices.
669  // -- Base nodes
670  for (auto& partition : _partitions)
671  {
672  for (std::size_t i = 0; i < partition.number_of_non_ghost_base_nodes;
673  i++)
674  {
675  _nodes_global_ids[partition.nodes[i]->getID()] =
676  node_global_id_offset;
677  node_global_id_offset++;
678  }
679  }
680 
681  if (!is_mixed_high_order_linear_elems)
682  {
683  return;
684  }
685 
686  // -- Nodes for high order elements.
687  for (auto& partition : _partitions)
688  {
689  const std::size_t end_id = partition.number_of_base_nodes +
690  partition.number_of_non_ghost_nodes -
691  partition.number_of_non_ghost_base_nodes;
692  for (std::size_t i = partition.number_of_base_nodes; i < end_id; i++)
693  {
694  _nodes_global_ids[partition.nodes[i]->getID()] =
695  node_global_id_offset;
696  node_global_id_offset++;
697  }
698  }
699 }
700 
701 template <typename T>
702 void writePropertyVectorValuesBinary(std::ostream& os,
703  MeshLib::PropertyVector<T> const& pv)
704 {
705  os.write(reinterpret_cast<const char*>(pv.data()), pv.size() * sizeof(T));
706 }
707 
708 template <typename T>
710  MeshLib::Properties const& partitioned_properties, std::string const& name,
711  std::ostream& out_val, std::ostream& out_meta)
712 {
713  if (!partitioned_properties.existsPropertyVector<T>(name))
714  {
715  return false;
716  }
717 
719  pvmd.property_name = name;
720  auto* pv = partitioned_properties.getPropertyVector<T>(name);
722  pvmd.number_of_components = pv->getNumberOfComponents();
723  pvmd.number_of_tuples = pv->getNumberOfTuples();
724  writePropertyVectorValuesBinary(out_val, *pv);
726  return true;
727 }
728 
729 void writePropertiesBinary(const std::string& file_name_base,
730  MeshLib::Properties const& partitioned_properties,
731  std::vector<Partition> const& partitions,
732  MeshLib::MeshItemType const mesh_item_type)
733 {
734  auto const& property_names =
735  partitioned_properties.getPropertyVectorNames(mesh_item_type);
736  if (property_names.empty())
737  {
738  return;
739  }
740 
741  auto const file_name_infix = toString(mesh_item_type);
742 
743  auto const file_name_cfg = file_name_base + "_partitioned_" +
744  file_name_infix + "_properties_cfg" +
745  std::to_string(partitions.size()) + ".bin";
746  std::ofstream out(file_name_cfg, std::ios::binary);
747  if (!out)
748  {
749  OGS_FATAL("Could not open file '%s' for output.",
750  file_name_cfg.c_str());
751  }
752 
753  auto const file_name_val = file_name_base + "_partitioned_" +
754  file_name_infix + "_properties_val" +
755  std::to_string(partitions.size()) + ".bin";
756  std::ofstream out_val(file_name_val, std::ios::binary);
757  if (!out_val)
758  {
759  OGS_FATAL("Could not open file '%s' for output.",
760  file_name_val.c_str());
761  }
762 
763  std::size_t const number_of_properties(property_names.size());
765 
766  applyToPropertyVectors(property_names,
767  [&](auto type, std::string const& name) {
768  return writePropertyVectorBinary<decltype(type)>(
769  partitioned_properties, name, out_val, out);
770  });
771 
772  unsigned long offset = 0;
773  for (const auto& partition : partitions)
774  {
776  offset, static_cast<unsigned long>(
777  partition.numberOfMeshItems(mesh_item_type))};
778  DBUG(
779  "Write meta data for node-based PropertyVector: global offset %d, "
780  "number of tuples %d",
781  pvpmd.offset, pvpmd.number_of_tuples);
783  offset += pvpmd.number_of_tuples;
784  }
785 }
786 
788 {
792 
793  std::ostream& writeConfigBinary(std::ostream& os) const;
794 };
795 
796 std::ostream& ConfigOffsets::writeConfigBinary(std::ostream& os) const
797 {
798  os.write(reinterpret_cast<const char*>(this), sizeof(ConfigOffsets));
799 
800  static long reserved = 0; // Value reserved in the binary format, not used
801  // in the partitioning process.
802  return os.write(reinterpret_cast<const char*>(&reserved), sizeof(long));
803 }
804 
806 {
807  long node;
810 };
811 
813 {
814  return {static_cast<long>(partition.nodes.size()),
815  static_cast<long>(partition.regular_elements.size() +
817  partition.regular_elements)),
818  static_cast<long>(partition.ghost_elements.size() +
820  partition.ghost_elements))};
821 }
822 
824  PartitionOffsets const& offsets)
825 {
826  return {
827  static_cast<long>(oldConfig.node_rank_offset +
828  offsets.node * sizeof(NodeStruct)),
829  // Offset the ending entry of the element integer variales of
830  // the non-ghost elements of this partition in the vector of elem_info.
831  static_cast<long>(oldConfig.element_rank_offset +
832  offsets.regular_elements * sizeof(long)),
833 
834  // Offset the ending entry of the element integer variales of
835  // the ghost elements of this partition in the vector of elem_info.
836  static_cast<long>(oldConfig.ghost_element_rank_offset +
837  offsets.ghost_elements * sizeof(long))};
838 }
839 
845 std::tuple<std::vector<long>, std::vector<long>> writeConfigDataBinary(
846  const std::string& file_name_base, std::vector<Partition> const& partitions)
847 {
848  auto const file_name_cfg = file_name_base + "_partitioned_msh_cfg" +
849  std::to_string(partitions.size()) + ".bin";
850  std::ofstream of_bin_cfg(file_name_cfg, std::ios::binary);
851  if (!of_bin_cfg)
852  {
853  OGS_FATAL("Could not open file '%s' for output.",
854  file_name_cfg.c_str());
855  }
856 
857  std::vector<long> num_elem_integers;
858  num_elem_integers.reserve(partitions.size());
859  std::vector<long> num_g_elem_integers;
860  num_g_elem_integers.reserve(partitions.size());
861 
862  ConfigOffsets config_offsets = {0, 0, 0}; // 0 for first partition.
863  for (const auto& partition : partitions)
864  {
865  partition.writeConfigBinary(of_bin_cfg);
866 
867  config_offsets.writeConfigBinary(of_bin_cfg);
868  auto const& new_offsets = computePartitionElementOffsets(partition);
869  config_offsets = incrementConfigOffsets(config_offsets, new_offsets);
870 
871  num_elem_integers.push_back(new_offsets.regular_elements);
872  num_g_elem_integers.push_back(new_offsets.ghost_elements);
873  }
874 
875  return std::make_tuple(num_elem_integers, num_g_elem_integers);
876 }
877 
886  const MeshLib::Element& elem,
887  const std::unordered_map<std::size_t, long>& local_node_ids,
888  std::vector<long>& elem_info,
889  long& counter)
890 {
891  unsigned mat_id = 0; // TODO: Material ID to be set from the mesh data
892  const long nn = elem.getNumberOfNodes();
893  elem_info[counter++] = mat_id;
894  elem_info[counter++] = static_cast<long>(elem.getCellType());
895  elem_info[counter++] = nn;
896 
897  for (long i = 0; i < nn; i++)
898  {
899  auto const& n = *elem.getNode(i);
900  elem_info[counter++] = local_node_ids.at(n.getID());
901  }
902 }
903 
905 std::unordered_map<std::size_t, long> enumerateLocalNodeIds(
906  std::vector<MeshLib::Node*> const& nodes)
907 {
908  std::unordered_map<std::size_t, long> local_ids;
909  local_ids.reserve(nodes.size());
910 
911  long local_node_id = 0;
912  for (const auto* node : nodes)
913  {
914  local_ids[node->getID()] = local_node_id++;
915  }
916  return local_ids;
917 }
918 
925 void writeElementsBinary(std::string const& file_name_base,
926  std::vector<Partition> const& partitions,
927  std::vector<long> const& num_elem_integers,
928  std::vector<long> const& num_g_elem_integers)
929 {
930  const std::string npartitions_str = std::to_string(partitions.size());
931 
932  auto const file_name_ele =
933  file_name_base + "_partitioned_msh_ele" + npartitions_str + ".bin";
934  std::ofstream element_info_os(file_name_ele, std::ios::binary);
935  if (!element_info_os)
936  {
937  OGS_FATAL("Could not open file '%s' for output.",
938  file_name_ele.c_str());
939  }
940 
941  auto const file_name_ele_g =
942  file_name_base + "_partitioned_msh_ele_g" + npartitions_str + ".bin";
943  std::ofstream ghost_element_info_os(file_name_ele_g, std::ios::binary);
944  if (!ghost_element_info_os)
945  {
946  OGS_FATAL("Could not open file '%s' for output.",
947  file_name_ele_g.c_str());
948  }
949 
950  for (std::size_t i = 0; i < partitions.size(); i++)
951  {
952  const auto& partition = partitions[i];
953  auto const local_node_ids = enumerateLocalNodeIds(partition.nodes);
954 
955  // A vector contians all element integer variables of
956  // the non-ghost elements of this partition
957  std::vector<long> ele_info(num_elem_integers[i]);
958 
959  // Non-ghost elements.
960  long counter = partition.regular_elements.size();
961 
962  for (std::size_t j = 0; j < partition.regular_elements.size(); j++)
963  {
964  const auto* elem = partition.regular_elements[j];
965  ele_info[j] = counter;
966  getElementIntegerVariables(*elem, local_node_ids, ele_info,
967  counter);
968  }
969  // Write vector data of non-ghost elements
970  element_info_os.write(reinterpret_cast<const char*>(ele_info.data()),
971  ele_info.size() * sizeof(long));
972 
973  // Ghost elements
974  ele_info.resize(num_g_elem_integers[i]);
975 
976  counter = partition.ghost_elements.size();
977 
978  for (std::size_t j = 0; j < partition.ghost_elements.size(); j++)
979  {
980  const auto* elem = partition.ghost_elements[j];
981  ele_info[j] = counter;
982  getElementIntegerVariables(*elem, local_node_ids, ele_info,
983  counter);
984  }
985  // Write vector data of ghost elements
986  ghost_element_info_os.write(
987  reinterpret_cast<const char*>(ele_info.data()),
988  ele_info.size() * sizeof(long));
989  }
990 }
991 
996 void writeNodesBinary(const std::string& file_name_base,
997  std::vector<Partition> const& partitions,
998  std::vector<std::size_t> const& global_node_ids)
999 {
1000  auto const file_name = file_name_base + "_partitioned_msh_nod" +
1001  std::to_string(partitions.size()) + ".bin";
1002  std::ofstream os(file_name, std::ios::binary);
1003  if (!os)
1004  {
1005  OGS_FATAL("Could not open file '%s' for output.", file_name.c_str());
1006  }
1007 
1008  for (const auto& partition : partitions)
1009  {
1010  partition.writeNodesBinary(os, global_node_ids);
1011  }
1012 }
1013 
1014 void NodeWiseMeshPartitioner::writeBinary(const std::string& file_name_base)
1015 {
1016  writePropertiesBinary(file_name_base, _partitioned_properties, _partitions,
1018  writePropertiesBinary(file_name_base, _partitioned_properties, _partitions,
1020 
1021  const auto elem_integers =
1022  writeConfigDataBinary(file_name_base, _partitions);
1023 
1024  const std::vector<IntegerType>& num_elem_integers =
1025  std::get<0>(elem_integers);
1026  const std::vector<IntegerType>& num_g_elem_integers =
1027  std::get<1>(elem_integers);
1028  writeElementsBinary(file_name_base, _partitions, num_elem_integers,
1029  num_g_elem_integers);
1030 
1031  writeNodesBinary(file_name_base, _partitions, _nodes_global_ids);
1032 }
1033 
1035  std::string const& output_filename_base,
1036  std::vector<Partition> const& partitions,
1037  MeshLib::Properties const& partitioned_properties) const
1038 {
1039  writeNodesBinary(output_filename_base, partitions, _nodes_global_ids);
1040 
1041  const auto elem_integers =
1042  writeConfigDataBinary(output_filename_base, partitions);
1043 
1044  const std::vector<IntegerType>& num_elem_integers =
1045  std::get<0>(elem_integers);
1046  const std::vector<IntegerType>& num_g_elem_integers =
1047  std::get<1>(elem_integers);
1048  writeElementsBinary(output_filename_base, partitions, num_elem_integers,
1049  num_g_elem_integers);
1050 
1051  writePropertiesBinary(output_filename_base, partitioned_properties,
1052  partitions, MeshLib::MeshItemType::Node);
1053  writePropertiesBinary(output_filename_base, partitioned_properties,
1054  partitions, MeshLib::MeshItemType::Cell);
1055 }
1056 
1058  const std::string& file_name_base)
1059 {
1060  const std::string fname = file_name_base + "_partitioned_cfg" +
1061  std::to_string(_npartitions) + ".msh";
1062  std::fstream os_subd_head(fname, std::ios::out | std::ios::trunc);
1063  const std::string mesh_info =
1064  "Subdomain mesh ("
1065  "Number of nodes; Number of base nodes;"
1066  " Number of regular elements; Number of ghost elements;"
1067  " Number of non-ghost base nodes; Number of non-ghost nodes"
1068  " Number of base nodes of the global mesh;"
1069  " Number of nodes of the global mesh;"
1070  " Number of integer variables to define non-ghost elements;"
1071  " Number of integer variables to define ghost elements.)";
1072  os_subd_head << mesh_info << "\n";
1073  os_subd_head << _npartitions << "\n";
1074 
1075  for (const auto& partition : _partitions)
1076  {
1077  os_subd_head << partition.nodes.size();
1078  os_subd_head << " " << partition.number_of_base_nodes;
1079  os_subd_head << " " << partition.regular_elements.size();
1080  os_subd_head << " " << partition.ghost_elements.size();
1081  os_subd_head << " " << partition.number_of_non_ghost_base_nodes;
1082  os_subd_head << " " << partition.number_of_non_ghost_nodes;
1083  os_subd_head << " " << _mesh->getNumberOfBaseNodes();
1084  os_subd_head << " " << _mesh->getNumberOfNodes();
1085  os_subd_head << " "
1087  partition.regular_elements);
1088  os_subd_head << " "
1090  partition.ghost_elements)
1091  << " 0\n";
1092  }
1093 }
1094 
1096  const std::string& file_name_base)
1097 {
1098  const std::string fname = file_name_base + "_partitioned_elems_" +
1099  std::to_string(_npartitions) + ".msh";
1100  std::fstream os_subd(fname, std::ios::out | std::ios::trunc);
1101  for (const auto& partition : _partitions)
1102  {
1103  // Set the local node indices of the current partition.
1104  IntegerType node_local_id_offset = 0;
1105  std::vector<IntegerType> nodes_local_ids(_mesh->getNumberOfNodes(), -1);
1106  for (const auto* node : partition.nodes)
1107  {
1108  nodes_local_ids[node->getID()] = node_local_id_offset;
1109  ++node_local_id_offset;
1110  }
1111 
1112  for (const auto* elem : partition.regular_elements)
1113  {
1114  writeLocalElementNodeIndices(os_subd, *elem, nodes_local_ids);
1115  }
1116  for (const auto* elem : partition.ghost_elements)
1117  {
1118  writeLocalElementNodeIndices(os_subd, *elem, nodes_local_ids);
1119  }
1120  os_subd << std::endl;
1121  }
1122 }
1123 
1124 void NodeWiseMeshPartitioner::writeNodesASCII(const std::string& file_name_base)
1125 {
1126  const std::string fname = file_name_base + "_partitioned_nodes_" +
1127  std::to_string(_npartitions) + ".msh";
1128  std::fstream os_subd_node(fname, std::ios::out | std::ios::trunc);
1129  os_subd_node.precision(std::numeric_limits<double>::digits10);
1130  os_subd_node.setf(std::ios::scientific);
1131 
1132  for (const auto& partition : _partitions)
1133  {
1134  for (const auto* node : partition.nodes)
1135  {
1136  double const* coords = node->getCoords();
1137  os_subd_node << _nodes_global_ids[node->getID()] << " " << coords[0]
1138  << " " << coords[1] << " " << coords[2] << "\n";
1139  }
1140  os_subd_node << std::endl;
1141  }
1142 }
1143 
1144 void NodeWiseMeshPartitioner::writeASCII(const std::string& file_name_base)
1145 {
1146  writeConfigDataASCII(file_name_base);
1147  writeElementsASCII(file_name_base);
1148  writeNodesASCII(file_name_base);
1149 }
1150 
1152  std::ostream& os,
1153  const MeshLib::Element& elem,
1154  const std::vector<IntegerType>& local_node_ids)
1155 {
1156  unsigned mat_id = 0; // TODO: Material ID to be set from the mesh data
1157  os << mat_id << " " << static_cast<unsigned>(elem.getCellType()) << " "
1158  << elem.getNumberOfNodes() << " ";
1159  for (unsigned i = 0; i < elem.getNumberOfNodes(); i++)
1160  {
1161  os << " " << local_node_ids[elem.getNodeIndex(i)];
1162  }
1163  os << "\n";
1164 }
1165 
1166 } // 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:64
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:176
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:156
void writeConfigDataASCII(const std::string &file_name_base)