MFEM  v4.3.0
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
datacollection.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2021, Lawrence Livermore National Security, LLC. Produced
2 // at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3 // LICENSE and NOTICE for details. LLNL-CODE-806117.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability visit https://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the BSD-3 license. We welcome feedback and contributions, see file
10 // CONTRIBUTING.md for details.
11 
12 #ifndef MFEM_DATACOLLECTION
13 #define MFEM_DATACOLLECTION
14 
15 #include "../config/config.hpp"
16 #include "gridfunc.hpp"
17 #ifdef MFEM_USE_MPI
18 #include "pgridfunc.hpp"
19 #endif
20 #include <string>
21 #include <map>
22 #include <fstream>
23 
24 namespace mfem
25 {
26 
27 /// Lightweight adaptor over an std::map from strings to pointer to T
28 template<typename T>
30 {
31 public:
32  typedef std::map<std::string, T*> MapType;
33  typedef typename MapType::iterator iterator;
34  typedef typename MapType::const_iterator const_iterator;
35 
36  /// Register field @a field with name @a fname
37  /** Replace existing field associated with @a fname (and optionally
38  delete associated pointer if @a own_data is true) */
39  void Register(const std::string& fname, T* field, bool own_data)
40  {
41  T*& ref = field_map[fname];
42  if (own_data)
43  {
44  delete ref; // if newly allocated -> ref is null -> OK
45  }
46  ref = field;
47  }
48 
49  /// Unregister association between field @a field and name @a fname
50  /** Optionally delete associated pointer if @a own_data is true */
51  void Deregister(const std::string& fname, bool own_data)
52  {
53  iterator it = field_map.find(fname);
54  if ( it != field_map.end() )
55  {
56  if (own_data)
57  {
58  delete it->second;
59  }
60  field_map.erase(it);
61  }
62  }
63 
64  /// Clear all associations between names and fields
65  /** Delete associated pointers when @a own_data is true */
66  void DeleteData(bool own_data)
67  {
68  for (iterator it = field_map.begin(); it != field_map.end(); ++it)
69  {
70  if (own_data)
71  {
72  delete it->second;
73  }
74  it->second = NULL;
75  }
76  }
77 
78  /// Predicate to check if a field is associated with name @a fname
79  bool Has(const std::string& fname) const
80  {
81  return field_map.find(fname) != field_map.end();
82  }
83 
84  /// Get a pointer to the field associated with name @a fname
85  /** @return Pointer to field associated with @a fname or NULL */
86  T* Get(const std::string& fname) const
87  {
88  const_iterator it = field_map.find(fname);
89  return it != field_map.end() ? it->second : NULL;
90  }
91 
92  /// Returns a const reference to the underlying map
93  const MapType& GetMap() const { return field_map; }
94 
95  /// Returns the number of registered fields
96  int NumFields() const { return field_map.size(); }
97 
98  /// Returns a begin iterator to the registered fields
99  iterator begin() { return field_map.begin(); }
100  /// Returns a begin const iterator to the registered fields
101  const_iterator begin() const { return field_map.begin(); }
102 
103  /// Returns an end iterator to the registered fields
104  iterator end() { return field_map.end(); }
105  /// Returns an end const iterator to the registered fields
106  const_iterator end() const { return field_map.end(); }
107 
108  /// Returns an iterator to the field @a fname
109  iterator find(const std::string& fname)
110  { return field_map.find(fname); }
111 
112  /// Returns a const iterator to the field @a fname
113  const_iterator find(const std::string& fname) const
114  { return field_map.find(fname); }
115 
116  /// Clears the map of registered fields without reclaiming memory
117  void clear() { field_map.clear(); }
118 
119 protected:
121 };
122 
123 
124 /** A class for collecting finite element data that is part of the same
125  simulation. Currently, this class groups together grid functions (fields),
126  quadrature functions (q-fields), and the mesh that they are defined on. */
128 {
129 private:
130  /// A collection of named GridFunctions
132 
133  /// A collection of named QuadratureFunctions
135 public:
139 
143 
144  /// Format constants to be used with SetFormat().
145  /** Derived classes can define their own format enumerations and override the
146  method SetFormat() to perform input validation. */
147  enum Format
148  {
149  SERIAL_FORMAT = 0, /**<
150  MFEM's serial ascii format, using the methods Mesh::Print() /
151  ParMesh::Print(), and GridFunction::Save() / ParGridFunction::Save().*/
153  MFEM's parallel ascii format, using the methods ParMesh::ParPrint() and
154  GridFunction::Save() / ParGridFunction::Save(). */
155  };
156 
157 protected:
158  /// Name of the collection, used as a directory name when saving
159  std::string name;
160 
161  /** @brief A path where the directory with results is saved.
162  If not empty, it has '/' at the end. */
163  std::string prefix_path;
164 
165  /** A FieldMap mapping registered field names to GridFunction pointers. */
167 
168  /** A FieldMap mapping registered names to QuadratureFunction pointers. */
170 
171  /// The (common) mesh for the collected fields
173 
174  /// Time cycle; for time-dependent simulations cycle >= 0, otherwise = -1.
175  /** When cycle >= 0, it is appended to directory names. */
176  int cycle;
177  /// Physical time (for time-dependent simulations)
178  double time;
179 
180  /// Time step i.e. delta_t (for time-dependent simulations)
181  double time_step;
182 
183  /// Serial or parallel run? False iff mesh is a ParMesh
184  bool serial;
185  /// Append rank to any output file names.
187 
188  /// MPI rank (in parallel)
189  int myid;
190  /// Number of MPI ranks (in parallel)
192 #ifdef MFEM_USE_MPI
193  /// Associated MPI communicator
194  MPI_Comm m_comm;
195 #endif
196 
197  /// Precision (number of digits) used for the text output of doubles
199  /// Number of digits used for the cycle and MPI rank in filenames
201 
202  /// Default value for precision
203  static const int precision_default = 6;
204  /// Default value for pad_digits_*
205  static const int pad_digits_default = 6;
206 
207  /// Output mesh format: see the #Format enumeration
208  int format;
210 
211  /// Should the collection delete its mesh and fields
212  bool own_data;
213 
214  /// Error state
215  int error;
216 
217  /// Delete data owned by the DataCollection keeping field information
218  void DeleteData();
219  /// Delete data owned by the DataCollection including field information
220  void DeleteAll();
221 
222  std::string GetMeshShortFileName() const;
223  std::string GetMeshFileName() const;
224  std::string GetFieldFileName(const std::string &field_name) const;
225 
226  /// Save one field to disk, assuming the collection directory exists
227  void SaveOneField(const FieldMapIterator &it);
228 
229  /// Save one q-field to disk, assuming the collection directory exists
230  void SaveOneQField(const QFieldMapIterator &it);
231 
232  // Helper method
233  static int create_directory(const std::string &dir_name,
234  const Mesh *mesh, int myid);
235 
236 public:
237  /// Initialize the collection with its name and Mesh.
238  /** When @a mesh_ is NULL, then the real mesh can be set with SetMesh(). */
239  explicit DataCollection(const std::string& collection_name,
240  Mesh *mesh_ = NULL);
241 
242  /// Add a grid function to the collection
243  virtual void RegisterField(const std::string& field_name, GridFunction *gf)
244  { field_map.Register(field_name, gf, own_data); }
245 
246  /// Remove a grid function from the collection
247  virtual void DeregisterField(const std::string& field_name)
248  { field_map.Deregister(field_name, own_data); }
249 
250  /// Add a QuadratureFunction to the collection.
251  virtual void RegisterQField(const std::string& q_field_name,
252  QuadratureFunction *qf)
253  { q_field_map.Register(q_field_name, qf, own_data); }
254 
255 
256  /// Remove a QuadratureFunction from the collection
257  virtual void DeregisterQField(const std::string& field_name)
258  { q_field_map.Deregister(field_name, own_data); }
259 
260  /// Check if a grid function is part of the collection
261  bool HasField(const std::string& name) const
262  { return field_map.Has(name); }
263 
264  /// Get a pointer to a grid function in the collection.
265  /** Returns NULL if @a field_name is not in the collection. */
266  GridFunction *GetField(const std::string& field_name)
267  { return field_map.Get(field_name); }
268 
269 #ifdef MFEM_USE_MPI
270  /// Return the associated MPI communicator or MPI_COMM_NULL.
271  MPI_Comm GetComm() const { return m_comm; }
272 
273  /// Get a pointer to a parallel grid function in the collection.
274  /** Returns NULL if @a field_name is not in the collection.
275  @note The GridFunction pointer stored in the collection is statically
276  cast to ParGridFunction pointer. */
277  ParGridFunction *GetParField(const std::string& field_name)
278  { return static_cast<ParGridFunction*>(GetField(field_name)); }
279 #endif
280 
281  /// Check if a QuadratureFunction with the given name is in the collection.
282  bool HasQField(const std::string& q_field_name) const
283  { return q_field_map.Has(q_field_name); }
284 
285  /// Get a pointer to a QuadratureFunction in the collection.
286  /** Returns NULL if @a field_name is not in the collection. */
287  QuadratureFunction *GetQField(const std::string& q_field_name)
288  { return q_field_map.Get(q_field_name); }
289 
290  /// Get a const reference to the internal field map.
291  /** The keys in the map are the field names and the values are pointers to
292  GridFunction%s. */
293  const FieldMapType &GetFieldMap() const
294  { return field_map.GetMap(); }
295 
296  /// Get a const reference to the internal q-field map.
297  /** The keys in the map are the q-field names and the values are pointers to
298  QuadratureFunction%s. */
300  { return q_field_map.GetMap(); }
301 
302  /// Get a pointer to the mesh in the collection
303  Mesh *GetMesh() { return mesh; }
304  /// Set/change the mesh associated with the collection
305  /** When passed a Mesh, assumes the serial case: MPI rank id is set to 0 and
306  MPI num_procs is set to 1. When passed a ParMesh, MPI info from the
307  ParMesh is used to set the DataCollection's MPI rank and num_procs. */
308  virtual void SetMesh(Mesh *new_mesh);
309 #ifdef MFEM_USE_MPI
310  /// Set/change the mesh associated with the collection.
311  /** For this case, @a comm is used to set the DataCollection's MPI rank id
312  and MPI num_procs, which influences the how files are saved for domain
313  decomposed meshes. */
314  virtual void SetMesh(MPI_Comm comm, Mesh *new_mesh);
315 #endif
316 
317  /// Set time cycle (for time-dependent simulations)
318  void SetCycle(int c) { cycle = c; }
319  /// Set physical time (for time-dependent simulations)
320  void SetTime(double t) { time = t; }
321 
322  /// Set the simulation time step (for time-dependent simulations)
323  void SetTimeStep(double ts) { time_step = ts; }
324 
325  /// Get time cycle (for time-dependent simulations)
326  int GetCycle() const { return cycle; }
327  /// Get physical time (for time-dependent simulations)
328  double GetTime() const { return time; }
329  /// Get the simulation time step (for time-dependent simulations)
330  double GetTimeStep() const { return time_step; }
331 
332  /// Get the name of the collection
333  const std::string& GetCollectionName() const { return name; }
334  /// Set the ownership of collection data
335  void SetOwnData(bool o) { own_data = o; }
336 
337  /// Set the precision (number of digits) used for the text output of doubles
338  void SetPrecision(int prec) { precision = prec; }
339  /// Set the number of digits used for both the cycle and the MPI rank
340  void SetPadDigits(int digits) { pad_digits_cycle=pad_digits_rank = digits; }
341  /// Set the number of digits used for the cycle
342  void SetPadDigitsCycle(int digits) { pad_digits_cycle = digits; }
343  /// Set the number of digits used for the MPI rank in filenames
344  void SetPadDigitsRank(int digits) { pad_digits_rank = digits; }
345  /// Set the desired output mesh and data format.
346  /** See the enumeration #Format for valid options. Derived classes can define
347  their own format enumerations and override this method to perform input
348  validation. */
349  virtual void SetFormat(int fmt);
350 
351  /// Set the flag for use of gz compressed files
352  virtual void SetCompression(bool comp);
353 
354  /// Set the path where the DataCollection will be saved.
355  void SetPrefixPath(const std::string &prefix);
356 
357  /// Get the path where the DataCollection will be saved.
358  const std::string &GetPrefixPath() const { return prefix_path; }
359 
360  /// Save the collection to disk.
361  /** By default, everything is saved in the "prefix_path" directory with
362  subdirectory name "collection_name" or "collection_name_cycle" for
363  time-dependent simulations. */
364  virtual void Save();
365  /// Save the mesh, creating the collection directory.
366  virtual void SaveMesh();
367  /// Save one field, assuming the collection directory already exists.
368  virtual void SaveField(const std::string &field_name);
369  /// Save one q-field, assuming the collection directory already exists.
370  virtual void SaveQField(const std::string &q_field_name);
371 
372  /// Load the collection. Not implemented in the base class DataCollection.
373  virtual void Load(int cycle_ = 0);
374 
375  /// Delete the mesh and fields if owned by the collection
376  virtual ~DataCollection();
377 
378  /// Errors returned by Error()
379  enum { NO_ERROR = 0, READ_ERROR = 1, WRITE_ERROR = 2 };
380 
381  /// Get the current error state
382  int Error() const { return error; }
383  /// Reset the error state
384  void ResetError(int err = NO_ERROR) { error = err; }
385 
386 #ifdef MFEM_USE_MPI
387  friend class ParMesh;
388 #endif
389 };
390 
391 
392 /// Helper class for VisIt visualization data
394 {
395 public:
396  std::string association;
398  int lod;
400  VisItFieldInfo(std::string association_, int num_components_, int lod_ = 1)
401  { association = association_; num_components = num_components_; lod =lod_;}
402 };
403 
404 /// Data collection with VisIt I/O routines
406 {
407 protected:
408  // Additional data needed in the VisIt root file, which describes the mesh
409  // and all the fields in the collection
413  std::map<std::string, VisItFieldInfo> field_info_map;
414  typedef std::map<std::string, VisItFieldInfo>::iterator FieldInfoMapIterator;
415 
416  /// Prepare the VisIt root file in JSON format for the current collection
417  std::string GetVisItRootString();
418  /// Read in a VisIt root file in JSON format
419  void ParseVisItRootString(const std::string& json);
420 
421  void UpdateMeshInfo();
422 
423  // Helper functions for Load()
424  void LoadVisItRootFile(const std::string& root_name);
425  void LoadMesh();
426  void LoadFields();
427 
428 public:
429  /// Constructor. The collection name is used when saving the data.
430  /** If @a mesh_ is NULL, then the mesh can be set later by calling either
431  SetMesh() or Load(). The latter works only in serial. */
432  VisItDataCollection(const std::string& collection_name, Mesh *mesh_ = NULL);
433 
434 #ifdef MFEM_USE_MPI
435  /// Construct a parallel VisItDataCollection to be loaded from files.
436  /** Before loading the collection with Load(), some parameters in the
437  collection can be adjusted, e.g. SetPadDigits(), SetPrefixPath(), etc. */
438  VisItDataCollection(MPI_Comm comm, const std::string& collection_name,
439  Mesh *mesh_ = NULL);
440 #endif
441 
442  /// Set/change the mesh associated with the collection
443  virtual void SetMesh(Mesh *new_mesh);
444 
445 #ifdef MFEM_USE_MPI
446  /// Set/change the mesh associated with the collection.
447  virtual void SetMesh(MPI_Comm comm, Mesh *new_mesh);
448 #endif
449 
450  /// Add a grid function to the collection and update the root file
451  virtual void RegisterField(const std::string& field_name, GridFunction *gf);
452 
453  /// Add a quadrature function to the collection and update the root file.
454  /** Visualization of quadrature function is not supported in VisIt(3.12).
455  A patch has been sent to VisIt developers in June 2020. */
456  virtual void RegisterQField(const std::string& q_field_name,
457  QuadratureFunction *qf);
458 
459  /// Set VisIt parameter: default levels of detail for the MultiresControl
460  void SetLevelsOfDetail(int levels_of_detail);
461 
462  /// Set VisIt parameter: maximum levels of detail for the MultiresControl
463  void SetMaxLevelsOfDetail(int max_levels_of_detail);
464 
465  /** @brief Delete all data owned by VisItDataCollection including field data
466  information. */
467  void DeleteAll();
468 
469  /// Save the collection and a VisIt root file
470  virtual void Save();
471 
472  /// Save a VisIt root file for the collection
473  void SaveRootFile();
474 
475  /// Load the collection based on its VisIt data (described in its root file)
476  virtual void Load(int cycle_ = 0);
477 
478  /// We will delete the mesh and fields if we own them
479  virtual ~VisItDataCollection() {}
480 };
481 
482 
483 /// Helper class for ParaView visualization data
485 {
486 private:
487  int levels_of_detail;
488  std::fstream pvd_stream;
489  VTKFormat pv_data_format;
490  bool high_order_output;
491  bool restart_mode;
492 
493 protected:
494  void SaveDataVTU(std::ostream &out, int ref);
495  void SaveGFieldVTU(std::ostream& out, int ref_, const FieldMapIterator& it);
496  void SaveQFieldVTU(std::ostream &out, int ref, const QFieldMapIterator& it);
497  const char *GetDataFormatString() const;
498  const char *GetDataTypeString() const;
499 
500  std::string GenerateCollectionPath();
501  std::string GenerateVTUFileName();
502  std::string GenerateVTUFileName(int rank);
503  std::string GenerateVTUPath();
504  std::string GeneratePVDFileName();
505  std::string GeneratePVTUFileName();
506  std::string GeneratePVTUPath();
507 
508 public:
509  /// Constructor. The collection name is used when saving the data.
510  /** If @a mesh_ is NULL, then the mesh can be set later by calling SetMesh().
511  Before saving the data collection, some parameters in the collection can
512  be adjusted, e.g. SetPadDigits(), SetPrefixPath(), etc. */
513  ParaViewDataCollection(const std::string& collection_name,
514  mfem::Mesh *mesh_ = NULL);
515 
516  /// Set refinement levels - every element is uniformly split based on
517  /// levels_of_detail_
518  void SetLevelsOfDetail(int levels_of_detail_);
519 
520  /// Save the collection - the directory name is constructed based on the
521  /// cycle value
522  virtual void Save() override;
523 
524  /// Set the data format for the ParaView output files. Possible options are
525  /// VTKFormat::ASCII, VTKFormat::BINARY, and VTKFormat::BINARY32.
526  /// The ASCII and BINARY options output double precision data, whereas the
527  /// BINARY32 option outputs single precision data.
528  void SetDataFormat(VTKFormat fmt);
529 
530  /// Set the zlib compression level. 0 indicates no compression, -1 indicates
531  /// the default compression level. Otherwise, specify a number between 1 and
532  /// 9, 1 being the fastest, and 9 being the best compression. Compression
533  /// only takes effect if the output format is BINARY or BINARY32. MFEM must
534  /// be compiled with MFEM_USE_ZLIB = YES.
535  void SetCompressionLevel(int compression_level_);
536 
537  /// Enable or disable zlib compression. If the input is true, use the default
538  /// zlib compression level (unless the compression level has previously been
539  /// set by calling SetCompressionLevel).
540  void SetCompression(bool compression_) override;
541 
542  /// Returns true if the output format is BINARY or BINARY32, false if ASCII.
543  bool IsBinaryFormat() const;
544 
545  /// Sets whether or not to output the data as high-order elements (false
546  /// by default). Reading high-order data requires ParaView 5.5 or later.
547  void SetHighOrderOutput(bool high_order_output_);
548 
549  /// Enable or disable restart mode. If restart is enabled, new writes will
550  /// preserve timestep metadata for any solutions prior to the currently
551  /// defined time.
552  void UseRestartMode(bool restart_mode_);
553 
554  /// Load the collection - not implemented in the ParaView writer
555  virtual void Load(int cycle_ = 0) override;
556 };
557 
558 }
559 #endif
virtual void SaveMesh()
Save the mesh, creating the collection directory.
void SaveOneField(const FieldMapIterator &it)
Save one field to disk, assuming the collection directory exists.
void SetPrecision(int prec)
Set the precision (number of digits) used for the text output of doubles.
Class for grid function - Vector with associated FE space.
Definition: gridfunc.hpp:30
const MapType & GetMap() const
Returns a const reference to the underlying map.
void SetCycle(int c)
Set time cycle (for time-dependent simulations)
int format
Output mesh format: see the Format enumeration.
int error
Error state.
void SetDataFormat(VTKFormat fmt)
double GetTime() const
Get physical time (for time-dependent simulations)
bool appendRankToFileName
Append rank to any output file names.
iterator begin()
Returns a begin iterator to the registered fields.
virtual void RegisterQField(const std::string &q_field_name, QuadratureFunction *qf)
Add a QuadratureFunction to the collection.
void DeleteAll()
Delete data owned by the DataCollection including field information.
const char * GetDataTypeString() const
bool HasQField(const std::string &q_field_name) const
Check if a QuadratureFunction with the given name is in the collection.
Helper class for ParaView visualization data.
const_iterator find(const std::string &fname) const
Returns a const iterator to the field fname.
iterator end()
Returns an end iterator to the registered fields.
int pad_digits_cycle
Number of digits used for the cycle and MPI rank in filenames.
GFieldMap::iterator FieldMapIterator
GridFunction * GetField(const std::string &field_name)
Get a pointer to a grid function in the collection.
Format
Format constants to be used with SetFormat().
const std::string & GetCollectionName() const
Get the name of the collection.
Helper class for VisIt visualization data.
bool HasField(const std::string &name) const
Check if a grid function is part of the collection.
virtual void Save()
Save the collection and a VisIt root file.
MPI_Comm GetComm() const
Return the associated MPI communicator or MPI_COMM_NULL.
std::string GetVisItRootString()
Prepare the VisIt root file in JSON format for the current collection.
QuadratureFunction * GetQField(const std::string &q_field_name)
Get a pointer to a QuadratureFunction in the collection.
QFieldMap::MapType QFieldMapType
bool Has(const std::string &fname) const
Predicate to check if a field is associated with name fname.
int Error() const
Get the current error state.
bool own_data
Should the collection delete its mesh and fields.
VisItFieldInfo(std::string association_, int num_components_, int lod_=1)
virtual void SetCompression(bool comp)
Set the flag for use of gz compressed files.
void ParseVisItRootString(const std::string &json)
Read in a VisIt root file in JSON format.
virtual void SaveQField(const std::string &q_field_name)
Save one q-field, assuming the collection directory already exists.
void SetPadDigitsCycle(int digits)
Set the number of digits used for the cycle.
double GetTimeStep() const
Get the simulation time step (for time-dependent simulations)
MPI_Comm m_comm
Associated MPI communicator.
virtual void Load(int cycle_=0) override
Load the collection - not implemented in the ParaView writer.
double time
Physical time (for time-dependent simulations)
MapType::iterator iterator
void ResetError(int err=NO_ERROR)
Reset the error state.
std::map< std::string, VisItFieldInfo >::iterator FieldInfoMapIterator
iterator find(const std::string &fname)
Returns an iterator to the field fname.
void SetCompressionLevel(int compression_level_)
std::string GetMeshShortFileName() const
MapType::const_iterator const_iterator
Mesh * mesh
The (common) mesh for the collected fields.
const QFieldMapType & GetQFieldMap() const
Get a const reference to the internal q-field map.
void SaveOneQField(const QFieldMapIterator &it)
Save one q-field to disk, assuming the collection directory exists.
void SetPadDigits(int digits)
Set the number of digits used for both the cycle and the MPI rank.
void SaveQFieldVTU(std::ostream &out, int ref, const QFieldMapIterator &it)
virtual void RegisterField(const std::string &field_name, GridFunction *gf)
Add a grid function to the collection.
bool serial
Serial or parallel run? False iff mesh is a ParMesh.
int cycle
Time cycle; for time-dependent simulations cycle &gt;= 0, otherwise = -1.
VisItDataCollection(const std::string &collection_name, Mesh *mesh_=NULL)
Constructor. The collection name is used when saving the data.
virtual void Save()
Save the collection to disk.
const FieldMapType & GetFieldMap() const
Get a const reference to the internal field map.
int GetCycle() const
Get time cycle (for time-dependent simulations)
Data collection with VisIt I/O routines.
void SaveDataVTU(std::ostream &out, int ref)
VTKFormat
Definition: vtk.hpp:61
void SaveGFieldVTU(std::ostream &out, int ref_, const FieldMapIterator &it)
virtual void SetFormat(int fmt)
Set the desired output mesh and data format.
static int create_directory(const std::string &dir_name, const Mesh *mesh, int myid)
void Register(const std::string &fname, T *field, bool own_data)
Register field field with name fname.
DataCollection(const std::string &collection_name, Mesh *mesh_=NULL)
Initialize the collection with its name and Mesh.
void DeleteData()
Delete data owned by the DataCollection keeping field information.
void SetLevelsOfDetail(int levels_of_detail)
Set VisIt parameter: default levels of detail for the MultiresControl.
void SetHighOrderOutput(bool high_order_output_)
static const int pad_digits_default
Default value for pad_digits_*.
void SetCompression(bool compression_) override
const char * GetDataFormatString() const
void SetTime(double t)
Set physical time (for time-dependent simulations)
const std::string & GetPrefixPath() const
Get the path where the DataCollection will be saved.
QFieldMap::const_iterator QFieldMapConstIterator
virtual void Load(int cycle_=0)
Load the collection based on its VisIt data (described in its root file)
void Deregister(const std::string &fname, bool own_data)
Unregister association between field field and name fname.
void SetOwnData(bool o)
Set the ownership of collection data.
virtual ~VisItDataCollection()
We will delete the mesh and fields if we own them.
void SetPadDigitsRank(int digits)
Set the number of digits used for the MPI rank in filenames.
std::map< std::string, VisItFieldInfo > field_info_map
ParGridFunction * GetParField(const std::string &field_name)
Get a pointer to a parallel grid function in the collection.
void UseRestartMode(bool restart_mode_)
int NumFields() const
Returns the number of registered fields.
QFieldMap::iterator QFieldMapIterator
int precision
Precision (number of digits) used for the text output of doubles.
OutStream err(std::cerr)
Global stream used by the library for standard error output. Initially it uses the same std::streambu...
Definition: globals.hpp:71
bool IsBinaryFormat() const
Returns true if the output format is BINARY or BINARY32, false if ASCII.
int myid
MPI rank (in parallel)
virtual void SaveField(const std::string &field_name)
Save one field, assuming the collection directory already exists.
static const int precision_default
Default value for precision.
const_iterator begin() const
Returns a begin const iterator to the registered fields.
GFieldMap::MapType FieldMapType
void SetMaxLevelsOfDetail(int max_levels_of_detail)
Set VisIt parameter: maximum levels of detail for the MultiresControl.
virtual void RegisterField(const std::string &field_name, GridFunction *gf)
Add a grid function to the collection and update the root file.
virtual ~DataCollection()
Delete the mesh and fields if owned by the collection.
void SetTimeStep(double ts)
Set the simulation time step (for time-dependent simulations)
std::string GetFieldFileName(const std::string &field_name) const
GFieldMap::const_iterator FieldMapConstIterator
void clear()
Clears the map of registered fields without reclaiming memory.
void SaveRootFile()
Save a VisIt root file for the collection.
virtual void DeregisterQField(const std::string &field_name)
Remove a QuadratureFunction from the collection.
const_iterator end() const
Returns an end const iterator to the registered fields.
ParaViewDataCollection(const std::string &collection_name, mfem::Mesh *mesh_=NULL)
Constructor. The collection name is used when saving the data.
std::string prefix_path
A path where the directory with results is saved. If not empty, it has &#39;/&#39; at the end...
Lightweight adaptor over an std::map from strings to pointer to T.
virtual void RegisterQField(const std::string &q_field_name, QuadratureFunction *qf)
Add a quadrature function to the collection and update the root file.
virtual void DeregisterField(const std::string &field_name)
Remove a grid function from the collection.
virtual void Load(int cycle_=0)
Load the collection. Not implemented in the base class DataCollection.
RefCoord t[3]
void SetLevelsOfDetail(int levels_of_detail_)
Mesh * GetMesh()
Get a pointer to the mesh in the collection.
std::string GetMeshFileName() const
void LoadVisItRootFile(const std::string &root_name)
T * Get(const std::string &fname) const
Get a pointer to the field associated with name fname.
virtual void SetMesh(Mesh *new_mesh)
Set/change the mesh associated with the collection.
virtual void Save() override
std::string name
Name of the collection, used as a directory name when saving.
virtual void SetMesh(Mesh *new_mesh)
Set/change the mesh associated with the collection.
void DeleteAll()
Delete all data owned by VisItDataCollection including field data information.
void DeleteData(bool own_data)
Clear all associations between names and fields.
Class for parallel grid function.
Definition: pgridfunc.hpp:32
double time_step
Time step i.e. delta_t (for time-dependent simulations)
OutStream out(std::cout)
Global stream used by the library for standard output. Initially it uses the same std::streambuf as s...
Definition: globals.hpp:66
int num_procs
Number of MPI ranks (in parallel)
Class for parallel meshes.
Definition: pmesh.hpp:32
Class representing a function through its values (scalar or vector) at quadrature points...
Definition: gridfunc.hpp:734
std::map< std::string, T * > MapType
void SetPrefixPath(const std::string &prefix)
Set the path where the DataCollection will be saved.