equation_systems_io.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2018 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 #include "libmesh/libmesh_common.h"
21 
22 
23 // C++ Includes
24 #include <cstdio> // for std::sprintf
25 #include <sstream>
26 
27 // Local Includes
30 #include "libmesh/mesh_base.h"
31 #include "libmesh/mesh_tools.h"
32 #include "libmesh/parallel.h"
33 #include "libmesh/xdr_cxx.h"
35 
36 namespace libMesh
37 {
38 
39 // Forward Declarations
40 
41 // Anonymous namespace for implementation details.
42 namespace {
43 std::string local_file_name (const unsigned int processor_id,
44  const std::string & name)
45 {
46  std::string basename(name);
47  char buf[256];
48 
49  if (basename.size() - basename.rfind(".bz2") == 4)
50  {
51  basename.erase(basename.end()-4, basename.end());
52  std::sprintf(buf, "%s.%04u.bz2", basename.c_str(), processor_id);
53  }
54  else if (basename.size() - basename.rfind(".gz") == 3)
55  {
56  basename.erase(basename.end()-3, basename.end());
57  std::sprintf(buf, "%s.%04u.gz", basename.c_str(), processor_id);
58  }
59  else
60  std::sprintf(buf, "%s.%04u", basename.c_str(), processor_id);
61 
62  return std::string(buf);
63 }
64 }
65 
66 
67 
68 
69 // ------------------------------------------------------------
70 // EquationSystem class implementation
71 template <typename InValType>
72 void EquationSystems::read (const std::string & name,
73  const unsigned int read_flags,
74  bool partition_agnostic)
75 {
76  XdrMODE mode = READ;
77  if (name.find(".xdr") != std::string::npos)
78  mode = DECODE;
79  this->read(name, mode, read_flags, partition_agnostic);
80 
81 #ifdef LIBMESH_ENABLE_AMR
82  MeshRefinement mesh_refine(_mesh);
83  mesh_refine.clean_refinement_flags();
84 #endif
85 }
86 
87 
88 
89 template <typename InValType>
90 void EquationSystems::read (const std::string & name,
91  const XdrMODE mode,
92  const unsigned int read_flags,
93  bool partition_agnostic)
94 {
95  // If we have exceptions enabled we can be considerate and try
96  // to read old restart files which contain infinite element
97  // information but do not have the " with infinite elements"
98  // string in the version information.
99 
100  // First try the read the user requested
101  libmesh_try
102  {
103  this->_read_impl<InValType> (name, mode, read_flags, partition_agnostic);
104  }
105 
106  // If that fails, try it again but explicitly request we look for infinite element info
107  libmesh_catch (...)
108  {
109  libMesh::out << "\n*********************************************************************\n"
110  << "READING THE FILE \"" << name << "\" FAILED.\n"
111  << "It is possible this file contains infinite element information,\n"
112  << "but the version string does not contain \" with infinite elements\"\n"
113  << "Let's try this again, but looking for infinite element information...\n"
114  << "*********************************************************************\n"
115  << std::endl;
116 
117  libmesh_try
118  {
119  this->_read_impl<InValType> (name, mode, read_flags | EquationSystems::TRY_READ_IFEMS, partition_agnostic);
120  }
121 
122  // If all that failed, we are out of ideas here...
123  libmesh_catch (...)
124  {
125  libMesh::out << "\n*********************************************************************\n"
126  << "Well, at least we tried!\n"
127  << "Good Luck!!\n"
128  << "*********************************************************************\n"
129  << std::endl;
130  LIBMESH_THROW();
131  }
132  }
133 
134 #ifdef LIBMESH_ENABLE_AMR
135  MeshRefinement mesh_refine(_mesh);
136  mesh_refine.clean_refinement_flags();
137 #endif
138 }
139 
140 
141 
142 template <typename InValType>
143 void EquationSystems::_read_impl (const std::string & name,
144  const XdrMODE mode,
145  const unsigned int read_flags,
146  bool partition_agnostic)
147 {
211  // Set booleans from the read_flags argument
212  const bool read_header = read_flags & EquationSystems::READ_HEADER;
213  const bool read_data = read_flags & EquationSystems::READ_DATA;
214  const bool read_additional_data = read_flags & EquationSystems::READ_ADDITIONAL_DATA;
215  const bool read_legacy_format = read_flags & EquationSystems::READ_LEGACY_FORMAT;
216  const bool try_read_ifems = read_flags & EquationSystems::TRY_READ_IFEMS;
217  const bool read_basic_only = read_flags & EquationSystems::READ_BASIC_ONLY;
218  bool read_parallel_files = false;
219 
220  std::vector<std::pair<std::string, System *>> xda_systems;
221 
222  // This will unzip a file with .bz2 as the extension, otherwise it
223  // simply returns the name if the file need not be unzipped.
224  Xdr io ((this->processor_id() == 0) ? name : "", mode);
225  libmesh_assert (io.reading());
226 
227  {
228  // 1.)
229  // Read the version header.
230  std::string version = "legacy";
231  if (!read_legacy_format)
232  {
233  if (this->processor_id() == 0) io.data(version);
234  this->comm().broadcast(version);
235 
236  // All processors have the version header, if it does not contain
237  // the libMesh_label string then it is a legacy file.
238  const std::string libMesh_label = "libMesh-";
239  std::string::size_type lm_pos = version.find(libMesh_label);
240  if (lm_pos==std::string::npos)
241  {
242  io.close();
243 
244  // Recursively call this read() function but with the
245  // EquationSystems::READ_LEGACY_FORMAT bit set.
246  this->read (name, mode, (read_flags | EquationSystems::READ_LEGACY_FORMAT), partition_agnostic);
247  return;
248  }
249 
250  // Figure out the libMesh version that created this file
251  std::istringstream iss(version.substr(lm_pos + libMesh_label.size()));
252  int ver_major = 0, ver_minor = 0, ver_patch = 0;
253  char dot;
254  iss >> ver_major >> dot >> ver_minor >> dot >> ver_patch;
255  io.set_version(LIBMESH_VERSION_ID(ver_major, ver_minor, ver_patch));
256 
257 
258  read_parallel_files = (version.rfind(" parallel") < version.size());
259 
260  // If requested that we try to read infinite element information,
261  // and the string " with infinite elements" is not in the version,
262  // then tack it on. This is for compatibility reading ifem
263  // files written prior to 11/10/2008 - BSK
264  if (try_read_ifems)
265  if (!(version.rfind(" with infinite elements") < version.size()))
266  version += " with infinite elements";
267 
268  }
269  else
270  libmesh_deprecated();
271 
272  START_LOG("read()","EquationSystems");
273 
274  // 2.)
275  // Read the number of equation systems
276  unsigned int n_sys=0;
277  if (this->processor_id() == 0) io.data (n_sys);
278  this->comm().broadcast(n_sys);
279 
280  for (unsigned int sys=0; sys<n_sys; sys++)
281  {
282  // 3.)
283  // Read the name of the sys-th equation system
284  std::string sys_name;
285  if (this->processor_id() == 0) io.data (sys_name);
286  this->comm().broadcast(sys_name);
287 
288  // 4.)
289  // Read the type of the sys-th equation system
290  std::string sys_type;
291  if (this->processor_id() == 0) io.data (sys_type);
292  this->comm().broadcast(sys_type);
293 
294  if (read_header)
295  this->add_system (sys_type, sys_name);
296 
297  // 5.) - 9.)
298  // Let System::read_header() do the job
299  System & new_system = this->get_system(sys_name);
300  new_system.read_header (io,
301  version,
302  read_header,
303  read_additional_data,
304  read_legacy_format);
305 
306  xda_systems.push_back(std::make_pair(sys_name, &new_system));
307 
308  // If we're only creating "basic" systems, we need to tell
309  // each system that before we call init() later.
310  if (read_basic_only)
311  new_system.set_basic_system_only();
312  }
313  }
314 
315 
316 
317  // Now we are ready to initialize the underlying data
318  // structures. This will initialize the vectors for
319  // storage, the dof_map, etc...
320  if (read_header)
321  this->init();
322 
323  // 10.) & 11.)
324  // Read and set the numeric vector values
325  if (read_data)
326  {
327  // the EquationSystems::read() method should look constant from the mesh
328  // perspective, but we need to assign a temporary numbering to the nodes
329  // and elements in the mesh, which requires that we abuse const_cast
330  if (!read_legacy_format && partition_agnostic)
331  {
332  MeshBase & mesh = const_cast<MeshBase &>(this->get_mesh());
334  }
335 
336  Xdr local_io (read_parallel_files ? local_file_name(this->processor_id(),name) : "", mode);
337 
338  for (auto & pr : xda_systems)
339  if (read_legacy_format)
340  {
341  libmesh_deprecated();
342 #ifdef LIBMESH_ENABLE_DEPRECATED
343  pr.second->read_legacy_data (io, read_additional_data);
344 #endif
345  }
346  else
347  if (read_parallel_files)
348  pr.second->read_parallel_data<InValType> (local_io, read_additional_data);
349  else
350  pr.second->read_serialized_data<InValType> (io, read_additional_data);
351 
352 
353  // Undo the temporary numbering.
354  if (!read_legacy_format && partition_agnostic)
356  }
357 
358  STOP_LOG("read()","EquationSystems");
359 
360  // Localize each system's data
361  this->update();
362 }
363 
364 
365 
366 void EquationSystems::write(const std::string & name,
367  const unsigned int write_flags,
368  bool partition_agnostic) const
369 {
370  XdrMODE mode = WRITE;
371  if (name.find(".xdr") != std::string::npos)
372  mode = ENCODE;
373  this->write(name, mode, write_flags, partition_agnostic);
374 }
375 
376 
377 
378 void EquationSystems::write(const std::string & name,
379  const XdrMODE mode,
380  const unsigned int write_flags,
381  bool partition_agnostic) const
382 {
446  // the EquationSystems::write() method should look constant,
447  // but we need to assign a temporary numbering to the nodes
448  // and elements in the mesh, which requires that we abuse const_cast
449  if (partition_agnostic)
450  {
451  MeshBase & mesh = const_cast<MeshBase &>(this->get_mesh());
453  }
454 
455  // set booleans from write_flags argument
456  const bool write_data = write_flags & EquationSystems::WRITE_DATA;
457  const bool write_additional_data = write_flags & EquationSystems::WRITE_ADDITIONAL_DATA;
458 
459  // always write parallel files if we're instructed to write in
460  // parallel
461  const bool write_parallel_files =
463  // Even if we're on a distributed mesh, we may or may not have a
464  // consistent way of reconstructing the same mesh partitioning
465  // later, but we need the same mesh partitioning if we want to
466  // reread the parallel solution safely, so let's write a serial file
467  // unless specifically requested not to.
468  // ||
469  // // but also write parallel files if we haven't been instructed to
470  // // write in serial and we're on a distributed mesh
471  // (!(write_flags & EquationSystems::WRITE_SERIAL_FILES) &&
472  // !this->get_mesh().is_serial())
473  ;
474 
475  // New scope so that io will close before we try to zip the file
476  {
477  Xdr io((this->processor_id()==0) ? name : "", mode);
478  libmesh_assert (io.writing());
479 
480  LOG_SCOPE("write()", "EquationSystems");
481 
482  const unsigned int proc_id = this->processor_id();
483 
484  unsigned int n_sys = 0;
485  for (auto & pr : _systems)
486  if (!pr.second->hide_output())
487  n_sys++;
488 
489  // set the version number in the Xdr object
490  io.set_version(LIBMESH_VERSION_ID(LIBMESH_MAJOR_VERSION,
491  LIBMESH_MINOR_VERSION,
492  LIBMESH_MICRO_VERSION));
493 
494  // Only write the header information
495  // if we are processor 0.
496  if (proc_id == 0)
497  {
498  std::string comment;
499  char buf[256];
500 
501  // 1.)
502  // Write the version header
503  std::string version("libMesh-" + libMesh::get_io_compatibility_version());
504  if (write_parallel_files) version += " parallel";
505 
506 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
507  version += " with infinite elements";
508 #endif
509  io.data (version, "# File Format Identifier");
510 
511  // 2.)
512  // Write the number of equation systems
513  io.data (n_sys, "# No. of Equation Systems");
514 
515  for (auto & pr : _systems)
516  {
517  // Ignore this system if it has been marked as hidden
518  if (pr.second->hide_output()) continue;
519 
520  // 3.)
521  // Write the name of the sys_num-th system
522  {
523  const unsigned int sys_num = pr.second->number();
524  std::string sys_name = pr.first;
525 
526  comment = "# Name, System No. ";
527  std::sprintf(buf, "%u", sys_num);
528  comment += buf;
529 
530  io.data (sys_name, comment.c_str());
531  }
532 
533  // 4.)
534  // Write the type of system handled
535  {
536  const unsigned int sys_num = pr.second->number();
537  std::string sys_type = pr.second->system_type();
538 
539  comment = "# Type, System No. ";
540  std::sprintf(buf, "%u", sys_num);
541  comment += buf;
542 
543  io.data (sys_type, comment.c_str());
544  }
545 
546  // 5.) - 9.)
547  // Let System::write_header() do the job
548  pr.second->write_header (io, version, write_additional_data);
549  }
550  }
551 
552  // Start from the first system, again,
553  // to write vectors to disk, if wanted
554  if (write_data)
555  {
556  // open a parallel buffer if warranted.
557  Xdr local_io (write_parallel_files ? local_file_name(this->processor_id(),name) : "", mode);
558 
559  for (auto & pr : _systems)
560  {
561  // Ignore this system if it has been marked as hidden
562  if (pr.second->hide_output()) continue;
563 
564  // 10.) + 11.)
565  if (write_parallel_files)
566  pr.second->write_parallel_data (local_io,write_additional_data);
567  else
568  pr.second->write_serialized_data (io,write_additional_data);
569  }
570  }
571  }
572 
573  // the EquationSystems::write() method should look constant,
574  // but we need to undo the temporary numbering of the nodes
575  // and elements in the mesh, which requires that we abuse const_cast
576  if (partition_agnostic)
577  const_cast<MeshBase &>(_mesh).fix_broken_node_and_element_numbering();
578 }
579 
580 
581 
582 // template specialization
583 
584 template void EquationSystems::read<Number> (const std::string & name, const unsigned int read_flags, bool partition_agnostic);
585 template void EquationSystems::read<Number> (const std::string & name, const XdrMODE mode, const unsigned int read_flags, bool partition_agnostic);
586 template void EquationSystems::_read_impl<Number> (const std::string & name, const XdrMODE mode, const unsigned int read_flags, bool partition_agnostic);
587 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
588 template void EquationSystems::read<Real> (const std::string & name, const unsigned int read_flags, bool partition_agnostic);
589 template void EquationSystems::read<Real> (const std::string & name, const XdrMODE mode, const unsigned int read_flags, bool partition_agnostic);
590 template void EquationSystems::_read_impl<Real> (const std::string & name, const XdrMODE mode, const unsigned int read_flags, bool partition_agnostic);
591 #endif
592 
593 } // namespace libMesh
std::string name(const ElemQuality q)
Definition: elem_quality.C:42
bool writing() const
Definition: xdr_cxx.h:117
void data(T &a, const char *comment="")
Definition: xdr_cxx.C:761
const T_sys & get_system(const std::string &name) const
void close()
Definition: xdr_cxx.C:273
MeshBase & mesh
const Parallel::Communicator & comm() const
virtual void fix_broken_node_and_element_numbering()=0
void read(const std::string &name, const XdrMODE, const unsigned int read_flags=(READ_HEADER|READ_DATA), bool partition_agnostic=true)
Base class for Mesh.
Definition: mesh_base.h:77
void write(const std::string &name, const XdrMODE, const unsigned int write_flags=(WRITE_DATA), bool partition_agnostic=true) const
void set_basic_system_only()
Definition: system.h:2097
Responsible for mesh refinement algorithms and data.
std::string get_io_compatibility_version()
void read_header(Xdr &io, const std::string &version, const bool read_header=true, const bool read_additional_data=true, const bool read_legacy_format=false)
Definition: system_io.C:115
Manages consistently variables, degrees of freedom, and coefficient vectors.
Definition: system.h:92
virtual System & add_system(const std::string &system_type, const std::string &name)
void globally_renumber_nodes_and_elements(MeshBase &)
Definition: mesh_tools.C:2424
C++ interface for the XDR (eXternal Data Representation) format.
Definition: xdr_cxx.h:66
void set_version(int ver)
Definition: xdr_cxx.h:159
const MeshBase & get_mesh() const
bool reading() const
Definition: xdr_cxx.h:111
void _read_impl(const std::string &name, const XdrMODE, const unsigned int read_flags, bool partition_agnostic=true)
processor_id_type processor_id() const
OStreamProxy out(std::cout)
void broadcast(T &data, const unsigned int root_id=0) const
std::map< std::string, System * > _systems