MFEM  v3.4
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, Lawrence Livermore National Security, LLC. Produced at
2 // the Lawrence Livermore National Laboratory. LLNL-CODE-443211. All Rights
3 // reserved. See file COPYRIGHT for details.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability see http://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the GNU Lesser General Public License (as published by the Free
10 // Software Foundation) version 2.1 dated February 1999.
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 
23 namespace mfem
24 {
25 
26 /// Lightweight adaptor over an std::map from strings to pointer to T
27 template<typename T>
29 {
30 public:
31  typedef std::map<std::string, T*> MapType;
32  typedef typename MapType::iterator iterator;
33  typedef typename MapType::const_iterator const_iterator;
34 
35  /// Register field @a field with name @a fname
36  /** Replace existing field associated with @a fname (and optionally
37  delete associated pointer if @a own_data is true) */
38  void Register(const std::string& fname, T* field, bool own_data)
39  {
40  T*& ref = field_map[fname];
41  if (own_data)
42  {
43  delete ref; // if newly allocated -> ref is null -> OK
44  }
45  ref = field;
46  }
47 
48  /// Unregister association between field @a field and name @a fname
49  /** Optionally delete associated pointer if @a own_data is true */
50  void Deregister(const std::string& fname, bool own_data)
51  {
52  iterator it = field_map.find(fname);
53  if ( it != field_map.end() )
54  {
55  if (own_data)
56  {
57  delete it->second;
58  }
59  field_map.erase(it);
60  }
61  }
62 
63  /// Clear all associations between names and fields
64  /** Delete associated pointers when @a own_data is true */
65  void DeleteData(bool own_data)
66  {
67  for (iterator it = field_map.begin(); it != field_map.end(); ++it)
68  {
69  if (own_data)
70  {
71  delete it->second;
72  }
73  it->second = NULL;
74  }
75  }
76 
77  /// Predicate to check if a field is associated with name @a fname
78  bool Has(const std::string& fname) const
79  {
80  return field_map.find(fname) != field_map.end();
81  }
82 
83  /// Get a pointer to the field associated with name @a fname
84  /** @return Pointer to field associated with @a fname or NULL */
85  T* Get(const std::string& fname) const
86  {
87  const_iterator it = field_map.find(fname);
88  return it != field_map.end() ? it->second : NULL;
89  }
90 
91  /// Returns a const reference to the underlying map
92  const MapType& GetMap() const { return field_map; }
93 
94  /// Returns the number of registered fields
95  int NumFields() const { return field_map.size(); }
96 
97  /// Returns a begin iterator to the registered fields
98  iterator begin() { return field_map.begin(); }
99  /// Returns a begin const iterator to the registered fields
100  const_iterator begin() const { return field_map.begin(); }
101 
102  /// Returns an end iterator to the registered fields
103  iterator end() { return field_map.end(); }
104  /// Returns an end const iterator to the registered fields
105  const_iterator end() const { return field_map.end(); }
106 
107  /// Returns an iterator to the field @a fname
108  iterator find(const std::string& fname)
109  { return field_map.find(fname); }
110 
111  /// Returns a const iterator to the field @a fname
112  const_iterator find(const std::string& fname) const
113  { return field_map.find(fname); }
114 
115  /// Clears the map of registered fields without reclaiming memory
116  void clear() { field_map.clear(); }
117 
118 protected:
120 };
121 
122 
123 /** A class for collecting finite element data that is part of the same
124  simulation. Currently, this class groups together grid functions (fields),
125  quadrature functions (q-fields), and the mesh that they are defined on. */
127 {
128 private:
129  /// A collection of named GridFunctions
131 
132  /// A collection of named QuadratureFunctions
134 public:
138 
142 
143  /// Format constants to be used with SetFormat().
144  /** Derived classes can define their own format enumerations and override the
145  method SetFormat() to perform input validation. */
146  enum Format
147  {
148  SERIAL_FORMAT = 0, /**<
149  MFEM's serial ascii format, using the methods Mesh::Print() /
150  ParMesh::Print(), and GridFunction::Save() / ParGridFunction::Save().*/
152  MFEM's parallel ascii format, using the methods ParMesh::ParPrint() and
153  GridFunction::Save() / ParGridFunction::Save(). */
154  };
155 
156 protected:
157  /// Name of the collection, used as a directory name when saving
158  std::string name;
159 
160  /** @brief A path where the directory with results is saved.
161  If not empty, it has '/' at the end. */
162  std::string prefix_path;
163 
164  /** A FieldMap mapping registered field names to GridFunction pointers. */
166 
167  /** A FieldMap mapping registered names to QuadratureFunction pointers. */
169 
170  /// The (common) mesh for the collected fields
172 
173  /// Time cycle; for time-dependent simulations cycle >= 0, otherwise = -1.
174  /** When cycle >= 0, it is appended to directory names. */
175  int cycle;
176  /// Physical time (for time-dependent simulations)
177  double time;
178 
179  /// Time step i.e. delta_t (for time-dependent simulations)
180  double time_step;
181 
182  /// Serial or parallel run? False iff mesh is a ParMesh
183  bool serial;
184  /// Append rank to any output file names.
186 
187  /// MPI rank (in parallel)
188  int myid;
189  /// Number of MPI ranks (in parallel)
191 #ifdef MFEM_USE_MPI
192  /// Associated MPI communicator
193  MPI_Comm m_comm;
194 #endif
195 
196  /// Precision (number of digits) used for the text output of doubles
198  /// Number of digits used for the cycle and MPI rank in filenames
200 
201  /// Default value for precision
202  static const int precision_default = 6;
203  /// Default value for pad_digits_*
204  static const int pad_digits_default = 6;
205 
206  /// Output mesh format: see the #Format enumeration
207  int format;
208 
209  /// Should the collection delete its mesh and fields
210  bool own_data;
211 
212  /// Error state
213  int error;
214 
215  /// Delete data owned by the DataCollection keeping field information
216  void DeleteData();
217  /// Delete data owned by the DataCollection including field information
218  void DeleteAll();
219 
220  std::string GetMeshShortFileName() const;
221  std::string GetMeshFileName() const;
222  std::string GetFieldFileName(const std::string &field_name) const;
223 
224  /// Save one field to disk, assuming the collection directory exists
225  void SaveOneField(const FieldMapIterator &it);
226 
227  /// Save one q-field to disk, assuming the collection directory exists
228  void SaveOneQField(const QFieldMapIterator &it);
229 
230  // Helper method
231  static int create_directory(const std::string &dir_name,
232  const Mesh *mesh, int myid);
233 
234 public:
235  /// Initialize the collection with its name and Mesh.
236  /** When @a mesh_ is NULL, then the real mesh can be set with SetMesh(). */
237  explicit DataCollection(const std::string& collection_name,
238  Mesh *mesh_ = NULL);
239 
240  /// Add a grid function to the collection
241  virtual void RegisterField(const std::string& field_name, GridFunction *gf)
242  { field_map.Register(field_name, gf, own_data); }
243 
244  /// Remove a grid function from the collection
245  virtual void DeregisterField(const std::string& field_name)
246  { field_map.Deregister(field_name, own_data); }
247 
248  /// Add a QuadratureFunction to the collection.
249  virtual void RegisterQField(const std::string& q_field_name,
250  QuadratureFunction *qf)
251  { q_field_map.Register(q_field_name, qf, own_data); }
252 
253 
254  /// Remove a QuadratureFunction from the collection
255  virtual void DeregisterQField(const std::string& field_name)
256  { q_field_map.Deregister(field_name, own_data); }
257 
258  /// Check if a grid function is part of the collection
259  bool HasField(const std::string& name) const
260  { return field_map.Has(name); }
261 
262  /// Get a pointer to a grid function in the collection.
263  /** Returns NULL if @a field_name is not in the collection. */
264  GridFunction *GetField(const std::string& field_name)
265  { return field_map.Get(field_name); }
266 
267 #ifdef MFEM_USE_MPI
268  /// Return the associated MPI communicator or MPI_COMM_NULL.
269  MPI_Comm GetComm() const { return m_comm; }
270 
271  /// Get a pointer to a parallel grid function in the collection.
272  /** Returns NULL if @a field_name is not in the collection.
273  @note The GridFunction pointer stored in the collection is statically
274  cast to ParGridFunction pointer. */
275  ParGridFunction *GetParField(const std::string& field_name)
276  { return static_cast<ParGridFunction*>(GetField(field_name)); }
277 #endif
278 
279  /// Check if a QuadratureFunction with the given name is in the collection.
280  bool HasQField(const std::string& q_field_name) const
281  { return q_field_map.Has(q_field_name); }
282 
283  /// Get a pointer to a QuadratureFunction in the collection.
284  /** Returns NULL if @a field_name is not in the collection. */
285  QuadratureFunction *GetQField(const std::string& q_field_name)
286  { return q_field_map.Get(q_field_name); }
287 
288  /// Get a const reference to the internal field map.
289  /** The keys in the map are the field names and the values are pointers to
290  GridFunction%s. */
291  const FieldMapType &GetFieldMap() const
292  { return field_map.GetMap(); }
293 
294  /// Get a const reference to the internal q-field map.
295  /** The keys in the map are the q-field names and the values are pointers to
296  QuadratureFunction%s. */
298  { return q_field_map.GetMap(); }
299 
300  /// Get a pointer to the mesh in the collection
301  Mesh *GetMesh() { return mesh; }
302  /// Set/change the mesh associated with the collection
303  /** When passed a Mesh, assumes the serial case: MPI rank id is set to 0 and
304  MPI num_procs is set to 1. When passed a ParMesh, MPI info from the
305  ParMesh is used to set the DataCollection's MPI rank and num_procs. */
306  virtual void SetMesh(Mesh *new_mesh);
307 #ifdef MFEM_USE_MPI
308  /// Set/change the mesh associated with the collection.
309  /** For this case, @a comm is used to set the DataCollection's MPI rank id
310  and MPI num_procs, which influences the how files are saved for domain
311  decomposed meshes. */
312  virtual void SetMesh(MPI_Comm comm, Mesh *new_mesh);
313 #endif
314 
315  /// Set time cycle (for time-dependent simulations)
316  void SetCycle(int c) { cycle = c; }
317  /// Set physical time (for time-dependent simulations)
318  void SetTime(double t) { time = t; }
319 
320  /// Set the simulation time step (for time-dependent simulations)
321  void SetTimeStep(double ts) { time_step = ts; }
322 
323  /// Get time cycle (for time-dependent simulations)
324  int GetCycle() const { return cycle; }
325  /// Get physical time (for time-dependent simulations)
326  double GetTime() const { return time; }
327  /// Get the simulation time step (for time-dependent simulations)
328  double GetTimeStep() const { return time_step; }
329 
330  /// Get the name of the collection
331  const std::string& GetCollectionName() const { return name; }
332  /// Set the ownership of collection data
333  void SetOwnData(bool o) { own_data = o; }
334 
335  /// Set the precision (number of digits) used for the text output of doubles
336  void SetPrecision(int prec) { precision = prec; }
337  /// Set the number of digits used for both the cycle and the MPI rank
338  void SetPadDigits(int digits) { pad_digits_cycle=pad_digits_rank = digits; }
339  /// Set the number of digits used for the cycle
340  void SetPadDigitsCycle(int digits) { pad_digits_cycle = digits; }
341  /// Set the number of digits used for the MPI rank in filenames
342  void SetPadDigitsRank(int digits) { pad_digits_rank = digits; }
343  /// Set the desired output mesh and data format.
344  /** See the enumeration #Format for valid options. Derived classes can define
345  their own format enumerations and override this method to perform input
346  validation. */
347  virtual void SetFormat(int fmt);
348 
349  /// Set the path where the DataCollection will be saved.
350  void SetPrefixPath(const std::string &prefix);
351 
352  /// Get the path where the DataCollection will be saved.
353  const std::string &GetPrefixPath() const { return prefix_path; }
354 
355  /// Save the collection to disk.
356  /** By default, everything is saved in the "prefix_path" directory with
357  subdirectory name "collection_name" or "collection_name_cycle" for
358  time-dependent simulations. */
359  virtual void Save();
360  /// Save the mesh, creating the collection directory.
361  virtual void SaveMesh();
362  /// Save one field, assuming the collection directory already exists.
363  virtual void SaveField(const std::string &field_name);
364  /// Save one q-field, assuming the collection directory already exists.
365  virtual void SaveQField(const std::string &q_field_name);
366 
367  /// Load the collection. Not implemented in the base class DataCollection.
368  virtual void Load(int cycle_ = 0);
369 
370  /// Delete the mesh and fields if owned by the collection
371  virtual ~DataCollection();
372 
373  /// Errors returned by Error()
374  enum { NO_ERROR = 0, READ_ERROR = 1, WRITE_ERROR = 2 };
375 
376  /// Get the current error state
377  int Error() const { return error; }
378  /// Reset the error state
379  void ResetError(int err = NO_ERROR) { error = err; }
380 };
381 
382 
383 /// Helper class for VisIt visualization data
385 {
386 public:
387  std::string association;
390  VisItFieldInfo(std::string _association, int _num_components)
391  { association = _association; num_components = _num_components; }
392 };
393 
394 /// Data collection with VisIt I/O routines
396 {
397 protected:
398  // Additional data needed in the VisIt root file, which describes the mesh
399  // and all the fields in the collection
402  std::map<std::string, VisItFieldInfo> field_info_map;
403  typedef std::map<std::string, VisItFieldInfo>::iterator FieldInfoMapIterator;
404 
405  /// Prepare the VisIt root file in JSON format for the current collection
406  std::string GetVisItRootString();
407  /// Read in a VisIt root file in JSON format
408  void ParseVisItRootString(const std::string& json);
409 
410  // Helper functions for Load()
411  void LoadVisItRootFile(const std::string& root_name);
412  void LoadMesh();
413  void LoadFields();
414 
415 public:
416  /// Constructor. The collection name is used when saving the data.
417  /** If @a mesh_ is NULL, then the mesh can be set later by calling either
418  SetMesh() or Load(). The latter works only in serial. */
419  VisItDataCollection(const std::string& collection_name, Mesh *mesh_ = NULL);
420 
421 #ifdef MFEM_USE_MPI
422  /// Construct a parallel VisItDataCollection to be loaded from files.
423  /** Before loading the collection with Load(), some parameters in the
424  collection can be adjusted, e.g. SetPadDigits(), SetPrefixPath(), etc. */
425  VisItDataCollection(MPI_Comm comm, const std::string& collection_name,
426  Mesh *mesh_ = NULL);
427 #endif
428 
429  /// Set/change the mesh associated with the collection
430  virtual void SetMesh(Mesh *new_mesh);
431 
432 #ifdef MFEM_USE_MPI
433  /// Set/change the mesh associated with the collection.
434  virtual void SetMesh(MPI_Comm comm, Mesh *new_mesh);
435 #endif
436 
437  /// Add a grid function to the collection and update the root file
438  virtual void RegisterField(const std::string& field_name, GridFunction *gf);
439 
440  /// Set VisIt parameter: maximum levels of detail for the MultiresControl
441  void SetMaxLevelsOfDetail(int max_levels_of_detail);
442 
443  /** @brief Delete all data owned by VisItDataCollection including field data
444  information. */
445  void DeleteAll();
446 
447  /// Save the collection and a VisIt root file
448  virtual void Save();
449 
450  /// Save a VisIt root file for the collection
451  void SaveRootFile();
452 
453  /// Load the collection based on its VisIt data (described in its root file)
454  virtual void Load(int cycle_ = 0);
455 
456  /// We will delete the mesh and fields if we own them
457  virtual ~VisItDataCollection() {}
458 };
459 
460 }
461 
462 #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:27
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.
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.
bool HasQField(const std::string &q_field_name) const
Check if a QuadratureFunction with the given name is in the collection.
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.
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.
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.
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.
virtual void RegisterField(const std::string &field_name, GridFunction *gf)
Add a grid function to the collection.
VisItFieldInfo(std::string _association, int _num_components)
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.
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.
static const int pad_digits_default
Default value for pad_digits_*.
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.
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:69
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.
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 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.
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.
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)
int num_procs
Number of MPI ranks (in parallel)
Class representing a function through its values (scalar or vector) at quadrature points...
Definition: gridfunc.hpp:401
std::map< std::string, T * > MapType
void SetPrefixPath(const std::string &prefix)
Set the path where the DataCollection will be saved.