libMesh
equation_systems_io.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 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"
20 #include "libmesh/libmesh_logging.h"
21 
22 
23 // C++ Includes
24 #include <cstdio> // for std::sprintf
25 #include <sstream>
26 
27 // Local Includes
28 #include "libmesh/libmesh_version.h"
29 #include "libmesh/equation_systems.h"
30 #include "libmesh/mesh_base.h"
31 #include "libmesh/mesh_tools.h"
32 #include "libmesh/parallel.h"
33 #include "libmesh/xdr_cxx.h"
34 #include "libmesh/mesh_refinement.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  std::vector<std::pair<std::string, System *>>::iterator
339  pos = xda_systems.begin();
340 
341  for (; pos != xda_systems.end(); ++pos)
342  if (read_legacy_format)
343  {
344  libmesh_deprecated();
345 #ifdef LIBMESH_ENABLE_DEPRECATED
346  pos->second->read_legacy_data (io, read_additional_data);
347 #endif
348  }
349  else
350  if (read_parallel_files)
351  pos->second->read_parallel_data<InValType> (local_io, read_additional_data);
352  else
353  pos->second->read_serialized_data<InValType> (io, read_additional_data);
354 
355 
356  // Undo the temporary numbering.
357  if (!read_legacy_format && partition_agnostic)
359  }
360 
361  STOP_LOG("read()","EquationSystems");
362 
363  // Localize each system's data
364  this->update();
365 }
366 
367 
368 
369 void EquationSystems::write(const std::string & name,
370  const unsigned int write_flags,
371  bool partition_agnostic) const
372 {
373  XdrMODE mode = WRITE;
374  if (name.find(".xdr") != std::string::npos)
375  mode = ENCODE;
376  this->write(name, mode, write_flags, partition_agnostic);
377 }
378 
379 
380 
381 void EquationSystems::write(const std::string & name,
382  const XdrMODE mode,
383  const unsigned int write_flags,
384  bool partition_agnostic) const
385 {
449  // the EquationSystems::write() method should look constant,
450  // but we need to assign a temporary numbering to the nodes
451  // and elements in the mesh, which requires that we abuse const_cast
452  if (partition_agnostic)
453  {
454  MeshBase & mesh = const_cast<MeshBase &>(this->get_mesh());
456  }
457 
458  // set booleans from write_flags argument
459  const bool write_data = write_flags & EquationSystems::WRITE_DATA;
460  const bool write_additional_data = write_flags & EquationSystems::WRITE_ADDITIONAL_DATA;
461 
462  // always write parallel files if we're instructed to write in
463  // parallel
464  const bool write_parallel_files =
466  // Even if we're on a distributed mesh, we may or may not have a
467  // consistent way of reconstructing the same mesh partitioning
468  // later, but we need the same mesh partitioning if we want to
469  // reread the parallel solution safely, so let's write a serial file
470  // unless specifically requested not to.
471  // ||
472  // // but also write parallel files if we haven't been instructed to
473  // // write in serial and we're on a distributed mesh
474  // (!(write_flags & EquationSystems::WRITE_SERIAL_FILES) &&
475  // !this->get_mesh().is_serial())
476  ;
477 
478  // New scope so that io will close before we try to zip the file
479  {
480  Xdr io((this->processor_id()==0) ? name : "", mode);
481  libmesh_assert (io.writing());
482 
483  LOG_SCOPE("write()", "EquationSystems");
484 
485  const unsigned int proc_id = this->processor_id();
486 
487  unsigned int n_sys = 0;
488  for (std::map<std::string, System *>::const_iterator pos = _systems.begin();
489  pos != _systems.end(); ++pos)
490  {
491  if (! pos->second->hide_output()) n_sys++;
492  }
493 
494  // set the version number in the Xdr object
495  io.set_version(LIBMESH_VERSION_ID(LIBMESH_MAJOR_VERSION,
496  LIBMESH_MINOR_VERSION,
497  LIBMESH_MICRO_VERSION));
498 
499  // Only write the header information
500  // if we are processor 0.
501  if (proc_id == 0)
502  {
503  std::string comment;
504  char buf[256];
505 
506  // 1.)
507  // Write the version header
508  std::string version("libMesh-" + libMesh::get_io_compatibility_version());
509  if (write_parallel_files) version += " parallel";
510 
511 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
512  version += " with infinite elements";
513 #endif
514  io.data (version, "# File Format Identifier");
515 
516  // 2.)
517  // Write the number of equation systems
518  io.data (n_sys, "# No. of Equation Systems");
519 
520  for (std::map<std::string, System *>::const_iterator pos = _systems.begin();
521  pos != _systems.end(); ++pos)
522  {
523  // Ignore this system if it has been marked as hidden
524  if (pos->second->hide_output()) continue;
525 
526  // 3.)
527  // Write the name of the sys_num-th system
528  {
529  const unsigned int sys_num = pos->second->number();
530  std::string sys_name = pos->first;
531 
532  comment = "# Name, System No. ";
533  std::sprintf(buf, "%u", sys_num);
534  comment += buf;
535 
536  io.data (sys_name, comment.c_str());
537  }
538 
539  // 4.)
540  // Write the type of system handled
541  {
542  const unsigned int sys_num = pos->second->number();
543  std::string sys_type = pos->second->system_type();
544 
545  comment = "# Type, System No. ";
546  std::sprintf(buf, "%u", sys_num);
547  comment += buf;
548 
549  io.data (sys_type, comment.c_str());
550  }
551 
552  // 5.) - 9.)
553  // Let System::write_header() do the job
554  pos->second->write_header (io, version, write_additional_data);
555  }
556  }
557 
558  // Start from the first system, again,
559  // to write vectors to disk, if wanted
560  if (write_data)
561  {
562  // open a parallel buffer if warranted.
563  Xdr local_io (write_parallel_files ? local_file_name(this->processor_id(),name) : "", mode);
564 
565  for (std::map<std::string, System *>::const_iterator pos = _systems.begin();
566  pos != _systems.end(); ++pos)
567  {
568  // Ignore this system if it has been marked as hidden
569  if (pos->second->hide_output()) continue;
570 
571  // 10.) + 11.)
572  if (write_parallel_files)
573  pos->second->write_parallel_data (local_io,write_additional_data);
574  else
575  pos->second->write_serialized_data (io,write_additional_data);
576  }
577  }
578  }
579 
580  // the EquationSystems::write() method should look constant,
581  // but we need to undo the temporary numbering of the nodes
582  // and elements in the mesh, which requires that we abuse const_cast
583  if (partition_agnostic)
584  const_cast<MeshBase &>(_mesh).fix_broken_node_and_element_numbering();
585 }
586 
587 
588 
589 // template specialization
590 
591 template void EquationSystems::read<Number> (const std::string & name, const unsigned int read_flags, bool partition_agnostic);
592 template void EquationSystems::read<Number> (const std::string & name, const XdrMODE mode, const unsigned int read_flags, bool partition_agnostic);
593 template void EquationSystems::_read_impl<Number> (const std::string & name, const XdrMODE mode, const unsigned int read_flags, bool partition_agnostic);
594 #ifdef LIBMESH_USE_COMPLEX_NUMBERS
595 template void EquationSystems::read<Real> (const std::string & name, const unsigned int read_flags, bool partition_agnostic);
596 template void EquationSystems::read<Real> (const std::string & name, const XdrMODE mode, const unsigned int read_flags, bool partition_agnostic);
597 template void EquationSystems::_read_impl<Real> (const std::string & name, const XdrMODE mode, const unsigned int read_flags, bool partition_agnostic);
598 #endif
599 
600 } // namespace libMesh
std::string name(const ElemQuality q)
This function returns a string containing some name for q.
Definition: elem_quality.C:39
void data(T &a, const char *comment="")
Inputs or outputs a single value.
Definition: xdr_cxx.C:761
void write(const std::string &name, const XdrMODE, const unsigned int write_flags=(WRITE_DATA), bool partition_agnostic=true) const
Write the systems to disk using the XDR data format.
ImplicitSystem & sys
void close()
Closes the file if it is open.
Definition: xdr_cxx.C:271
MeshBase & mesh
void clean_refinement_flags()
Sets the refinement flag to Elem::DO_NOTHING for each element in the mesh.
virtual void fix_broken_node_and_element_numbering()=0
There is no reason for a user to ever call this function.
void read(const std::string &name, const XdrMODE, const unsigned int read_flags=(READ_HEADER|READ_DATA), bool partition_agnostic=true)
Read & initialize the systems from disk using the XDR data format.
The libMesh namespace provides an interface to certain functionality in the library.
This is the MeshBase class.
Definition: mesh_base.h:68
MeshBase & _mesh
The mesh data structure.
libmesh_assert(j)
void set_basic_system_only()
Sets the system to be "basic only": i.e.
Definition: system.h:2078
This is the MeshRefinement class.
std::string get_io_compatibility_version()
Specifier for I/O file compatibility features.
XdrMODE
Defines an enum for read/write mode in Xdr format.
Definition: enum_xdr_mode.h:32
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)
Reads the basic data header for this System.
Definition: system_io.C:115
This is the base class for classes which contain information related to any physical process that mig...
Definition: system.h:76
bool reading() const
Definition: xdr_cxx.h:113
bool writing() const
Definition: xdr_cxx.h:119
virtual System & add_system(const std::string &system_type, const std::string &name)
Add the system of type system_type named name to the systems array.
void broadcast(T &data, const unsigned int root_id=0) const
Take a local value and broadcast it to all processors.
void globally_renumber_nodes_and_elements(MeshBase &)
There is no reason for a user to ever call this function.
Definition: mesh_tools.C:1968
This class implements a C++ interface to the XDR (eXternal Data Representation) format.
Definition: xdr_cxx.h:68
void set_version(int ver)
Sets the version of the file that is being read.
Definition: xdr_cxx.h:161
const Parallel::Communicator & comm() const
OStreamProxy out
const MeshBase & get_mesh() const
const T_sys & get_system(const std::string &name) const
virtual void init()
Initialize all the systems.
void _read_impl(const std::string &name, const XdrMODE, const unsigned int read_flags, bool partition_agnostic=true)
Actual read implementation.
processor_id_type processor_id()
Definition: libmesh_base.h:96
void update()
Updates local values for all the systems.
processor_id_type processor_id() const
std::map< std::string, System * > _systems
Data structure holding the systems.