OGS 6.2.0-97-g4a610c866
NodePartitionedMeshReader.cpp
Go to the documentation of this file.
1 
17 
18 #include <logog/include/logog.hpp>
19 
20 #ifdef USE_PETSC
21 #include <mpi.h>
22 #endif
23 
24 #include "BaseLib/FileTools.h"
25 #include "BaseLib/RunTime.h"
26 
28 #include "MeshLib/Properties.h"
29 
30 // Check if the value can by converted to given type without overflow.
31 template <typename VALUE, typename TYPE>
32 bool
33 is_safely_convertable(VALUE const& value)
34 {
35  bool const result = value <= std::numeric_limits<TYPE>::max();
36  if (!result)
37  {
38  ERR("The value %d is too large for conversion.", value);
39  ERR("Maximum available size is %d.", std::numeric_limits<TYPE>::max());
40  }
41  return result;
42 }
43 
44 namespace MeshLib
45 {
46 namespace IO
47 {
49  : _mpi_comm(comm)
50 {
51  MPI_Comm_size(_mpi_comm, &_mpi_comm_size);
52  MPI_Comm_rank(_mpi_comm, &_mpi_rank);
53 
55 }
56 
58 {
59  MPI_Type_free(&_mpi_node_type);
60 }
61 
63 {
64  int const count = 2;
65  int blocks[count] = {1, 3};
66  MPI_Datatype types[count] = {MPI_UNSIGNED_LONG, MPI_DOUBLE};
67  MPI_Aint displacements[count] = {0, sizeof(NodeData::index)};
68 
69  MPI_Type_create_struct(count, blocks, displacements, types, &_mpi_node_type);
70  MPI_Type_commit(&_mpi_node_type);
71 }
72 
74  const std::string &file_name_base)
75 {
76  BaseLib::RunTime timer;
77  timer.start();
78 
79  MeshLib::NodePartitionedMesh* mesh = nullptr;
80 
81  // Always try binary file first
82  std::string const fname_new = file_name_base + "_partitioned_msh_cfg" +
83  std::to_string(_mpi_comm_size) + ".bin";
84 
85  if(!BaseLib::IsFileExisting(fname_new)) // doesn't exist binary file.
86  {
87  INFO("Reading ASCII mesh file ...");
88 
89  mesh = readASCII(file_name_base);
90  }
91  else
92  {
93  INFO("Reading binary mesh file ...");
94 
95  mesh = readBinary(file_name_base);
96  }
97 
98  INFO("[time] Reading the mesh took %f s.", timer.elapsed());
99 
100  MPI_Barrier(_mpi_comm);
101 
102  return mesh;
103 }
104 
105 template <typename DATA>
106 bool
108  MPI_Offset offset, MPI_Datatype type, DATA& data) const
109 {
110  // Check container size
111  if (!is_safely_convertable<std::size_t, int>(data.size()))
112  {
113  ERR("The container size is too large for MPI_File_read() call.");
114  return false;
115  }
116 
117  // Open file
118  MPI_File file;
119 
120  char* filename_char = const_cast<char*>(filename.data());
121  int const file_status = MPI_File_open(_mpi_comm, filename_char,
122  MPI_MODE_RDONLY, MPI_INFO_NULL, &file);
123 
124  if(file_status != 0)
125  {
126  ERR("Error opening file %s. MPI error code %d", filename.c_str(), file_status);
127  return false;
128  }
129 
130  // Read data
131  char file_mode[] = "native";
132  MPI_File_set_view(file, offset, type, type, file_mode, MPI_INFO_NULL);
133  // The static cast is checked above.
134  MPI_File_read(file, data.data(), static_cast<int>(data.size()), type,
135  MPI_STATUS_IGNORE);
136  MPI_File_close(&file);
137 
138  return true;
139 }
140 
142  const std::string &file_name_base)
143 {
144  //----------------------------------------------------------------------------------
145  // Read headers
146  const std::string fname_header = file_name_base + "_partitioned_msh_";
147  const std::string fname_num_p_ext = std::to_string(_mpi_comm_size) + ".bin";
148 
149  if (!readBinaryDataFromFile(fname_header + "cfg" + fname_num_p_ext,
150  static_cast<MPI_Offset>(
151  static_cast<unsigned>(_mpi_rank) * sizeof(_mesh_info)),
152  MPI_LONG, _mesh_info))
153  return nullptr;
154 
155  //----------------------------------------------------------------------------------
156  // Read Nodes
157  std::vector<NodeData> nodes(_mesh_info.nodes);
158 
159  if (!readBinaryDataFromFile(fname_header + "nod" + fname_num_p_ext,
160  static_cast<MPI_Offset>(_mesh_info.offset[2]), _mpi_node_type, nodes))
161  return nullptr;
162 
163  std::vector<MeshLib::Node*> mesh_nodes;
164  std::vector<unsigned long> glb_node_ids;
165  setNodes(nodes, mesh_nodes, glb_node_ids);
166 
167  //----------------------------------------------------------------------------------
168  // Read non-ghost elements
169  std::vector<unsigned long> elem_data(
171  if (!readBinaryDataFromFile(fname_header + "ele" + fname_num_p_ext,
172  static_cast<MPI_Offset>(_mesh_info.offset[3]), MPI_LONG, elem_data))
173  return nullptr;
174 
175  std::vector<MeshLib::Element*> mesh_elems(
177  setElements(mesh_nodes, elem_data, mesh_elems);
178 
179  //----------------------------------------------------------------------------------
180  //Read ghost element
181  std::vector<unsigned long> ghost_elem_data(
183 
184  if (!readBinaryDataFromFile(fname_header + "ele_g" + fname_num_p_ext,
185  static_cast<MPI_Offset>(_mesh_info.offset[4]), MPI_LONG, ghost_elem_data))
186  return nullptr;
187 
188  const bool process_ghost = true;
189  setElements(mesh_nodes, ghost_elem_data, mesh_elems, process_ghost);
190 
191  //----------------------------------------------------------------------------------
192  // read the properties
193  MeshLib::Properties p(readPropertiesBinary(file_name_base));
194 
195  return newMesh(BaseLib::extractBaseName(file_name_base), mesh_nodes,
196  glb_node_ids, mesh_elems, p);
197 }
198 
200  const std::string& file_name_base) const
201 {
205  return p;
206 }
207 
209  const std::string& file_name_base, MeshLib::MeshItemType t,
210  MeshLib::Properties& p) const
211 {
212  std::string const item_type =
213  t == MeshLib::MeshItemType::Node ? "node" : "cell";
214  const std::string fname_cfg = file_name_base + "_partitioned_" + item_type +
215  "_properties_cfg" +
216  std::to_string(_mpi_comm_size) + ".bin";
217  std::ifstream is(fname_cfg.c_str(), std::ios::binary | std::ios::in);
218  if (!is)
219  {
220  WARN("Could not open file '%s'.\n"
221  "\tYou can ignore this warning if the mesh does not contain %s-"
222  "wise property data.", fname_cfg.c_str(), item_type.data());
223  return;
224  }
225  std::size_t number_of_properties = 0;
226  is.read(reinterpret_cast<char*>(&number_of_properties), sizeof(std::size_t));
227  std::vector<boost::optional<MeshLib::IO::PropertyVectorMetaData>> vec_pvmd(
228  number_of_properties);
229  for (std::size_t i(0); i < number_of_properties; ++i)
230  {
231  vec_pvmd[i] = MeshLib::IO::readPropertyVectorMetaData(is);
232  if (!vec_pvmd[i])
233  {
234  OGS_FATAL(
235  "Error in NodePartitionedMeshReader::readPropertiesBinary: "
236  "Could not read the meta data for the PropertyVector %d",
237  i);
238  }
239  }
240  for (std::size_t i(0); i < number_of_properties; ++i)
241  {
242  DBUG("[%d] +++++++++++++", _mpi_rank);
244  DBUG("[%d] +++++++++++++", _mpi_rank);
245  }
246  auto pos = is.tellg();
247  auto offset =
248  static_cast<long>(pos) +
249  static_cast<long>(_mpi_rank *
251  is.seekg(offset);
252  boost::optional<MeshLib::IO::PropertyVectorPartitionMetaData> pvpmd(
254  bool pvpmd_read_ok = static_cast<bool>(pvpmd);
255  bool all_pvpmd_read_ok;
256  MPI_Allreduce(&pvpmd_read_ok, &all_pvpmd_read_ok, 1, MPI_C_BOOL, MPI_LOR,
257  _mpi_comm);
258  if (!all_pvpmd_read_ok)
259  {
260  OGS_FATAL(
261  "Error in NodePartitionedMeshReader::readPropertiesBinary: "
262  "Could not read the partition meta data for the mpi process %d",
263  _mpi_rank);
264  }
265  DBUG("[%d] offset in the PropertyVector: %d", _mpi_rank, pvpmd->offset);
266  DBUG("[%d] %d tuples in partition.", _mpi_rank, pvpmd->number_of_tuples);
267  is.close();
268 
269  const std::string fname_val = file_name_base + "_partitioned_" + item_type +
270  "_properties_val" +
271  std::to_string(_mpi_comm_size) + ".bin";
272  is.open(fname_val.c_str(), std::ios::binary | std::ios::in);
273  if (!is)
274  {
275  ERR("Could not open file '%s'\n."
276  "\tYou can ignore this warning if the mesh does not contain %s-"
277  "wise property data.", fname_val.c_str(), item_type.data());
278  }
279 
280  readDomainSpecificPartOfPropertyVectors(vec_pvmd, *pvpmd, t, is, p);
281 }
282 
284  std::vector<boost::optional<MeshLib::IO::PropertyVectorMetaData>> const&
285  vec_pvmd,
288  std::istream& is,
289  MeshLib::Properties& p) const
290 {
291  unsigned long global_offset = 0;
292  std::size_t const number_of_properties = vec_pvmd.size();
293  for (std::size_t i(0); i < number_of_properties; ++i)
294  {
295  DBUG("[%d] global offset: %d, offset within the PropertyVector: %d.",
296  _mpi_rank, global_offset,
297  global_offset +
298  pvpmd.offset * vec_pvmd[i]->data_type_size_in_bytes);
299  if (vec_pvmd[i]->is_int_type)
300  {
301  if (vec_pvmd[i]->is_data_type_signed)
302  {
303  if (vec_pvmd[i]->data_type_size_in_bytes == sizeof(int))
304  createPropertyVectorPart<int>(is, *vec_pvmd[i], pvpmd, t,
305  global_offset, p);
306  if (vec_pvmd[i]->data_type_size_in_bytes == sizeof(long))
307  createPropertyVectorPart<long>(is, *vec_pvmd[i], pvpmd, t,
308  global_offset, p);
309  }
310  else
311  {
312  if (vec_pvmd[i]->data_type_size_in_bytes ==
313  sizeof(unsigned int))
314  createPropertyVectorPart<unsigned int>(
315  is, *vec_pvmd[i], pvpmd, t, global_offset, p);
316  if (vec_pvmd[i]->data_type_size_in_bytes ==
317  sizeof(unsigned long))
318  createPropertyVectorPart<unsigned long>(
319  is, *vec_pvmd[i], pvpmd, t, global_offset, p);
320  }
321  }
322  else
323  {
324  if (vec_pvmd[i]->data_type_size_in_bytes == sizeof(float))
325  createPropertyVectorPart<float>(is, *vec_pvmd[i], pvpmd, t,
326  global_offset, p);
327  if (vec_pvmd[i]->data_type_size_in_bytes == sizeof(double))
328  createPropertyVectorPart<double>(is, *vec_pvmd[i], pvpmd, t,
329  global_offset, p);
330  }
331  global_offset += vec_pvmd[i]->data_type_size_in_bytes *
332  vec_pvmd[i]->number_of_tuples *
333  vec_pvmd[i]->number_of_components;
334  }
335 }
336 
338  std::string const& file_name_base, std::ifstream& is_cfg,
339  std::ifstream& is_node, std::ifstream& is_elem) const
340 {
341  const std::string fname_header = file_name_base + "_partitioned_";
342  const std::string fname_num_p_ext = std::to_string(_mpi_comm_size) + ".msh";
343 
344  { // Configuration.
345  std::string const filename = fname_header + "cfg" + fname_num_p_ext;
346  is_cfg.open(filename);
347 
348  if( !is_cfg.good() )
349  {
350  ERR("Error opening file %s for input.", filename.c_str());
351  return false;
352  }
353 
354  std::string tmp_line;
355  std::getline(is_cfg, tmp_line);
356  int num_parts;
357  is_cfg >> num_parts >> std::ws;
358 
359  if(num_parts != _mpi_comm_size)
360  {
361  ERR("Aborting computation because of number of cores"
362  "/ subdomains mismatch.");
363  return false;
364  }
365  }
366 
367  { // Nodes.
368  std::string const filename = fname_header + "nodes_" + fname_num_p_ext;
369  is_node.open(filename);
370  if( !is_node.good() )
371  {
372  ERR("Error opening file %s for input.", filename.c_str());
373  return false;
374  }
375  }
376 
377  { // Elements.
378  std::string const filename = fname_header + "elems_" + fname_num_p_ext;
379  is_elem.open(filename);
380  if( !is_elem.good() )
381  {
382  ERR("Error opening file %s for input.", filename.c_str());
383  return false;
384  }
385  }
386 
387  return true;
388 }
389 
391  const int part_id, std::vector<MeshLib::Node*> &mesh_nodes,
392  std::vector<unsigned long> &glb_node_ids) const
393 {
394  int const message_tag = 0;
395 
396  // MPI_Send/Recv can handle only int. Check for overflow.
397  if (!is_safely_convertable<unsigned long, int>(_mesh_info.nodes))
398  {
399  ERR("Too large number of nodes to read.");
400  return false;
401  }
402 
403  std::vector<NodeData> nodes(_mesh_info.nodes);
404 
405  if(_mpi_rank == 0)
406  {
407  for(unsigned long k = 0; k < _mesh_info.nodes; k++)
408  {
409  NodeData &node = nodes[k];
410  is_node >> node.index >> node.x >> node.y >> node.z >> std::ws;
411  }
412 
413  if(part_id == 0)
414  setNodes(nodes, mesh_nodes, glb_node_ids);
415  else
416  MPI_Send(nodes.data(), static_cast<int>(_mesh_info.nodes),
417  _mpi_node_type, part_id, message_tag, _mpi_comm);
418  }
419  else if(_mpi_rank == part_id)
420  {
421  MPI_Recv(nodes.data(), static_cast<int>(_mesh_info.nodes),
422  _mpi_node_type, 0, message_tag, _mpi_comm, MPI_STATUS_IGNORE);
423  setNodes(nodes, mesh_nodes, glb_node_ids);
424  }
425 
426  return true;
427 }
428 
430  const int part_id, const std::size_t data_size, const bool process_ghost,
431  const std::vector<MeshLib::Node*> &mesh_nodes,
432  std::vector<MeshLib::Element*> &mesh_elems) const
433 {
434  int const message_tag = 0;
435 
436  // MPI_Send/Recv can handle only int. Check for overflow.
437  if (!is_safely_convertable<std::size_t, int>(data_size))
438  {
439  ERR("Too large number of elements to read.");
440  return false;
441  }
442 
443  std::vector<unsigned long> elem_data(data_size);
444  if(_mpi_rank == 0)
445  {
446  readElementASCII(is_elem, elem_data, process_ghost);
447 
448  if(part_id == 0)
449  {
450  if(!process_ghost)
451  mesh_elems.resize(
453  setElements(mesh_nodes, elem_data, mesh_elems, process_ghost);
454  }
455  else
456  MPI_Send(elem_data.data(), static_cast<int>(data_size), MPI_LONG,
457  part_id, message_tag, _mpi_comm);
458  }
459  else if(_mpi_rank == part_id)
460  {
461  MPI_Recv(elem_data.data(), static_cast<int>(data_size), MPI_LONG,
462  0, message_tag, _mpi_comm, MPI_STATUS_IGNORE);
463 
464  if(!process_ghost)
465  mesh_elems.resize(
467  setElements(mesh_nodes, elem_data, mesh_elems, process_ghost);
468  }
469 
470  return true;
471 }
472 
474  const std::string &file_name_base)
475 {
476  std::ifstream is_cfg;
477  std::ifstream is_node;
478  std::ifstream is_elem;
479 
480  bool file_opened = false;
481  if(_mpi_rank == 0)
482  {
483  file_opened = openASCIIFiles(file_name_base, is_cfg, is_node, is_elem);
484 
485  if(!file_opened)
486  {
487  is_cfg.close();
488  is_node.close();
489  is_elem.close();
490  }
491  }
492 
493  MPI_Bcast(&file_opened, 1, MPI_INT, 0, _mpi_comm);
494 
495  if(!file_opened)
496  return nullptr;
497 
498  MeshLib::NodePartitionedMesh* np_mesh = nullptr;
499  std::vector<MeshLib::Node*> mesh_nodes;
500  std::vector<unsigned long> glb_node_ids;
501  std::vector<MeshLib::Element*> mesh_elems;
502 
503  for(int i = 0; i < _mpi_comm_size; i++)
504  {
505  if(_mpi_rank == 0)
506  {
507  // Read first part into _mesh_info which is equal with the binary
508  // structure.
509  for(unsigned long j = 0; j < 10; ++j)
510  is_cfg >> *(_mesh_info.data() + j);
511  // The last positon is the extra_flag.
512  is_cfg >> _mesh_info.extra_flag;
513  is_cfg >> std::ws;
514  }
515 
516  MPI_Bcast(_mesh_info.data(), static_cast<int>(_mesh_info.size()),
517  MPI_LONG, 0, _mpi_comm);
518 
519  //---------------------------------------------------------------------
520  // Read Nodes
521  if (!readCastNodesASCII(is_node, i, mesh_nodes, glb_node_ids))
522  break;
523 
524  //---------------------------------------------------------------------
525  // Read elements
526  if (!readCastElemsASCII(is_elem, i,
528  false, mesh_nodes, mesh_elems))
529  break;
530 
531  //---------------------------------------------------------------------
532  // Ghost elements
533  if (!readCastElemsASCII(is_elem, i,
535  true, mesh_nodes, mesh_elems))
536  break;
537 
538  if(_mpi_rank == i) {
539  // reading ascii properties is not implemented
540  MeshLib::Properties properties;
541  np_mesh = newMesh(BaseLib::extractBaseName(file_name_base),
542  mesh_nodes, glb_node_ids, mesh_elems, properties);
543  }
544  }
545 
546  if(_mpi_rank == 0)
547  {
548  is_cfg.close();
549  is_node.close();
550  is_elem.close();
551  }
552 
553  MPI_Barrier(_mpi_comm);
554 
555  return np_mesh;
556 }
557 
559  std::string const& mesh_name,
560  std::vector<MeshLib::Node*> const& mesh_nodes,
561  std::vector<unsigned long> const& glb_node_ids,
562  std::vector<MeshLib::Element*> const& mesh_elems,
563  MeshLib::Properties const& properties) const
564 {
565  return new MeshLib::NodePartitionedMesh(
566  mesh_name, mesh_nodes, glb_node_ids, mesh_elems, properties,
570 }
571 
573  std::vector<unsigned long>& elem_data, const bool ghost) const
574 {
575  // Set number of elements.
576  unsigned long const ne =
578  unsigned long id_offset_elem = ne;
579  for(unsigned long j = 0; j < ne; j++)
580  {
581  elem_data[j] = id_offset_elem;
582  ins >> elem_data[id_offset_elem++]; //mat. idx
583  ins >> elem_data[id_offset_elem++]; //type
584  ins >> elem_data[id_offset_elem]; //nnodes
585  unsigned long const nn_e = elem_data[id_offset_elem++];
586  for(unsigned long k = 0; k < nn_e; k++)
587  ins >> elem_data[id_offset_elem++];
588  }
589 }
590 
591 void NodePartitionedMeshReader::setNodes(const std::vector<NodeData> &node_data,
592  std::vector<MeshLib::Node*> &mesh_node,
593  std::vector<unsigned long> &glb_node_ids) const
594 {
595  mesh_node.resize(_mesh_info.nodes);
596  glb_node_ids.resize(_mesh_info.nodes);
597 
598  for(std::size_t i = 0; i < mesh_node.size(); i++)
599  {
600  NodeData const& nd = node_data[i];
601  glb_node_ids[i] = nd.index;
602  mesh_node[i] = new MeshLib::Node(nd.x, nd.y, nd.z, i);
603  }
604 }
605 
607  const std::vector<MeshLib::Node*> &mesh_nodes,
608  const std::vector<unsigned long> &elem_data,
609  std::vector<MeshLib::Element*> &mesh_elems, const bool ghost) const
610 {
611  // Number of elements, ether ghost or regular
612  unsigned long const ne =
614  unsigned long const id_offset_ghost =
615  ghost ? _mesh_info.regular_elements : 0;
616 
617  for(unsigned long i = 0; i < ne; i++)
618  {
619  unsigned long id_offset_elem = elem_data[i];
620 
621  const unsigned mat_idx = static_cast<unsigned>( elem_data[id_offset_elem++] );
622  const unsigned long e_type = elem_data[id_offset_elem++];
623  unsigned long const nnodes = elem_data[id_offset_elem++];
624 
625  MeshLib::Node** elem_nodes = new MeshLib::Node*[nnodes];
626  for(unsigned long k = 0; k < nnodes; k++)
627  elem_nodes[k] = mesh_nodes[ elem_data[id_offset_elem++] ];
628 
629  // The element types below are defined by the mesh_partition tool
630  // available at https://github.com/ufz/mesh_partition .
631  switch (e_type)
632  {
633  case 1:
634  mesh_elems[i + id_offset_ghost] =
635  new MeshLib::Point(elem_nodes, mat_idx);
636  break;
637  case 2:
638  mesh_elems[i + id_offset_ghost] =
639  new MeshLib::Line(elem_nodes, mat_idx);
640  break;
641  case 6:
642  mesh_elems[i + id_offset_ghost] =
643  new MeshLib::Quad(elem_nodes, mat_idx);
644  break;
645  case 11:
646  mesh_elems[i + id_offset_ghost] =
647  new MeshLib::Hex(elem_nodes, mat_idx);
648  break;
649  case 4:
650  mesh_elems[i + id_offset_ghost] =
651  new MeshLib::Tri(elem_nodes, mat_idx);
652  break;
653  case 9:
654  mesh_elems[i + id_offset_ghost] =
655  new MeshLib::Tet(elem_nodes, mat_idx);
656  break;
657  case 14:
658  mesh_elems[i + id_offset_ghost] =
659  new MeshLib::Prism(elem_nodes, mat_idx);
660  break;
661  case 17:
662  mesh_elems[i + id_offset_ghost] =
663  new MeshLib::Pyramid(elem_nodes, mat_idx);
664  break;
665  default:
666  OGS_FATAL(
667  "NodePartitionedMeshReader: construction of element type "
668  "%d is not implemented.",
669  e_type);
670  }
671  }
672 }
673 } // namespace IO
674 } // namespace MeshLib
TemplateElement< MeshLib::TriRule3 > Tri
Definition: Tri.h:26
TemplateElement< MeshLib::QuadRule4 > Quad
Definition: Quad.h:28
bool readBinaryDataFromFile(std::string const &filename, MPI_Offset offset, MPI_Datatype type, DATA &data) const
Parallel reading of a binary file via MPI_File_read, and it is called by readBinary to read files of ...
unsigned long ghost_elements
3: Number of ghost element of a partition,
void readDomainSpecificPartOfPropertyVectors(std::vector< boost::optional< MeshLib::IO::PropertyVectorMetaData >> const &vec_pvmd, MeshLib::IO::PropertyVectorPartitionMetaData const &pvpmd, MeshLib::MeshItemType t, std::istream &is, MeshLib::Properties &p) const
std::string extractBaseName(std::string const &pathname)
Definition: FileTools.cpp:116
MPI_Datatype _mpi_node_type
MPI data type for struct NodeData.
MPI_Comm _mpi_comm
Pointer to MPI communicator.
boost::optional< PropertyVectorMetaData > readPropertyVectorMetaData(std::istream &is)
void readElementASCII(std::ifstream &ins, std::vector< unsigned long > &elem_data, const bool ghost=false) const
Read elements data from ASCII file.
int _mpi_comm_size
Number of processes in the communicator: _mpi_comm.
TemplateElement< MeshLib::TetRule4 > Tet
Definition: Tet.h:25
MeshLib::NodePartitionedMesh * readASCII(const std::string &file_name_base)
Create a NodePartitionedMesh object, read ASCII mesh data, and return a pointer to it...
bool readCastNodesASCII(std::ifstream &is_node, const int part_id, std::vector< MeshLib::Node *> &mesh_nodes, std::vector< unsigned long > &glb_node_ids) const
Read mesh nodes from an ASCII file and cast to the corresponding rank.
void writePropertyVectorMetaData(PropertyVectorMetaData const &pvmd)
bool is_safely_convertable(VALUE const &value)
void setNodes(const std::vector< NodeData > &node_data, std::vector< MeshLib::Node *> &mesh_node, std::vector< unsigned long > &glb_node_ids) const
Set mesh nodes from a tempory array containing node data read from file.
bool openASCIIFiles(std::string const &file_name_base, std::ifstream &is_cfg, std::ifstream &is_node, std::ifstream &is_elem) const
Open ASCII files of node partitioned mesh data.
TemplateElement< MeshLib::HexRule8 > Hex
Definition: Hex.h:25
unsigned long global_nodes
7: Number of all nodes of global mesh,
Property manager on mesh items. Class Properties manages scalar, vector or matrix properties...
Definition: Properties.h:37
bool IsFileExisting(const std::string &strFilename)
Returns true if given file exists. From http://www.techbytes.ca/techbyte103.html. ...
Definition: FileTools.cpp:36
Definition of the RunTime class.
boost::optional< PropertyVectorPartitionMetaData > readPropertyVectorPartitionMetaData(std::istream &is)
unsigned long nodes
0: Number of all nodes of a partition,
MeshLib::NodePartitionedMesh * readBinary(const std::string &file_name_base)
Create a NodePartitionedMesh object, read binary mesh data in the manner of parallel, and return a pointer to it. Four binary files have to been read in this function named as: file_name_base+_partitioned_msh_cfg[number of partitions].bin file_name_base+_partitioned_msh_nod[number of partitions].bin file_name_base+_partitioned_msh_ele[number of partitions].bin file_name_base+_partitioned_msh_ele_g[number of partitions].bin in which, the first file contains an array of integers for the PartitionMeshInfo for all partitions.
void setElements(const std::vector< MeshLib::Node *> &mesh_nodes, const std::vector< unsigned long > &elem_data, std::vector< MeshLib::Element *> &mesh_elems, const bool ghost=false) const
Set mesh elements from a tempory array containing node data read from file.
MeshLib::NodePartitionedMesh * read(const std::string &file_name_base)
Create a NodePartitionedMesh object, read data to it, and return a pointer to it. Data files are eith...
MeshItemType
Definition: Location.h:21
TemplateElement< MeshLib::PyramidRule5 > Pyramid
Definition: Pyramid.h:25
Interface for heuristic search length strategy.
Definition: ProjectData.h:29
MeshLib::Properties readPropertiesBinary(const std::string &file_name_base) const
void registerNodeDataMpiType()
Define MPI data type for NodeData struct.
Declare a class to read node-wise partitioned mesh with MPI functions.
Count the running time.
Definition: RunTime.h:32
MeshLib::NodePartitionedMesh * newMesh(std::string const &mesh_name, std::vector< MeshLib::Node *> const &mesh_nodes, std::vector< unsigned long > const &glb_node_ids, std::vector< MeshLib::Element *> const &mesh_elems, MeshLib::Properties const &properties) const
Create a new mesh of NodePartitionedMesh after reading and processing the data.
TemplateElement< PointRule1 > Point
Definition: Point.h:19
void start()
Start the timer.
Definition: RunTime.h:36
TemplateElement< MeshLib::PrismRule6 > Prism
Definition: Prism.h:25
Definition of the class Properties that implements a container of properties.
bool readCastElemsASCII(std::ifstream &is_elem, const int part_id, const std::size_t data_size, const bool process_ghost, const std::vector< MeshLib::Node *> &mesh_nodes, std::vector< MeshLib::Element *> &mesh_elems) const
Read mesh elements from an ASCII file and cast to the corresponding rank.
#define OGS_FATAL(fmt,...)
Definition: Error.h:63
Filename manipulation routines.
struct MeshLib::IO::NodePartitionedMeshReader::PartitionedMeshInfo _mesh_info
double elapsed() const
Get the elapsed time after started.
Definition: RunTime.h:52
unsigned long active_nodes
5: Number of all active nodes a parition,
TemplateElement< MeshLib::LineRule2 > Line
Definition: Line.h:25