Chombo + EB + MF  3.2
CH_HDF5.H
Go to the documentation of this file.
1 #ifdef CH_LANG_CC
2 /*
3  * _______ __
4  * / ___/ / ___ __ _ / / ___
5  * / /__/ _ \/ _ \/ V \/ _ \/ _ \
6  * \___/_//_/\___/_/_/_/_.__/\___/
7  * Please refer to Copyright.txt, in Chombo's root directory.
8  */
9 #endif
10 
11 #ifndef _CH_HDF5_H_
12 #define _CH_HDF5_H_
13 
14 #define CHOFFSET(object, member) (int)((char*)&(object.member) - (char*)&object)
15 
16 #ifdef CH_USE_HDF5 // if you don't have CH_USE_HDF5, then this file is useless
17 
18 #include <iostream>
19 using std::cout;
20 using std::endl;
21 
22 #ifdef CH_MPI
23 #include "mpi.h"
24 #endif
25 
26 #include "LevelData.H"
27 #include "HDF5Portable.H"
28 // hdf5 #defines inline.... duh!
29 #undef inline
30 #include <string>
31 #include <map>
32 #include "RealVect.H"
33 #include "CH_Timer.H"
34 #include "LoadBalance.H"
35 #include "LayoutIterator.H"
36 #include "Vector.H"
37 #include "memtrack.H"
38 #include "FluxBox.H"
39 
40 #ifdef CH_MULTIDIM
41 #include "BoxTools_ExternC_Mangler.H" // Generated by lib/utils/multidim/mangle_externs.sh
42 #endif
43 #include "NamespaceHeader.H"
44 #ifdef H5_USE_16_API
45 #define H516
46 #endif
47 
48 template <class T>
49 void read(T& item, Vector<Vector<char> >& a_allocatedBuffers, const Box& box,
50  const Interval& comps)
51 {
52  Vector<void*> b(a_allocatedBuffers.size());
53  for (int i=0; i<b.size(); ++i) b[i]= &(a_allocatedBuffers[i][0]);
54  read(item, b, box, comps);
55 }
56 
57 namespace CH_HDF5
58 {
59  enum IOPolicy
60  {
62  IOPolicyMultiDimHyperslab = (1<<0), // HDF5 will internally linearize the
63  // data
64  IOPolicyCollectiveWrite = (1<<1) // HDF5 will collectively write data
65  };
66 }
67 
68 using std::map;
69 
70 class HDF5Handle;
71 
72 // CH_HDF5.H
73 // ============
74 
75 /// user-friendly function to write out data on a AMR level
76 /**
77  all data is IN data.
78  a_handle: handle open and ready.
79  a_data: data, what else do you want to know ?
80  a_dx: the grid spacing at this level
81  a_dt: the timestep size that was last completed
82  a_time: the time of this level (might not be the same as other levels)
83  a_domain: the problem domain, represented at this level of refinement
84  a_refRatio:the refinement of a_level+1 wrt a_level. for vis systems it
85  would probably help if you use 1 for a_level==max_level.
86 */
87 template <class T>
88 int writeLevel(HDF5Handle& a_handle,
89  const int& a_level,
90  const T& a_data,
91  const Real& a_dx,
92  const Real& a_dt,
93  const Real& a_time,
94  const Box& a_domain,
95  const int& a_refRatio,
96  const IntVect& outputGhost = IntVect::Zero,
97  const Interval& comps = Interval());
98 
99 template <class T>
100 int readLevel(HDF5Handle& a_handle,
101  const int& a_level,
102  LevelData<T>& a_data,
103  Real& a_dx,
104  Real& a_dt,
105  Real& a_time,
106  Box& a_domain,
107  int& a_refRatio,
108  const Interval& a_comps = Interval(),
109  bool setGhost = false);
110 
111 template <class T>
112 int readLevel(HDF5Handle& a_handle,
113  const int& a_level,
114  LevelData<T>& a_data,
115  RealVect& a_dx,
116  Real& a_dt,
117  Real& a_time,
118  Box& a_domain,
119  IntVect& a_refRatio,
120  const Interval& a_comps = Interval(),
121  bool setGhost = false);
122 
123 // More basic HDF5 functions. These can be used at a persons own
124 // discretion. Refer to the User's Guide for a description of how these
125 // functions interact with the convenience functions writeLevel, etc.
126 
127 /// writes BoxLayout to HDF5 file.
128 /**
129  Writes BoxLayout to HDF5 file. Only one BoxLayout per group is permitted, this operation overwrites
130  previous entries.
131  This operation assumes boxes are cell-centered.\\
132  returns: success: 0\\
133  HDF5 error: negative error code.\\
134 */
135 int write(HDF5Handle& a_handle,
136  const BoxLayout& a_layout,
137  const std::string& name = "boxes");
138 
139 /// writes a BoxLayoutData<T> to an HDF5 file.
140 /**
141  writes a BoxLayoutData<T> to an HDF5 file.\\
142  returns: success: 0\\
143  HDF5 error: negative error code.\\
144 */
145 template <class T>
146 int write(HDF5Handle& a_handle,
147  const BoxLayoutData<T>& a_data,
148  const std::string& a_name,
149  IntVect outputGhost = IntVect::Zero,
150  const Interval& comps = Interval(),
151  bool newForm = false);
152 
153 /// writes a LevelData<T> to an HDF5 file.
154 /**
155  Writes a LevelData<T> to an HDF5 file.
156  the DisjointBoxLayout is not written out with the data, the user is required to
157  handle that object seperately. (see the "read" interface below).\\
158  returns: success: 0\\
159  HDF5 error: negative error code.\\
160 */
161 template <class T>
162 int write(HDF5Handle& a_handle,
163  const LevelData<T>& a_data,
164  const std::string& a_name,
165  const IntVect& outputGhost = IntVect::Zero,
166  const Interval& comps = Interval());
167 
168 /// reads Vector<Box> from location specified by a_handle.
169 /**
170  Reads BoxLayout from the group specified by a_handle.
171  Only one BoxLayout per group is permitted, this operation overwrites.
172  This operation assumes boxes are cell-centered.\\
173  returns: success: 0\\
174  HDF5 error: negative error code.\\
175 
176  Arg name refers to the name of an HDF5 dataset in which the boxes you want
177  are stored in a particular format, which is the one used to output DataLayouts
178  to HDF5. Here is an example of that format (as shown by h5dump):
179 
180  DATASET "datalayout" {
181  DATATYPE H5T_COMPOUND {
182  H5T_STD_I32LE "lo_i";
183  H5T_STD_I32LE "lo_j";
184  H5T_STD_I32LE "hi_i";
185  H5T_STD_I32LE "hi_j";
186  }
187  DATASPACE SIMPLE { ( 2 ) / ( 2 ) }
188  DATA {
189  (0): {
190  0,
191  0,
192  1,
193  0
194  },
195  (1): {
196  2,
197  0,
198  2,
199  1
200  }
201  }
202  }
203 
204 */
205 int read(HDF5Handle& a_handle,
206  Vector<Box>& boxes,
207  const std::string& name = "boxes");
208 
209 /// reads the set of Boxes out from the level_* groups of a Chombo HDF5 AMR file
210 /**
211  goes to all groups named level_n for level n = 0 to numLevel and fills the
212  vector of boxes. CH_HDF5.Handle must be set to the root group (the default location when the
213  file is just opened.
214 
215  A handy *get everything* version of int read(HDF5Handle& a_handle, Vector<Box>& boxes)
216 */
217 int readBoxes(HDF5Handle& a_handle,
218  Vector<Vector<Box> >& boxes);
219 
220 /// FArrayBox-at-a-time read function. FArrayBox gets redefined in the function. Reads data field named by a_dataName.
221 /** FArrayBox gets redefined in the function.
222  it will be sized to the box at the [level,boxNumber] position and have components
223  the go from [0,a_components.size()-1]
224 
225  some meta-data has to be read again and again internally to do this function, but
226  it should make some out-of-core tools possible.
227 */
228 
229 int readFArrayBox(HDF5Handle& a_handle,
230  FArrayBox& a_fab,
231  int a_level,
232  int a_boxNumber,
233  const Interval& a_components,
234  const std::string& a_dataName = "data" );
235 
236 /// read BoxLayoutData named a_name from location specified by a_handle.
237 /**
238  Read BoxLayoutData named a_name from location specified by a_handle. User must supply the correct BoxLayout for this function if redefineData == true. \\
239  returns: success: 0\\
240  bad location: 1\\
241  HDF5 error: negative error code.\\
242 */
243 template <class T>
244 int read(HDF5Handle& a_handle,
245  BoxLayoutData<T>& a_data,
246  const std::string& a_name,
247  const BoxLayout& a_layout,
248  const Interval& a_comps = Interval(),
249  bool redefineData = true);
250 
251 /// read LevelData named a_name from location specified by a_handle.
252 /**
253 
254  Read LevelData named a_name from location specified by a_handle.
255  User must supply the correct BoxLayout for this function.
256 
257  Arg a_name is significant: the HDF5 group to which a_handle is set is
258  is assumed to contain a dataset called <a_name>:datatype=<some integer>,
259  a dataset called <a_name>:offsets=<some integer>, and a subgroup named
260  <a_name>_attributes. You will have all these items if you dumped your
261  LevelData out using the corresponding write() function defined here.
262 
263  If arg redefineData==false, then the user must pass in a valid LevelData.
264  Otherwise, this function figures out how many components and ghost cells there
265  are, and allocates the correct amount of space. The actual FArray data held
266  by the LevelData gets filled in here, regardless of redefineData; "redefine"
267  alludes to the family of define() functions.
268 
269  returns: success: 0\\
270  bad location: 1\\
271  HDF5 error: negative error code.\\
272 */
273 template <class T>
274 int read(HDF5Handle& a_handle,
275  LevelData<T>& a_data,
276  const std::string& a_name,
277  const DisjointBoxLayout& a_layout,
278  const Interval& a_comps = Interval(),
279  bool redefineData = true);
280 
281 /// Handle to a particular group in an HDF file.
282 /**
283  HDF5Handle is a handle to a particular group in an HDF file. Upon
284  construction, group is defined to be the root. All data is
285  written and read assuming the native representations for data on
286  the architecture it is running on. When a file is opened, these
287  settings are checked and an error is flagged when things don't
288  match up. It is the USER'S responsibility to close() this object
289  when it is no longer needed.
290 
291 */
293 {
294 public:
295  ///
296  /**
297  Enumeration of opening modes for HDF5 files. \\
298 
299  CREATE: file is created if it didn't exist, or an existing file of
300  the same name is clobbered.\\
301 
302  CREATE_SERIAL: in serial execution this is equivalent to CREATE. In parallel this opens
303  a file just on this particular calling processor. Working with LevelData and this kind of
304  file is difficult to get right and most users will not use this mode. \\
305 
306  OPEN_RDONLY: existing file is opened in read-only mode. If the
307  file doesn't already exist then open fails and isOpen() returns
308  false.\\
309 
310  OPEN_RDWR: existing file is opened in read-write mode. If the file
311  doesn't already exist then open fails and isOpen() returns false.\\
312 
313  */
314  enum mode
315  {
320  };
321 
322  /// {\bf constructor}
323 
324  ///
325  /**
326  Default constructor. User must call open() prior to using
327  constructed object.
328  */
329  HDF5Handle();
330 
331  ///
332  /** Opens file and sets the current group to the root "/" group. \\
333 
334  if mode == CREATE, then file is created if it didn't exist, or an
335  existing file of the same name is clobbered.\\
336 
337  if mode == OPEN_*, then existing file is opened, if the file doesn't
338  already exist then open fails and isOpen() returns false.\\
339 
340  Writes basic information such as SpaceDim and testReal to a global
341  group. The global group is named Chombo_global in the call signature
342  where no globalGroupName is passed.\\
343  */
344  HDF5Handle(
345  const std::string& a_filename,
346  mode a_mode,
347  const char *a_globalGroupName="Chombo_global");
348 
349  // HDF5Handle(const std::string& a_filename, mode a_mode);
350 
351  ~HDF5Handle();
352 
353  /// {\bf File functions}
354 
355  ///
356  /**
357  Opens file and sets the current group of this HDF5Handle to the
358  root "/" group. File that this HDF5Handle previously pointed at is
359  NOT closed, that is the users responsibility.\\
360 
361  if mode == CREATE, then file is created if it didn't exist, or an
362  existing file of the same name is clobbered.\\
363 
364  if mode == OPEN_*, then existing file is opened, if the file doesn't
365  already exist then open fails and isOpen() returns false.\\
366 
367  Writes basic information such as SpaceDim and testReal to a global
368  group. The global group is named Chombo_global in the call signature
369  where no globalGroupName is passed.\\
370 
371  returns:\\
372  0 on success\\
373  negative number if file open failed (return code from HDF5)\\
374  1 file does not appear to contain datacheck info, probably not a Chombo file\\
375  2 on data bit size differences between code and file.\\
376 
377  aborts on SpaceDim not matching between code and file\\
378 
379  */
380  int open(
381  const std::string& a_filename,
382  mode a_mode,
383  const char *a_globalGroupName="Chombo_global");
384 
385  // int open(const std::string& a_filename, mode a_mode);
386 
387  ///
388  /**
389  A NULL or failed constructed HDF5Handle will return false.
390  */
391  bool isOpen() const;
392 
393  ///
394  /**
395  Closes the file. Must be called to close file. Files are not
396  automatically closed.
397 
398  */
399  void close();
400 
401  /// {\bf Group functions}
402 
403  ///
404  /**
405  Sets the current group to be "/level_x" where x=a_level.
406  */
407  void setGroupToLevel(int a_level);
408 
409  ///
410  /**
411  Set group to users choice, referenced from file root.
412  groupAbsPath will look like a Unix file path:
413  "/mySpecialData/group1/" "/" is the root group of the
414  file. returns a negative value on failure
415 
416  */
417  int setGroup(const std::string& groupAbsPath);
418 
419  ///
420  /**
421  Add the indicated string to the group path. For example, if
422  getGroup() returns "/foo" then after pushGroup("bar"), getGroup()
423  will return "/foo/bar".
424  Return value is whatever the internal setGroup returns.
425  */
426  int pushGroup( const std::string& grp );
427 
428  ///
429  /**
430  Pop off the last element of the group path. For example, "/foo/bar"
431  becomes "/foo". It's an error to call this function if the group,
432  going in, is "/".
433  Return value is whatever the internal setGroup returns, when we call
434  it to reset the group's path.
435  */
436  int popGroup();
437 
438  ///
439  /**
440  Returns name of current group. groupAbsPath will look like a Unix
441  file path: "/mySpecialData/group1/" "/" is the root group of the
442  file.
443 
444  */
445  const std::string& getGroup() const;
446 
447  HDF5Handle::mode openMode() const {return m_mode;}
448  const hid_t& fileID() const;
449  const hid_t& groupID() const;
450  static hid_t box_id;
451  static hid_t intvect_id;
452  static hid_t realvect_id;
453  static map<std::string, std::string> groups;
454 
455 private:
456 
457  HDF5Handle(const HDF5Handle&);
459 
460  hid_t m_fileID;
463  bool m_isOpen;
464  std::string m_filename; // keep around for debugging
465  std::string m_group;
466  int m_level;
467 
468  // static hid_t file_access;
469  static bool initialized;
470  static void initialize();
471 
472 };
473 
474 /// data to be added to HDF5 files.
475 /**
476  HDF5HeaderData is a wrapper for some data maps to be added to HDF5
477  files. instead of an overdose of access functions, the maps are
478  made public and they can be manipulated by the user at will. They
479  maintain type safety.
480 
481  to add a Real data entry, a user can simply program as follows:
482 
483  <PRE>
484  Real dx;
485  .
486  .
487  HDF5HeaderData metaData;
488  metaData.m_real["dx"] = dx;
489  </PRE>
490 
491  If "dx" already existed, then it is overwritten, otherwise an entry is
492  created and added with the new value;
493 
494  To search for entries, the user does the following:
495  <PRE>
496 
497  HDF5HeaderData metaData;
498  HDF5Handle currentStep(filename);
499  currentStep.setGroupToLevel(0);
500  metaData.readFromFile(currentStep);
501  if (metaData.m_intvect.find("ghost") != metaData.m_intvect.end())
502  ghost = metaData.m_intvect["ghost"];
503  else
504  ghost = defaultGhostIntVect;
505 
506  </PRE>
507 
508  A user can skip the check for existence if they have reason to "know" the
509  data will be there. It is just good coding practice.
510 
511  To erase an entry, you can use:
512  <PRE>
513  metaData.m_real.erase("dx");
514  </PRE>
515 
516 */
518 {
519 public:
520 
521  ///
522  /**
523  Writes this HDF5HeaderData's current attribute list to the
524  current group in 'file.' Returns 0 on success, returns the
525  error code from HDF5 on failure.
526 
527  */
528  int writeToFile(HDF5Handle& file) const;
529 
530  ///
531  /**
532  Reads into this HDF5HeaderData's attribute list from file. Read
533  process is add/change, does not remove key-value pairs. Reads
534  from current group. Returns 0 on success, positive number if a
535  particular member of group caused an error, negative on general
536  error.
537 
538  */
539  int readFromFile(HDF5Handle& file);
540 
541  ///
542  void clear();
543 
544  ///
545  map<std::string, Real> m_real;
546 
547  ///
548  map<std::string, int> m_int;
549 
550  ///
551  map<std::string, std::string> m_string;
552 
553  ///
554  map<std::string, IntVect> m_intvect;
555 
556  ///
557  map<std::string, Box> m_box;
558 
559  ///
560  map<std::string, RealVect> m_realvect;
561 
562  //users should not need these functions in general
563 
564  int writeToLocation(hid_t loc_id) const;
565  int readFromLocation(hid_t loc_id);
566 
567  /// useful for debugging. dumps contents to std::cout
568  void dump() const;
569 
570 private:
571  static herr_t attributeScan(hid_t loc_id, const char *name, void *opdata);
572 };
573 
574 extern "C"
575 {
576 #ifdef H516
577  herr_t HDF5HeaderDataattributeScan(hid_t loc_id, const char *name, void *opdata);
578 #else
579  herr_t HDF5HeaderDataattributeScan(hid_t loc_id, const char *name, const H5A_info_t* info, void *opdata);
580 #endif
581 }
582 
583 std::ostream& operator<<(std::ostream& os, const HDF5HeaderData& data);
584 
585 /// Methods for writing multiple LevelData to an HDF5 file.
586 /**
587  While the write functions can be used to write a single LevelData, they
588  do not support writing intervals from multiple LevelData. This is often
589  the case when you want to add diagnostic information to the output file.
590  The purpose of this class is to handle writing multiple LevelData to
591  a single file. The dataset for the LevelData is created during
592  construction, when the boxes are written. Any number of intervals from
593  LevelData can then be written.
594 
595  Special behavior can be selected with the policy flags:
596  <ul>
597  <li> CH_HDF5::IOPolicyMultiDimHyperslab - The memory dataspace will be
598  set up so that HDF5 linearizes the data. Using this requires
599  T.dataPtr() and contiguous memory (so things like T=FluxBox will
600  not work).
601  <li> CH_HDF5::IOPolicyCollectiveWrite - The write for each dataset will
602  be collective. The write is always parallel but otherwise will
603  be done independently. May or may not provide a speedup but it
604  should work. Rumor has it the BG requires this.
605  </ul>
606 
607  With MPI-everywhere on x86 architectures, it is highly likely that
608  the default policy CH_HDF5::IOPolicyDefault gives the best performance
609  (i.e., use internal linearization and independent I/O). In some limited
610  testing (2 procs, 20 cores each):
611  IOPolicyDefault - 6.77
612  IOPolicyCollectiveWrite - 33.33
613  IOPolicyMultiDimHyperslab - 168.14
614  IOPolicyCollectiveWrite | IOPolicyMultiDimHyperslab - 30.94
615  which is difficult to make sense of.
616 
617  \note
618  <ul>
619  <li> Users are resposible for matching the total number of components
620  used to allocate the dataset match the number of levelData
621  added to the dataset
622  <li> Only a single datatype is supported per file
623  <li> Since the dataset must be created in advance, T can only be
624  statically allocatable (T::preAllocatable() == 0)
625  <li> API 1.6 is not supported
626  </ul>
627 */
628 template <class T>
630 {
631 
632 public:
633 
634  /// Constructor writes boxes and allocates dataset for LevelData
635  WriteMultiData(HDF5Handle& a_handle,
636  const BoxLayout& a_layout,
637  const int a_numIntv,
638  const string& a_name,
639  const int a_policyFlags = CH_HDF5::IOPolicyDefault,
640  const IntVect& a_outputGhost = IntVect::Zero,
641  const bool a_newForm = false);
642 
643  /// Destructor
645  {
646  H5Dclose(m_dataSet);
647  }
648 
649 private:
650 
651  // Copy and assignment not allowed
654 
655 public:
656 
657  /// Write an interval of LevelData to the dataset
658  int writeData(const BoxLayoutData<T>& a_data,
659  const Interval& a_intvMem,
660  const Interval& a_intvFile);
661 
662 protected:
663 
664  const int m_policyFlags; ///< Policies
665  const IntVect m_outputGhost; ///< Number of ghost cells written
666  const bool m_newForm; ///< ?
667 
668  char m_dataname[128]; ///< Name for level data dataset
669  hid_t m_dataSet; ///< Dataset for level data
670  Interval m_allIntvFile; ///< Interval for components in file
671  long m_maxBoxPerProc; ///< Maximum boxes written by any proc
672 
673  // Both of these vectors are size 1 in the outermost dimension since
674  // only 1 type is allowed.
676  ///< Offset in file for each process to
677  ///< write to
678  Vector<hid_t> m_types; ///< Type of data written
679 };
680 
681 //=============================================================================
682 //
683 // end of declarations.
684 //
685 //=============================================================================
686 
687 #if ( H5_VERS_MAJOR == 1 && H5_VERS_MINOR > 6 )
688 typedef hsize_t ch_offset_t;
689 #else
690 #if ( H5_VERS_MAJOR == 1 && H5_VERS_MINOR == 6 && H5_VERS_RELEASE >= 4 )
691 typedef hsize_t ch_offset_t;
692 #else
693 typedef hssize_t ch_offset_t;
694 #endif
695 #endif
696 
697 template<class T>
698 hid_t H5Type(const T* dummy);
699 
700 template< >
701 hid_t H5Type(const int* dummy);
702 
703 template< >
704 hid_t H5Type(const long long* dummy);
705 
706 template< >
707 hid_t H5Type(const float* dummy);
708 
709 template< >
710 hid_t H5Type(const double* dummy);
711 
712 template< >
713 hid_t H5Type(const Box* dummy);
714 
715 template< >
716 hid_t H5Type(const RealVect* dummy);
717 
718 template< >
719 hid_t H5Type(const IntVect* dummy);
720 
721 template<class T>
722 hid_t H5Type(const T* dummy)
723 {
724  // no such definition;
725  MayDay::Error(" H5Type(const T* dummy)");
726  return -4;
727 }
728 
729 void createData(hid_t& a_dataset,
730  hid_t& a_dataspace,
731  HDF5Handle& handle,
732  const std::string& name,
733  hid_t type,
734  hsize_t size);
735 
736 template <class T>
737 void createDataset(hid_t& a_dataset,
738  hid_t& a_dataspace,
739  HDF5Handle& handle,
740  const std::string& name,
741  const T* dummy,
742  hsize_t size)
743 {
744  createData(a_dataset, a_dataspace, handle, name, H5Type(dummy), size);
745 
746 }
747 
748 void writeDataset(hid_t a_dataset,
749  hid_t a_dataspace,
750  const void* start,
751  ch_offset_t off,
752  hsize_t count);
753 
754 void readDataset(hid_t a_dataset,
755  hid_t a_dataspace,
756  void* start,
757  ch_offset_t off,
758  hsize_t count);
759 
760 // non-user code used in implementation of communication
761 
763 {
766  void operator=(const OffsetBuffer& rhs);
767 };
768 
769 ostream& operator<<(ostream& os, const OffsetBuffer& ob);
770 
771 #include "NamespaceFooter.H"
772 
773 #include "BaseNamespaceHeader.H"
774 
775 #include "NamespaceVar.H"
776 
777 //OffsetBuffer specialization of linearSize
778 template < >
779 int linearSize(const CH_XDIR::OffsetBuffer& a_input);
780 
781 //OffsetBuffer specialization of linearIn
782 template < >
783 void linearIn(CH_XDIR::OffsetBuffer& a_outputT, const void* const a_inBuf);
784 
785 //OffsetBuffer specialization of linearOut
786 template < >
787 void linearOut(void* const a_outBuf, const CH_XDIR::OffsetBuffer& a_inputT);
788 
789 template < > int linearSize(const Vector<CH_XDIR::OffsetBuffer>& a_input);
790 template < > void linearIn(Vector<CH_XDIR::OffsetBuffer>& a_outputT, const void* const inBuf);
791 template < > void linearOut(void* const a_outBuf, const Vector<CH_XDIR::OffsetBuffer>& a_inputT);
792 
793 #include "BaseNamespaceFooter.H"
794 #include "NamespaceHeader.H"
795 
796 // First, template specializations for read/write for FArrayBox.
797 
798 template <>
799 inline void dataTypes(Vector<hid_t>& a_types, const BaseFab<int>& dummy)
800 {
801  a_types.resize(1);
802  a_types[0] = H5T_NATIVE_INT;
803 }
804 
805 template <>
806 inline void dataTypes(Vector<hid_t>& a_types, const BaseFab<char>& dummy)
807 {
808  a_types.resize(1);
809  a_types[0] = H5T_NATIVE_CHAR;
810 }
811 
812 /* since many compilers choke on the proper template specialization syntax
813  I am forced to pass in a dummy specialization argument
814 */
815 template <>
816 inline void dataTypes(Vector<hid_t>& a_types, const FArrayBox& dummy)
817 {
818  a_types.resize(1);
819  a_types[0] = H5T_NATIVE_REAL;
820 }
821 
822 template <>
823 inline void dataTypes(Vector<hid_t>& a_types, const FluxBox& dummy)
824 {
825  a_types.resize(1);
826  a_types[0] = H5T_NATIVE_REAL;
827 }
828 
829 template <>
830 inline void dataSize(const BaseFab<int>& item, Vector<int>& a_sizes,
831  const Box& box, const Interval& comps)
832 {
833  a_sizes[0] = box.numPts() * comps.size();
834 }
835 
836 template <>
837 inline void dataSize(const BaseFab<char>& item, Vector<int>& a_sizes,
838  const Box& box, const Interval& comps)
839 {
840  a_sizes[0] = box.numPts() * comps.size();
841 }
842 
843 template <>
844 inline void dataSize(const FArrayBox& item, Vector<int>& a_sizes,
845  const Box& box, const Interval& comps)
846 {
847  a_sizes[0] = box.numPts() * comps.size();
848 }
849 
850 template <>
851 inline void dataSize(const FluxBox& item, Vector<int>& a_sizes,
852  const Box& box, const Interval& comps)
853 {
854  int& t = a_sizes[0];
855  t = 0;
856  for (int dir =0; dir<CH_SPACEDIM; dir++)
857  {
858  Box edgeBox(surroundingNodes(box,dir));
859  t += edgeBox.numPts()*comps.size();
860  }
861 }
862 
863 template <>
864 inline const char* name(const FArrayBox& a_dummySpecializationArg)
865 {
866  // Attempt to get rid of warnings on IBM...
867  //static const char* name = "FArrayBox";
868  const char* name = "FArrayBox";
869  return name;
870 }
871 
872 template <>
873 inline const char* name(const BaseFab<int>& a_dummySpecializationArg)
874 {
875  // Attempt to get rid of warnings on IBM...
876  //static const char* name = "BaseFab<int>";
877  const char* name = "BaseFab<int>";
878  return name;
879 }
880 
881 template <>
882 inline const char* name(const BaseFab<char>& a_dummySpecializationArg)
883 {
884  // Attempt to get rid of warnings on IBM...
885  //static const char* name = "BaseFab<int>";
886  const char* name = "BaseFab<char>";
887  return name;
888 }
889 //
890 // now, generic, non-binary portable version of template functions
891 // for people who just want ot rely on linearIn/linearOut
892 
893 template <class T>
894 inline void dataTypes(Vector<hid_t>& a_types, const T& dummy)
895 {
896  a_types.resize(1);
897  a_types[0] = H5T_NATIVE_CHAR;
898 }
899 
900 template <class T>
901 inline void dataSize(const T& item, Vector<int>& a_sizes,
902  const Box& box, const Interval& comps)
903 {
904  a_sizes[0] = item.size(box, comps);
905 }
906 
907 template <class T>
908 inline void write(const T& item, Vector<void*>& a_allocatedBuffers,
909  const Box& box, const Interval& comps)
910 {
911  item.linearOut(a_allocatedBuffers[0], box, comps);
912 }
913 
914 template <class T>
915 inline void read(T& item, Vector<void*>& a_allocatedBuffers,
916  const Box& box, const Interval& comps)
917 {
918  item.linearIn(a_allocatedBuffers[0], box, comps);
919 }
920 
921 template <class T>
922 inline const char* name(const T& a_dummySpecializationArg)
923 {
924  static const char* name = "unknown";
925  return name;
926 }
927 
928 template <class T>
929 void getOffsets(Vector<Vector<long long> >& offsets, const BoxLayoutData<T>& a_data,
930  int types, const Interval& comps, const IntVect& outputGhost)
931 {
932  CH_TIME("getOffsets");
933  const BoxLayout& layout = a_data.boxLayout();
934  {
935  CH_TIMELEAF("offsets.resize");
936  offsets.resize(types, Vector<long long>(layout.size()+1));
937  // offsets.resize(layout.size() + 1, Vector<long long>(types));
938  }
939  for (int t=0; t<types; t++) offsets[t][0] = 0;
940  Vector<int> thisSize(types);
941  if (T::preAllocatable()==0)
942  { // static preAllocatable
943  T dummy;
944  unsigned int index = 1;
945  for (LayoutIterator it(layout.layoutIterator()); it.ok(); ++it)
946  {
947  Box region = layout[it()];
948  region.grow(outputGhost);
949  {
950  CH_TIMELEAF("dataSize");
951  dataSize(dummy, thisSize, region, comps);
952  }
953  for (unsigned int i=0; i<thisSize.size(); ++i)
954  {
955  CH_TIMELEAF("offsets calc");
956  //offsets[index][i] = offsets[index-1][i] + thisSize[i];
957  offsets[i][index] = offsets[i][index-1] + thisSize[i];
958  }
959  ++index;
960  }
961  }
962  else
963  { // symmetric and dynamic preallocatable need two pass I/O
964  OffsetBuffer buff;
965  //int index = 0;
966  for (DataIterator dit(a_data.dataIterator()); dit.ok(); ++dit)
967  {
968  int index = a_data.boxLayout().index(dit());
969  //int index = dit().intCode();
970  buff.index.push_back(index);
971  Box region = layout[dit()];
972  region.grow(outputGhost);
973  {
974  CH_TIMELEAF("dataSize");
975  dataSize(a_data[dit()], thisSize, region, comps);
976  }
977  buff.offsets.push_back(thisSize);
978  }
979  Vector<OffsetBuffer> gathering(numProc());
980  {
981  CH_TIMELEAF("gather");
982  gather(gathering, buff, uniqueProc(SerialTask::compute));
983  }
984  {
985  CH_TIMELEAF("broadcast");
987  }
988  // pout() << gathering<<endl;
989  for (int i=0; i<numProc(); ++i)
990  {
991  OffsetBuffer& offbuf = gathering[i];
992  for (int num=0; num<offbuf.index.size(); num++)
993  {
994  int index = offbuf.index[num];
995  for (unsigned int j=0; j<types; ++j)
996  {
997  CH_TIMELEAF("offsets calc");
998  //offsets[index+1][j] = offbuf.offsets[num][j];
999  offsets[j][index+1] = offbuf.offsets[num][j];
1000  }
1001  }
1002  }
1003  for (int i=0; i<layout.size(); i++)
1004  {
1005  for (unsigned int j=0; j<types; ++j)
1006  {
1007  CH_TIMELEAF("offsets calc");
1008  //offsets[i+1][j] += offsets[i][j];
1009  offsets[j][i+1] += offsets[j][i];
1010  }
1011  }
1012  }
1013 
1014  // pout() << offsets<<endl;
1015 }
1016 
1017 //==================================================================
1018 //
1019 // Statically pre-allocatable version only requires BoxLayout, not
1020 // the data
1021 //
1022 template <class T>
1024  const BoxLayout a_layout,
1025  int a_numTypes,
1026  const Interval& a_comps,
1027  const IntVect& a_outputGhost)
1028 {
1029  CH_TIME("getOffsets (prealloc)");
1030  CH_assert(T::preAllocatable() == 0);
1031 
1032  a_offsets.resize(a_numTypes, Vector<long long>(a_layout.size()+1));
1033  for (int t = 0; t != a_numTypes; ++t)
1034  {
1035  a_offsets[t][0] = 0;
1036  }
1037  Vector<int> thisSize(a_numTypes);
1038  T dummy;
1039  unsigned int index = 1;
1040  for (LayoutIterator it(a_layout.layoutIterator()); it.ok(); ++it)
1041  {
1042  Box region = a_layout[it()];
1043  region.grow(a_outputGhost);
1044  dataSize(dummy, thisSize, region, a_comps);
1045  for (int i = 0; i != thisSize.size(); ++i) // Loop over (1) types
1046  {
1047  //offsets[index][i] = offsets[index-1][i] + thisSize[i];
1048  a_offsets[i][index] = a_offsets[i][index-1] + thisSize[i];
1049  }
1050  ++index;
1051  }
1052 }
1053 
1054 //==================================================================
1055 //
1056 // Now, linear IO routines for a BoxLayoutData of T
1057 //
1058 template <class T>
1059 int write(HDF5Handle& a_handle, const BoxLayoutData<T>& a_data,
1060  const std::string& a_name, IntVect outputGhost,
1061  const Interval& in_comps, bool newForm)
1062 {
1063  CH_TIME("write_Level");
1064  int ret = 0;
1065 
1066  Interval comps(in_comps);
1067  if ( comps.size() == 0) comps = a_data.interval();
1068  T dummy; // used for preallocatable methods for dumb compilers.
1069  Vector<hid_t> types;
1070  dataTypes(types, dummy);
1071 
1072  Vector<Vector<long long> > offsets;
1073  Vector<long long> bufferCapacity(types.size(), 1); // noel (was 0)
1074  Vector<void*> buffers(types.size(), NULL);
1075 
1076  getOffsets(offsets, a_data, types.size(), comps, outputGhost);
1077 
1078  // create datasets collectively.
1079  hsize_t flatdims[1];
1080  char dataname[100];
1081  Vector<hid_t> dataspace(types.size());
1082  Vector<hid_t> dataset(types.size());
1083 
1084  herr_t err;
1085  hsize_t count[1];
1086  ch_offset_t offset[1];
1087  CH_assert(!(newForm && types.size() != 1));
1088 
1089  for (unsigned int i=0; i<types.size(); ++i)
1090  {
1091  flatdims[0] = offsets[i][offsets[i].size()-1];
1092  if (newForm)
1093  {
1094  if (a_name == "M")
1095  {
1096  strcpy(dataname, "Mask");
1097  }
1098  else
1099  {
1100  sprintf(dataname, "%sRegular", a_name.c_str());
1101  }
1102  }
1103  else
1104  {
1105  sprintf(dataname, "%s:datatype=%i",a_name.c_str(), i);
1106  }
1107  {
1108  CH_TIME("H5Screate");
1109  dataspace[i] = H5Screate_simple(1, flatdims, NULL);
1110  }
1111  CH_assert(dataspace[i] >=0);
1112  {
1113  CH_TIME("H5Dcreate");
1114 #ifdef H516
1115  dataset[i] = H5Dcreate(a_handle.groupID(), dataname,
1116  types[i],
1117  dataspace[i], H5P_DEFAULT);
1118 #else
1119  dataset[i] = H5Dcreate2(a_handle.groupID(), dataname,
1120  types[i],
1121  dataspace[i], H5P_DEFAULT,
1122  H5P_DEFAULT, H5P_DEFAULT);
1123 #endif
1124  }
1125  CH_assert(dataset[i] >= 0);
1126  }
1127 
1128  hid_t offsetspace, offsetData;
1129  for (unsigned int i=0; i<types.size(); ++i)
1130  {
1131  flatdims[0] = offsets[i].size();
1132  if (newForm)
1133  {
1134  if (a_name == "M")
1135  {
1136  strcpy(dataname, "MaskOffsets");
1137  }
1138  else
1139  {
1140  sprintf(dataname, "%sOffsets",a_name.c_str());
1141  }
1142  }
1143  else
1144  {
1145  sprintf(dataname, "%s:offsets=%i",a_name.c_str(), i);
1146  }
1147  {
1148  CH_TIME("H5S_H5D_offsets_create");
1149  offsetspace = H5Screate_simple(1, flatdims, NULL);
1150  CH_assert(offsetspace >= 0);
1151 #ifdef H516
1152  offsetData = H5Dcreate(a_handle.groupID(), dataname,
1153  H5T_NATIVE_LLONG, offsetspace,
1154  H5P_DEFAULT);
1155 #else
1156  offsetData = H5Dcreate2(a_handle.groupID(), dataname,
1157  H5T_NATIVE_LLONG, offsetspace,
1158  H5P_DEFAULT,
1159  H5P_DEFAULT,H5P_DEFAULT);
1160 #endif
1161  CH_assert(offsetData >= 0);
1162  }
1163  if (procID() == 0)
1164  {
1165  CH_TIME("WriteOffsets");
1166  hid_t memdataspace = H5Screate_simple(1, flatdims, NULL);
1167  CH_assert(memdataspace >= 0);
1168  err = H5Dwrite(offsetData, H5T_NATIVE_LLONG, memdataspace, offsetspace,
1169  H5P_DEFAULT, &(offsets[i][0]));
1170  CH_assert(err >= 0);
1171  H5Sclose(memdataspace);
1172  }
1173  {
1174  CH_TIME("H5S_H5D_offsets_close");
1175  H5Sclose(offsetspace);
1176  H5Dclose(offsetData);
1177  }
1178  }
1179 
1180  // write BoxLayoutData attributes into Dataset[0]
1181  if (!newForm)
1182  {
1183  CH_TIME("WriteAttributes");
1184  HDF5HeaderData info;
1185  info.m_int["comps"] = comps.size();
1186  info.m_string["objectType"] = name(dummy);
1187  info.m_intvect["outputGhost"] = outputGhost;
1188  info.m_intvect["ghost"] = outputGhost;
1189  std::string group = a_handle.getGroup();
1190  a_handle.setGroup(group+"/"+a_name+"_attributes");
1191  info.writeToFile(a_handle);
1192  a_handle.setGroup(group);
1193  }
1194 
1195  // collective operations finished, now perform parallel writes
1196  // to specified hyperslabs.
1197 
1198  Vector<size_t> type_size(types.size());
1199  for (unsigned int i=0; i<types.size(); ++i)
1200  {
1201  type_size[i] = H5Tget_size(types[i]);
1202  }
1203 
1204  Vector<int> thisSize(types.size());
1205 
1206  // Hooks for aggregated collective buffering (ACB).
1207  // Devendran, et al. Collective I/O Optimizations for Adaptive Mesh
1208  // Refinement Data Writes on Lustre File System, CUG 2016 for more on
1209  // ACB.
1210  // Comment out next line if you don't want ACB.
1211  // Also, comment out the line "#define TRY_MPI_COLLECTIVES_ 0" below if
1212  // you don't want MPI-IO collective buffering.
1213  // In general, you should turn on collective buffering if you use
1214  // ACB.
1215 #ifdef CH_MPI
1216 //#define AGGREGATE_BOXES_ 1
1217 #endif
1218 #ifdef AGGREGATE_BOXES_
1219  // pout() << "Turning on aggregated collective buffering (ACB)." << endl;
1220  // For each type types[i], aggCount[i] is the total number of data points for all
1221  // the boxes on this processor.
1222  Vector<hsize_t> aggCount(types.size(), 0);
1223  // size of aggregated buffers; initialize to 0 so we can add up
1224  // the box sizes to compute the size of the buffers
1225  Vector<long long> aggBufferSize(types.size(), 0);
1226  Vector<void*> aggBuffers(types.size(), NULL);
1227 #endif
1228 
1229  // step 1, create buffer big enough to hold the biggest linearized T
1230  // that I will have to output.
1231  // pout()<<"offsets ";
1232  // for (int i=0; i<offsets[0].size(); i++)
1233  // {
1234  // pout()<<" "<<offsets[0][i];
1235  // }
1236  // pout()<<"\n";
1237  {
1238  CH_TIME("ComputeBufferCapacity");
1239  for (DataIterator it = a_data.dataIterator(); it.ok(); ++it)
1240  {
1241  unsigned int index = a_data.boxLayout().index(it());
1242  for (unsigned int i=0; i<types.size(); ++i)
1243  {
1244  long long size = (offsets[i][index+1] - offsets[i][index])
1245  * type_size[i];
1246  // pout()<<"index:size "<<index<<" "<<size<<"\n";
1247  CH_assert(size >= 0);
1248  if (size > bufferCapacity[i]) // grow buffer if necessary.....
1249  {
1250  bufferCapacity[i] = size;
1251  }
1252 #ifdef AGGREGATE_BOXES_
1253  aggBufferSize[i] += size;
1254 #endif
1255  }
1256  }
1257  }
1258  // CH_assert(bufferCapacity[0] > 1);
1259  for (unsigned int i=0; i<types.size(); ++i)
1260  {
1261  CH_TIME("mallocMT_buffers");
1262  buffers[i] = mallocMT(bufferCapacity[i]);
1263  if (buffers[i] == NULL)
1264  {
1265  pout() << " i=" << i
1266  << " types.size() = " << types.size()
1267  << " bufferCapacity[i] = " << (int)bufferCapacity[i]
1268  << endl;
1269  MayDay::Error("memory error in buffer allocation write");
1270  }
1271 #ifdef AGGREGATE_BOXES_
1272  aggBuffers[i] = mallocMT(aggBufferSize[i]);
1273  if (aggBuffers[i] == NULL)
1274  {
1275  MayDay::Error("memory error in aggregate buffer allocation in write");
1276  }
1277 #endif
1278  }
1279 
1280 #ifdef CH_MPI
1281 //#define TRY_MPI_COLLECTIVES_ 1
1282 #endif
1283 #ifdef TRY_MPI_COLLECTIVES_
1284  // pout() << "Turned on MPI-IO collective buffering." << endl;
1285  // MPI collective buffering requires all processes call H5Dwrite collectively.
1286  // In particular, the processes must issue the same number of H5Dwrite calls,
1287  // or the application will hang when collective is turned on.
1288  // In the case where we do a separate write for each box, gather the
1289  // maximum number of boxes assigned to any process, so we know how many
1290  // H5Dwrites to do.
1291  DataIterator dit = a_data.dataIterator();
1292  int maxNumBoxes = dit.size();
1293  int sendBuf = maxNumBoxes;
1294  int result = MPI_Allreduce(&sendBuf, &maxNumBoxes, 1, MPI_INT,MPI_MAX, Chombo_MPI::comm);
1295  if (result != MPI_SUCCESS)
1296  {
1297  MayDay::Error("Couldn't collect maximum number of boxes!");
1298  }
1299  // Set dataset transfer property to collective mode. This is how we turn
1300  // on MPI-IO collective in HDF5.
1301  hid_t DXPL = H5Pcreate(H5P_DATASET_XFER);
1302  if(!(a_handle.openMode()==HDF5Handle::CREATE_SERIAL)) // can't set MPI collective if file was created for serial IO
1303  H5Pset_dxpl_mpio(DXPL, H5FD_MPIO_COLLECTIVE);
1304 #endif /*end TRY_MPI_COLLECTIVES_ */
1305 
1306  // Step 2. actually a) write each of my T objects into the
1307  // buffer, then b) write that buffered data out to the
1308  // write position in the data file using hdf5 hyperslab functions.
1309  {
1310  CH_TIME("linearize_H5Dwrite");
1311 #ifdef AGGREGATE_BOXES_
1312  // the first non-empty hyperslab needs to be handled separately from
1313  // the others, so keep track of the first non-empty hyperslab
1314  bool isFirstHyperslab = true;
1315 
1316  // bufferLoc keeps track of last location written to in the aggregrated
1317  // buffer
1318  Vector<void*> bufferLoc(types.size(), NULL);
1319  // Set the bufferLoc to the beginning of the aggregated buffer
1320  for(unsigned int i=0; i<types.size(); ++i)
1321  {
1322  bufferLoc[i] = aggBuffers[i];
1323  }
1324 #endif
1325  for (DataIterator it = a_data.dataIterator(); it.ok(); ++it)
1326  {
1327  const T& data = a_data[it()];
1328  unsigned int index = a_data.boxLayout().index(it());
1329  Box box = a_data.box(it());
1330  box.grow(outputGhost);
1331  // First, linearize the box, and put data into the buffer
1332  {
1333  CH_TIMELEAF("linearize");
1334 #ifdef AGGREGATE_BOXES_
1335  // Under aggregation, we need to be careful to write to
1336  // the buffer starting at bufferLoc.
1337  for(unsigned int i=0; i<types.size(); i++)
1338  {
1339  data.linearOut(bufferLoc[i], box, comps);
1340  char* tempLoc = ((char*)bufferLoc[i]) + data.size(box,comps);
1341  bufferLoc[i] = (void*)tempLoc;
1342  // bufferLoc[i] += data.size(box, comps);
1343  }
1344 #else
1345  write(data, buffers, box, comps); //write T to buffer
1346 #endif
1347  }
1348  // Next select HDF5 hyperslabs to specify where to write in HDF5 file
1349  // BUG: The hyperslab union generation is wrong if types.size() > 1
1350  // because isFirstHyperslab isn't independent for different sizes.
1351  // We don't make a fix because types.size() is almost always 1 (and
1352  // will probably be set to 1 in future redesigns).
1353  for (unsigned int i=0; i<types.size(); ++i)
1354  {
1355  offset[0] = offsets[i][index];
1356  count[0] = offsets[i][index+1] - offset[0];
1357 #ifdef AGGREGATE_BOXES_
1358  // Under aggregation, we may be sending multiple boxes to HDF5.
1359  // Because boxes on a process are not necessarily written to
1360  // contiguous locations in the HDF5 file, we need to specify
1361  // multiple file locations to HDF5. This is done using a
1362  // union of hyperslabs.
1363 
1364  aggCount[i] += count[0];
1365  if (isFirstHyperslab)
1366  {
1367  // If the box is non-empty, select the first hyperslab
1368  // and set isFirstHyperslab to false. Otherwise,
1369  // we haven't encountered the first non-empty
1370  // hyperslab, so keep isFirstHyperslab equal to true.
1371  if(count[0] > 0)
1372  {
1373  err = H5Sselect_hyperslab(dataspace[i], H5S_SELECT_SET,
1374  offset, NULL,
1375  count, NULL);
1376  CH_assert(err >= 0);
1377  isFirstHyperslab = false;
1378  }
1379  else // must explicitly tell HDF5 we are doing an empty write
1380  {
1381  H5Sselect_none(dataspace[i]);
1382  }
1383  }
1384  else
1385  {
1386  if(count[0] > 0)
1387  {
1388  // H5S_SELECT_OR creates a union of the selected
1389  // hyperslab with the already selected hyperslabs
1390  // in dataspace
1391  err = H5Sselect_hyperslab(dataspace[i], H5S_SELECT_OR,
1392  offset, NULL,
1393  count, NULL);
1394  CH_assert(err >= 0);
1395  }
1396  // else don't do H5Sselect_none in case it overrides
1397  // the already existing union of hyperslabs.
1398  }
1399 #else
1400  // Without aggregation, we create a simple hyperslab for the box.
1401  hid_t memdataspace=0;
1402  if (count[0] > 0)
1403  {
1404  err = H5Sselect_hyperslab(dataspace[i], H5S_SELECT_SET,
1405  offset, NULL,
1406  count, NULL);
1407  CH_assert(err >= 0);
1408  memdataspace = H5Screate_simple(1, count, NULL);
1409  CH_assert(memdataspace >= 0);
1410  }
1411  else // must explicitly tell HDF5 we are doing an empty write
1412  {
1413  H5Sselect_none(dataspace[i]);
1414  H5Sselect_none(memdataspace);
1415  }
1416 
1417  // Write out box if we are NOT performing aggregation.
1418  // (Under aggregation, one single write call is issued at the
1419  // end of the function.)
1420  {
1421  CH_TIMELEAF("H5Dwrite");
1422 #ifdef TRY_MPI_COLLECTIVES_
1423  err = H5Dwrite(dataset[i], types[i], memdataspace, dataspace[i],
1424  DXPL, buffers[i]);
1425 #else
1426  err = H5Dwrite(dataset[i], types[i], memdataspace, dataspace[i],
1427  H5P_DEFAULT, buffers[i]);
1428 #endif
1429  }
1430  CH_assert(err >= 0);
1431  H5Sclose(memdataspace);
1432  if (err < 0)
1433  {
1434  ret = err;
1435  pout() << "Before goto cleanup" << endl;
1436  goto cleanup;
1437  }
1438 #endif // end of ifdef AGGREGATE_BOXES_
1439  } // end of loop over types
1440  } // end of loop over data iterator
1441 
1442  // Under aggregation, now we issue one write call to send all boxes to HDF5
1443  // If aggregation is turned off, but MPI collective is turned on, we may
1444  // have to do empty writes to make sure all processes call H5Dwrite
1445  // the same number of times.
1446 #ifdef AGGREGATE_BOXES_
1447  for(unsigned int i=0; i<types.size(); ++i)
1448  {
1449  if(aggCount[i] > 0)
1450  {
1451  hid_t memdataspace = 0;
1452  memdataspace = H5Screate_simple(1, &(aggCount[i]), NULL);
1453  CH_assert(memdataspace >= 0);
1454  {
1455  CH_TIMELEAF("H5Dwrite");
1456 #ifdef TRY_MPI_COLLECTIVES_
1457  err = H5Dwrite(dataset[i], types[i], memdataspace, dataspace[i],
1458  DXPL, aggBuffers[i]);
1459 #else
1460  err = H5Dwrite(dataset[i], types[i], memdataspace, dataspace[i],
1461  H5P_DEFAULT, aggBuffers[i]);
1462 #endif
1463  }
1464  H5Sclose(memdataspace);
1465  if (err < 0)
1466  {
1467  ret = err;
1468  pout() << "Error! goto cleanup" << endl;
1469  goto cleanup;
1470  }
1471  }
1472  //else aggCount[i] is 0, and this processor has no data to write.
1473  // For MPI collectives, still have to do an empty write call
1474 #ifdef TRY_MPI_COLLECTIVES_
1475  else
1476  {
1477  hid_t memdataspace = 0;
1478  memdataspace = H5Screate_simple(1, &(aggCount[i]), NULL);
1479  H5Sselect_none(memdataspace);
1480  H5Sselect_none(dataspace[i]);
1481  err = H5Dwrite(dataset[i], types[i], memdataspace, dataspace[i],
1482  DXPL, aggBuffers[i]);
1483  if (err < 0)
1484  {
1485  ret = err;
1486  pout() << "Before goto cleanup" << endl;
1487  goto cleanup;
1488  }
1489  }
1490 #endif
1491  } // end for loop over types
1492 #else // not using aggregated collective buffering
1493 #ifdef TRY_MPI_COLLECTIVES_
1494  // MPI collectives expects all processes to make the same number of H5Dwrite calls,
1495  // or it will hang. So, call H5Dwrite with empty data
1496  // First create memdataspace according to example in
1497  // https://www.hdfgroup.org/ftp/HDF5/examples/parallel/coll_test.c
1498  hid_t memdataspace = 0;
1499  // Setting first argument to 1 b/c that's the value used for non-empty writes.
1500  // (See H5Sselect_hyperslab code above.)
1501  memdataspace = H5Screate_simple(1, count, NULL);
1502  H5Sselect_none(memdataspace);
1503  int nBoxes = a_data.dataIterator().size();
1504  for(int iwrite = nBoxes; iwrite < maxNumBoxes; iwrite++)
1505  {
1506  for (unsigned int i=0; i<types.size(); ++i)
1507  {
1508  H5Sselect_none(dataspace[i]);
1509  // for debugging: buffers has junk data in it (but none of that data should
1510  // be written out here)
1511  err = H5Dwrite(dataset[i], types[i], memdataspace, dataspace[i],
1512  DXPL, buffers[i]);
1513  if (err < 0)
1514  {
1515  ret = err;
1516  pout() << "Before goto cleanup" << endl;
1517  goto cleanup;
1518  }
1519  }
1520  }
1521  H5Sclose(memdataspace);
1522 
1523 #endif // end of #ifdef TRY_MPI_COLLECTIVES_
1524 #endif // end of #ifdef AGGREGATE_BOXES_
1525 
1526  } // end of region for CH_TIME("linearize_H5Dwrite")
1527 
1528 #ifdef TRY_MPI_COLLECTIVES_
1529  H5Pclose(DXPL);
1530 #endif
1531 
1532  // OK, clean up data structures
1533 
1534  cleanup:
1535  for (unsigned int i=0; i<types.size(); ++i)
1536  {
1537  {
1538  CH_TIME("freeMT");
1539  freeMT(buffers[i]);
1540  #ifdef AGGREGATE_BOXES_
1541  freeMT(aggBuffers[i]);
1542  #endif
1543  }
1544  {
1545  CH_TIME("H5Sclose");
1546  H5Sclose(dataspace[i]);
1547  }
1548  {
1549  CH_TIME("H5Dclose");
1550  H5Dclose(dataset[i]);
1551  }
1552  }
1553  return ret;
1554 
1555 }
1556 
1557 template <class T>
1558 int write(HDF5Handle& a_handle, const LevelData<T>& a_data,
1559  const std::string& a_name, const IntVect& outputGhost, const Interval& in_comps)
1560 {
1561  CH_TIMERS("Write Level");
1562  CH_TIMER("calc minimum in outputGhost",t1);
1563  CH_TIMER("writeToFile",t2);
1564  CH_TIMER("setGroup",t3);
1565  HDF5HeaderData info;
1566  info.m_intvect["ghost"] = a_data.ghostVect();
1567  IntVect og(outputGhost);
1568  CH_START(t1);
1569  og.min(a_data.ghostVect());
1570  CH_STOP(t1);
1571  info.m_intvect["outputGhost"] = og;
1572  std::string group = a_handle.getGroup();
1573  a_handle.setGroup(group+"/"+a_name+"_attributes");
1574  CH_START(t2);
1575  info.writeToFile(a_handle);
1576  CH_STOP(t2);
1577  CH_START(t3);
1578  a_handle.setGroup(group);
1579  CH_STOP(t3);
1580  return write(a_handle, (const BoxLayoutData<T>&)a_data, a_name, og, in_comps);
1581 }
1582 
1583 template <class T>
1584 int read(HDF5Handle& a_handle, LevelData<T>& a_data, const std::string& a_name,
1585  const DisjointBoxLayout& a_layout, const Interval& a_comps, bool a_redefineData)
1586 {
1587  if (a_redefineData)
1588  {
1589  HDF5HeaderData info;
1590  std::string group = a_handle.getGroup();
1591  if (a_handle.setGroup(group+"/"+a_name+"_attributes"))
1592  {
1593  std::string message = "error opening "
1594  +a_handle.getGroup()+"/"+a_name+"_attributes" ;
1595  MayDay::Warning(message.c_str());
1596  return 1;
1597  }
1598  info.readFromFile(a_handle);
1599  a_handle.setGroup(group);
1600  int ncomp = info.m_int["comps"];
1601  IntVect ghost = info.m_intvect["ghost"];
1602  if (a_comps.end() > 0 && ncomp < a_comps.end())
1603  {
1604  MayDay::Error("attempt to read component interval that is not available");
1605  }
1606  if (a_comps.size() == 0)
1607  a_data.define(a_layout, ncomp, ghost);
1608  else
1609  a_data.define(a_layout, a_comps.size(), ghost);
1610  }
1611  return read(a_handle, (BoxLayoutData<T>&)a_data, a_name, a_layout, a_comps, false);
1612 
1613 }
1614 
1615 template <class T>
1616 int read(HDF5Handle& a_handle, BoxLayoutData<T>& a_data, const std::string& a_name,
1617  const BoxLayout& a_layout, const Interval& a_comps, bool a_redefineData)
1618 {
1619  CH_TIME("read");
1620  int ret = 0; // return value;
1621 
1622  herr_t err;
1623 
1624  char dataname[100];
1625  hsize_t count[1];
1626  ch_offset_t offset[1];
1627  Vector<Vector<long long> > offsets;
1628 
1629  T dummy;
1630  Vector<hid_t> types;
1631  dataTypes(types, dummy);
1632  Vector<hid_t> dataspace(types.size());
1633  Vector<hid_t> dataset(types.size());
1634  offsets.resize(types.size(), Vector<long long>(a_layout.size() +1));
1635 
1636  //Vector<int> bufferCapacity(types.size(), 500);
1637  //Vector<void*> buffers(types.size(), NULL);
1638  Vector<Vector<char> > buffers(types.size(), 500);
1639 
1640  for (unsigned int i=0; i<types.size(); ++i)
1641  {
1642  sprintf(dataname, "%s:datatype=%i",a_name.c_str(), i);
1643 #ifdef H516
1644  dataset[i] = H5Dopen(a_handle.groupID(), dataname);
1645 #else
1646  dataset[i] = H5Dopen2(a_handle.groupID(), dataname, H5P_DEFAULT);
1647 #endif
1648  if (dataset[i] < 0)
1649  {
1650  MayDay::Warning("dataset open failure"); return dataset[i];
1651  }
1652  dataspace[i] = H5Dget_space(dataset[i]);
1653  if (dataspace[i] < 0)
1654  {
1655  MayDay::Warning("dataspace open failure"); return dataspace[i];
1656  }
1657  }
1658 
1659  hid_t offsetspace, offsetData;
1660  hsize_t flatdims[1];
1661  for (unsigned int i=0; i<types.size(); ++i)
1662  {
1663  flatdims[0] = offsets[i].size();
1664  sprintf(dataname, "%s:offsets=%i",a_name.c_str(), i);
1665  offsetspace = H5Screate_simple(1, flatdims, NULL);
1666  CH_assert(offsetspace >= 0);
1667 #ifdef H516
1668  offsetData = H5Dopen(a_handle.groupID(), dataname);
1669 #else
1670  offsetData = H5Dopen2(a_handle.groupID(), dataname,H5P_DEFAULT);
1671 #endif
1672  CH_assert(offsetData >= 0);
1673  hid_t memdataspace = H5Screate_simple(1, flatdims, NULL);
1674  CH_assert(memdataspace >= 0);
1675  err = H5Dread(offsetData, H5T_NATIVE_LLONG, memdataspace, offsetspace,
1676  H5P_DEFAULT, &(offsets[i][0]));
1677  CH_assert(err >=0);
1678  H5Sclose(memdataspace);
1679  H5Sclose(offsetspace);
1680  H5Dclose(offsetData);
1681  }
1682 
1683  HDF5HeaderData info;
1684  std::string group = a_handle.getGroup();
1685  if (a_handle.setGroup(a_handle.getGroup()+"/"+a_name+"_attributes"))
1686  {
1687  std::string message = "error opening "+a_handle.getGroup()+"/"+a_name ;
1688  MayDay::Warning(message.c_str());
1689  return 1;
1690  }
1691 
1692  info.readFromFile(a_handle);
1693  a_handle.setGroup(group);
1694  int ncomps = info.m_int["comps"];
1695  IntVect outputGhost(IntVect::Zero); // backwards file compatible mode.
1696  if (info.m_intvect.find("outputGhost") != info.m_intvect.end())
1697  {
1698  outputGhost = info.m_intvect["outputGhost"];
1699  }
1700  if (ncomps <= 0)
1701  {
1702  MayDay::Warning("ncomps <= 0 in read");
1703  return ncomps;
1704  }
1705 
1706  if (a_redefineData)
1707  {
1708  if (a_comps.size() != 0)
1709  {
1710  a_data.define(a_layout, a_comps.size());
1711  }
1712  else
1713  {
1714  a_data.define(a_layout, ncomps);
1715  }
1716  }
1717 
1718  Interval comps(0, ncomps-1);
1719  //huh?
1720  // if (a_comps.size() != 0) comps = Interval(0, a_comps.size());
1721 
1722  // getOffsets(offsets, a_data, types.size(), comps);
1723 
1724  Vector<size_t> type_size(types.size());
1725  for (unsigned int i=0; i<types.size(); ++i)
1726  {
1727  type_size[i] = H5Tget_size(types[i]);
1728 
1729  }
1730 
1731  Vector<int> thisSize(types.size());
1732  DataIterator it = a_data.dataIterator();
1733  for ( ; it.ok(); ++it)
1734  {
1735  CH_TIMELEAF("H5Dread");
1736  T& data = a_data[it()];
1737  unsigned int index = a_data.boxLayout().index(it());
1738  Box box = a_data.box(it());
1739 
1740  for (unsigned int i=0; i<types.size(); ++i)
1741  {
1742 
1743  offset[0] = offsets[i][index];
1744  count[0] = offsets[i][index+1] - offset[0];
1745  if (count[0] > 0)
1746  {
1747  size_t size = count[0] * type_size[i];
1748  while (size > buffers[i].size())
1749  {
1750  buffers[i].resize(2*buffers[i].size());
1751  }
1752 
1753  err = H5Sselect_hyperslab(dataspace[i], H5S_SELECT_SET,
1754  offset, NULL,
1755  count, NULL);
1756  CH_assert(err >= 0);
1757  hid_t memdataspace = H5Screate_simple(1, count, NULL);
1758  CH_assert(memdataspace >= 0);
1759  err = H5Dread(dataset[i], types[i], memdataspace, dataspace[i],
1760  H5P_DEFAULT, &(buffers[i][0]));
1761  CH_assert(err >= 0);
1762  H5Sclose(memdataspace);
1763  if (err < 0)
1764  {
1765  ret = err;
1766  goto cleanup;
1767  }
1768  }
1769  }
1770  box.grow(outputGhost);
1771  read(data, buffers, box, comps);
1772  }
1773 // if (it.size()==0)
1774 // {
1775 // // try doing an empty H5Dread operation to make the collective system happier.
1776 // for (unsigned int i=0; i<types.size(); ++i)
1777 // {
1778 // hid_t filespace, memspace;
1779 // H5Sselect_none(filespace);
1780 // H5Sselect_none(memspace);
1781 // err = H5Dread(dataset[i], types[i], memspace, filespace,
1782 // H5P_DEFAULT, 0);
1783 // H5Sclose(filespace);
1784 // H5Sclose(memspace);
1785 // }
1786 // }
1787 
1788  cleanup:
1789  for (unsigned int i=0; i<types.size(); ++i)
1790  {
1791  // freeMT(buffers[i]);
1792  H5Sclose(dataspace[i]);
1793  H5Dclose(dataset[i]);
1794  }
1795  return ret;
1796 }
1797 
1798 template <class T>
1799 int writeLevel(HDF5Handle& a_handle,
1800  const int& a_level,
1801  const T& a_data,
1802  const Real& a_dx,
1803  const Real& a_dt,
1804  const Real& a_time,
1805  const Box& a_domain,
1806  const int& a_refRatio,
1807  const IntVect& outputGhost,
1808  const Interval& comps)
1809 {
1810  int error;
1811  char levelName[10];
1812  std::string currentGroup = a_handle.getGroup();
1813  sprintf(levelName, "/level_%i",a_level);
1814  error = a_handle.setGroup(currentGroup + levelName);
1815  if (error != 0) return 1;
1816 
1817  HDF5HeaderData meta;
1818  meta.m_real["dx"] = a_dx;
1819  meta.m_real["dt"] = a_dt;
1820  meta.m_real["time"] = a_time;
1821  meta.m_box["prob_domain"] = a_domain;
1822  meta.m_int["ref_ratio"] = a_refRatio;
1823 
1824  error = meta.writeToFile(a_handle);
1825  if (error != 0) return 2;
1826 
1827  error = write(a_handle, a_data.boxLayout());
1828  if (error != 0) return 3;
1829 
1830  error = write(a_handle, a_data, "data", outputGhost, comps);
1831  if (error != 0) return 4;
1832 
1833  a_handle.setGroup(currentGroup);
1834 
1835  return 0;
1836 }
1837 
1838 template <class T>
1839 int writeLevel(HDF5Handle& a_handle,
1840  const int& a_level,
1841  const T& a_data,
1842  const RealVect& a_dx, // dx for each direction
1843  const Real& a_dt,
1844  const Real& a_time,
1845  const Box& a_domain,
1846  const IntVect& a_refRatios, // ref ratio for each direction
1847  const IntVect& outputGhost,
1848  const Interval& comps)
1849 {
1850  int error;
1851  char levelName[10];
1852  std::string currentGroup = a_handle.getGroup();
1853  sprintf(levelName, "/level_%i",a_level);
1854  error = a_handle.setGroup(currentGroup + levelName);
1855  if (error != 0) return 1;
1856 
1857  HDF5HeaderData meta;
1858  meta.m_realvect["vec_dx"] = a_dx;
1859  meta.m_real["dt"] = a_dt;
1860  meta.m_real["time"] = a_time;
1861  meta.m_box["prob_domain"] = a_domain;
1862  meta.m_intvect["vec_ref_ratio"] = a_refRatios;
1863 
1864  error = meta.writeToFile(a_handle);
1865  if (error != 0) return 2;
1866 
1867  error = write(a_handle, a_data.boxLayout());
1868  if (error != 0) return 3;
1869 
1870  error = write(a_handle, a_data, "data", outputGhost, comps);
1871  if (error != 0) return 4;
1872 
1873  a_handle.setGroup(currentGroup);
1874 
1875  return 0;
1876 }
1877 
1878 template <class T>
1879 int readLevel(HDF5Handle& a_handle,
1880  const int& a_level,
1881  LevelData<T>& a_data,
1882  Real& a_dx,
1883  Real& a_dt,
1884  Real& a_time,
1885  Box& a_domain,
1886  int& a_refRatio,
1887  const Interval& a_comps,
1888  bool setGhost)
1889 {
1890  HDF5HeaderData header;
1891  header.readFromFile(a_handle);
1892  //unused
1893  // int nComp = header.m_int["num_components"];
1894 
1895  int error;
1896  char levelName[10];
1897  std::string currentGroup = a_handle.getGroup();
1898  sprintf(levelName, "/level_%i",a_level);
1899  error = a_handle.setGroup(currentGroup + levelName);
1900  if (error != 0) return 1;
1901 
1902  HDF5HeaderData meta;
1903  error = meta.readFromFile(a_handle);
1904  if (error != 0) return 2;
1905  a_dx = meta.m_real["dx"];
1906  a_dt = meta.m_real["dt"];
1907  a_time = meta.m_real["time"];
1908  a_domain = meta.m_box["prob_domain"];
1909  a_refRatio = meta.m_int["ref_ratio"];
1910  Vector<Box> boxes;
1911  error = read(a_handle, boxes);
1912  Vector<int> procIDs;
1913  LoadBalance(procIDs, boxes);
1914 
1915  DisjointBoxLayout layout(boxes, procIDs, a_domain);
1916 
1917  layout.close();
1918  if (error != 0) return 3;
1919 
1920  error = read(a_handle, a_data, "data", layout, a_comps, true);
1921 
1922  if (error != 0) return 4;
1923 
1924  a_handle.setGroup(currentGroup);
1925 
1926  return 0;
1927 }
1928 
1929 template <class T>
1930 int readLevel(HDF5Handle& a_handle,
1931  const int& a_level,
1932  LevelData<T>& a_data,
1933  RealVect& a_dx,
1934  Real& a_dt,
1935  Real& a_time,
1936  Box& a_domain,
1937  IntVect& a_refRatio,
1938  const Interval& a_comps,
1939  bool setGhost)
1940 {
1941  HDF5HeaderData header;
1942  header.readFromFile(a_handle);
1943  //unused
1944  // int nComp = header.m_int["num_components"];
1945 
1946  int error;
1947  char levelName[10];
1948  std::string currentGroup = a_handle.getGroup();
1949  sprintf(levelName, "/level_%i",a_level);
1950  error = a_handle.setGroup(currentGroup + levelName);
1951  if (error != 0) return 1;
1952 
1953  HDF5HeaderData meta;
1954  error = meta.readFromFile(a_handle);
1955  if (error != 0) return 2;
1956  // a_dx = meta.m_realvect["vec_dx"];
1957  // Allow for vec_dx to be absent, as it will be in isotropic files.
1958  if (meta.m_realvect.find("vec_dx") != meta.m_realvect.end())
1959  {
1960  a_dx = meta.m_realvect["vec_dx"];
1961  }
1962  else
1963  { // vec_dx is not present, so get dx.
1964  Real dxScalar = meta.m_real["dx"];
1965  a_dx = dxScalar * RealVect::Unit;
1966  }
1967  a_dt = meta.m_real["dt"];
1968  a_time = meta.m_real["time"];
1969  a_domain = meta.m_box["prob_domain"];
1970  // a_refRatio = meta.m_intvect["vec_ref_ratio"];
1971  // Allow for vec_ref_ratio to be absent, as it will be in isotropic files.
1972  if (meta.m_intvect.find("vec_ref_ratio") != meta.m_intvect.end())
1973  {
1974  a_refRatio = meta.m_intvect["vec_ref_ratio"];
1975  }
1976  else
1977  { // vec_ref_ratio is not present, so get ref_ratio.
1978  int refRatioScalar = meta.m_int["ref_ratio"];
1979  a_refRatio = refRatioScalar * IntVect::Unit;
1980  }
1981  Vector<Box> boxes;
1982  error = read(a_handle, boxes);
1983  Vector<int> procIDs;
1984  LoadBalance(procIDs, boxes);
1985 
1986  DisjointBoxLayout layout(boxes, procIDs, a_domain);
1987 
1988  layout.close();
1989  if (error != 0) return 3;
1990 
1991  error = read(a_handle, a_data, "data", layout, a_comps, true);
1992 
1993  if (error != 0) return 4;
1994 
1995  a_handle.setGroup(currentGroup);
1996 
1997  return 0;
1998 }
1999 
2000 
2001 /*******************************************************************************
2002  *
2003  * Class WriteMultiData: member definitions
2004  *
2005  ******************************************************************************/
2006 
2007 /*--------------------------------------------------------------------*/
2008 // Constructor writes boxes and allocates dataset for LevelData
2009 /** Write boxes, processors, attributes and creates a dataset for all
2010  * the level data to be written. The dataset is closed on
2011  * destruction.
2012  * \param[in] a_handle
2013  * Chombo HDF5 handle holding open file for
2014  * writing
2015  * \param[in] a_layout
2016  * Box layout for all data to be written
2017  * \param[in] a_numIntv
2018  * Total number of components that will be
2019  * written from all BoxLayoutData or LevelData
2020  * \param[in] a_name Name of the dataset for the level data.
2021  * \param[in] a_policyFlags
2022  * Either
2023  * CH_HDF5::IOPolicyDefault - Use internal
2024  * T.linearOut(...) to linearize data and
2025  * processes independently write to file,
2026  * but still in parallel.
2027  * or a union (|) of the following flags:
2028  * CH_HDF5::IOPolicyMultiDimHyperslab - The
2029  * memory dataspace will be set up so that
2030  * HDF5 linearizes the data. Using this
2031  * requires T.dataPtr() and contiguous memory
2032  * (so things like T=FluxBox will not work).
2033  * CH_HDF5::IOPolicyCollectiveWrite - The write
2034  * for each dataset will be collective.
2035  * \param[in] a_outputGhost
2036  * Number of ghost cells that will be written to
2037  * the data file. Any data written must have
2038  * this many ghosts
2039  * \param[in] a_newForm
2040  * ?
2041  *//*-----------------------------------------------------------------*/
2042 
2043 template <class T>
2045  const BoxLayout& a_layout,
2046  const int a_numIntv,
2047  const string& a_name,
2048  const int a_policyFlags,
2049  const IntVect& a_outputGhost,
2050  const bool a_newForm)
2051  :
2052  m_policyFlags(a_policyFlags),
2053  m_outputGhost(a_outputGhost),
2054  m_newForm(a_newForm)
2055 {
2056  CH_TIME("WriteMultiData::constructor");
2057  CH_assert(T::preAllocatable() == 0);
2058 #ifdef H516
2059  CH_assert(false); // API 1.6 not supported
2060 #endif
2061 
2062 //--Write the boxes
2063 
2064  write(a_handle, a_layout);
2065 
2066 //--Create the dataset for the level data
2067 
2068  T dummy; // Used for preallocatable methods for dumb compilers.
2069  dataTypes<T>(m_types, dummy);
2070  // Only allow one type
2071  CH_assert(m_types.size() == 1);
2072 
2073  m_allIntvFile.define(0, a_numIntv-1);
2074  getOffsets<T>(m_offsets,
2075  a_layout,
2076  m_types.size(),
2077  m_allIntvFile,
2078  m_outputGhost);
2079 
2080  hsize_t flatdims[1];
2081  // Create dataset for data
2082  {
2083  hid_t dataspace;
2084  flatdims[0] = m_offsets[0].back();
2085  if (m_newForm)
2086  {
2087  if (a_name == "M")
2088  {
2089  strcpy(m_dataname, "Mask");
2090  }
2091  else
2092  {
2093  sprintf(m_dataname, "%sRegular", a_name.c_str());
2094  }
2095  }
2096  else
2097  {
2098  sprintf(m_dataname, "%s:datatype=%i", a_name.c_str(), 0);
2099  }
2100  {
2101  CH_TIME("H5Screate");
2102  dataspace = H5Screate_simple(1, flatdims, NULL);
2103  }
2104  CH_assert(dataspace >=0);
2105  {
2106  CH_TIME("H5Dcreate");
2107  m_dataSet = H5Dcreate2(a_handle.groupID(), m_dataname,
2108  m_types[0],
2109  dataspace, H5P_DEFAULT,
2110  H5P_DEFAULT, H5P_DEFAULT);
2111  }
2112  CH_assert(m_dataSet >= 0);
2113  {
2114  CH_TIME("H5S_H5D_data_close");
2115  H5Sclose(dataspace);
2116  }
2117  }
2118 
2119  // Create dataset for offsets
2120  {
2121  hid_t offsetspace, offsetdataset;
2122  char offsetname[128];
2123  flatdims[0] = m_offsets[0].size(); // Number of boxes
2124  if (m_newForm)
2125  {
2126  if (a_name == "M")
2127  {
2128  strcpy(offsetname, "MaskOffsets");
2129  }
2130  else
2131  {
2132  sprintf(offsetname, "%sOffsets", a_name.c_str());
2133  }
2134  }
2135  else
2136  {
2137  sprintf(offsetname, "%s:offsets=%i", a_name.c_str(), 0);
2138  }
2139  {
2140  CH_TIME("H5S_H5D_offsets_create");
2141  offsetspace = H5Screate_simple(1, flatdims, NULL);
2142  CH_assert(offsetspace >= 0);
2143  offsetdataset = H5Dcreate2(a_handle.groupID(), offsetname,
2144  H5T_NATIVE_LLONG, offsetspace,
2145  H5P_DEFAULT,
2146  H5P_DEFAULT,H5P_DEFAULT);
2147  CH_assert(offsetdataset >= 0);
2148  }
2149  if (procID() == 0)
2150  // Write the offsets
2151  {
2152  herr_t err;
2153  CH_TIME("WriteOffsets");
2154  hid_t memdataspace = H5Screate_simple(1, flatdims, NULL);
2155  CH_assert(memdataspace >= 0);
2156  err = H5Dwrite(offsetdataset, H5T_NATIVE_LLONG, memdataspace,
2157  offsetspace, H5P_DEFAULT, &(m_offsets[0][0]));
2158  CH_assert(err >= 0);
2159  H5Sclose(memdataspace);
2160  }
2161  {
2162  CH_TIME("H5S_H5D_offsets_close");
2163  H5Sclose(offsetspace);
2164  H5Dclose(offsetdataset);
2165  }
2166  }
2167 
2168  // Write attributes
2169  {
2170  CH_TIME("WriteAttributes");
2171  HDF5HeaderData info;
2172  info.m_intvect["ghost"] = m_outputGhost;
2173  info.m_intvect["outputGhost"] = m_outputGhost;
2174  if (!m_newForm)
2175  {
2176  info.m_int["comps"] = m_allIntvFile.size();
2177  info.m_string["objectType"] = name(dummy);
2178  }
2179  std::string group = a_handle.getGroup();
2180  a_handle.setGroup(group+"/"+a_name+"_attributes");
2181  info.writeToFile(a_handle);
2182  a_handle.setGroup(group);
2183  }
2184 
2185  // Get the maximum boxes per process
2186  {
2187  m_maxBoxPerProc = a_layout.dataIterator().size();
2188 #ifdef CH_MPI
2189  long myCountBoxes = m_maxBoxPerProc;
2190  MPI_Allreduce(&myCountBoxes,
2191  &m_maxBoxPerProc,
2192  1, MPI_LONG, MPI_MAX, Chombo_MPI::comm);
2193 #endif
2194  }
2195 }
2196 
2197 /*--------------------------------------------------------------------*/
2198 // Write an interval of LevelData to the dataset
2199 /** Call this as many times as you want to write level data. User is
2200  * resposible for making sure 'a_intvFile' has relevance to
2201  * 'a_numIntv' specified during construction.
2202  * \param[in] a_data Level data to write
2203  * \param[in] a_intvMem
2204  * Interval of a_data to read
2205  * \param[in] a_intvFile
2206  * Interval to write to file
2207  *//*-----------------------------------------------------------------*/
2208 
2209 template <class T>
2210 int
2212  const Interval& a_intvMem,
2213  const Interval& a_intvFile)
2214 {
2215  CH_TIME("WriteMultiData::writeData");
2216  CH_assert(a_intvMem.size() > 0);
2217  CH_assert(a_intvMem.size() == a_intvFile.size());
2218  CH_assert(a_intvMem.begin() >= a_data.interval().begin() &&
2219  a_intvMem.end() <= a_data.interval().end());
2220  CH_assert(a_intvFile.begin() >= m_allIntvFile.begin() &&
2221  a_intvFile.end() <= m_allIntvFile.end());
2222 
2223  int ret = 0;
2224 
2225  hsize_t fileCount[1], memCount[SpaceDim+1];
2226  ch_offset_t fileOffset[1], memOffset[SpaceDim+1]; // Same type as hsize_t
2227  herr_t err;
2228 
2229 //--Open the file dataspace
2230 
2231  hid_t fileDataSpace = H5Dget_space(m_dataSet);
2232 
2233 //--Data transfer policies (independent or collective)
2234 
2235  hid_t DXPL = H5Pcreate(H5P_DATASET_XFER);
2236 #ifdef CH_MPI
2237  if (m_policyFlags & CH_HDF5::IOPolicyCollectiveWrite)
2238  {
2239  H5Pset_dxpl_mpio(DXPL, H5FD_MPIO_COLLECTIVE);
2240  }
2241 #endif
2242 
2243  if (m_policyFlags & CH_HDF5::IOPolicyMultiDimHyperslab)
2244 
2245 //--Perform parallel writes from specified hyperslabs (when defining hyperslabs
2246 //--in memory, remember that HDF5 assumes row-ordering).
2247 
2248  {
2249  CH_TIME("hyperslab_H5Dwrite");
2250  DataIterator dit = a_data.dataIterator();
2251  const void* buffer;
2252  for (int iBox = 0; iBox != m_maxBoxPerProc; ++iBox)
2253  {
2254  hid_t memDataSpace;
2255  if (dit.ok())
2256  {
2257  const T& data = a_data[dit];
2258  unsigned globalIdx = a_data.boxLayout().index(dit());
2259  Box dbox = data.box(); // Data box
2260  Box rbox = a_data.box(dit()); // Read box
2261  rbox.grow(m_outputGhost);
2262  // Create the dataspace in memory
2263  memCount[0] = a_data.nComp();
2264  D_TERM6(memCount[SpaceDim-0] = dbox.size(0);,
2265  memCount[SpaceDim-1] = dbox.size(1);,
2266  memCount[SpaceDim-2] = dbox.size(2);,
2267  memCount[SpaceDim-3] = dbox.size(3);,
2268  memCount[SpaceDim-4] = dbox.size(4);,
2269  memCount[SpaceDim-5] = dbox.size(5);)
2270  memDataSpace = H5Screate_simple(SpaceDim+1, memCount, NULL);
2271  // Select the hyperslab from the memory dataspace
2272  memCount[0] = a_intvMem.size();
2273  D_TERM6(memCount[SpaceDim-0] = rbox.size(0);,
2274  memCount[SpaceDim-1] = rbox.size(1);,
2275  memCount[SpaceDim-2] = rbox.size(2);,
2276  memCount[SpaceDim-3] = rbox.size(3);,
2277  memCount[SpaceDim-4] = rbox.size(4);,
2278  memCount[SpaceDim-5] = rbox.size(5);)
2279  memOffset[0] = a_intvMem.begin();
2280  D_TERM6(
2281  memOffset[SpaceDim-0] = rbox.smallEnd(0) - dbox.smallEnd(0);,
2282  memOffset[SpaceDim-1] = rbox.smallEnd(1) - dbox.smallEnd(1);,
2283  memOffset[SpaceDim-2] = rbox.smallEnd(2) - dbox.smallEnd(2);,
2284  memOffset[SpaceDim-3] = rbox.smallEnd(3) - dbox.smallEnd(3);,
2285  memOffset[SpaceDim-4] = rbox.smallEnd(4) - dbox.smallEnd(4);,
2286  memOffset[SpaceDim-5] = rbox.smallEnd(5) - dbox.smallEnd(5);)
2287  err = H5Sselect_hyperslab(memDataSpace, H5S_SELECT_SET,
2288  memOffset, NULL, memCount, NULL);
2289  CH_assert(err >= 0);
2290  // Create the hyperslab in the file dataspace
2291  fileOffset[0] = m_offsets[0][globalIdx];
2292  fileCount[0] = m_offsets[0][globalIdx + 1] - fileOffset[0];
2293  if (fileCount[0] > 0) // Else catches more processes than boxes
2294  {
2295  // Revise offsets based on selection of interval
2296  const hsize_t dpnts = rbox.numPts(); // Points per comp
2297  fileOffset[0] += dpnts*a_intvFile.begin();
2298  fileCount[0] = dpnts*a_intvFile.size();
2299  CH_assert(fileOffset[0] + fileCount[0] <=
2300  m_offsets[0][globalIdx+1]);
2301  err = H5Sselect_hyperslab(fileDataSpace, H5S_SELECT_SET,
2302  fileOffset, NULL, fileCount, NULL);
2303  CH_assert(err >= 0);
2304  }
2305  else // More processes than boxes
2306  {
2307  H5Sselect_none(memDataSpace);
2308  H5Sselect_none(fileDataSpace);
2309  }
2310  buffer = data.dataPtr();
2311  ++dit;
2312  }
2313  else // Does not have a box to participate in collective :(
2314  {
2315  std::memset(memCount, 0, (SpaceDim+1)*sizeof(hsize_t));
2316  memDataSpace = H5Screate_simple(SpaceDim+1, memCount, NULL);
2317  H5Sselect_none(memDataSpace);
2318  H5Sselect_none(fileDataSpace);
2319  buffer = 0;
2320  }
2321  // Write
2322  err = H5Dwrite(m_dataSet, m_types[0], memDataSpace, fileDataSpace,
2323  DXPL, buffer);
2324  CH_assert(err >= 0);
2325  H5Sclose(memDataSpace);
2326  if (err < 0)
2327  {
2328  ret = err;
2329  goto cleanup;
2330  }
2331  }
2332  }
2333  else
2334 
2335 //--Perform parallel writes with 1-D hyperslabs and using internal mechanisms
2336 //--for linearization.
2337 
2338  {
2339  CH_TIME("linearize_H5Dwrite");
2340  void* linearBuffer;
2341 
2342  // Step 1: create a buffer to linearize T
2343  {
2344  CH_TIMELEAF("allocateLinearBuffer");
2345  long long bufferSize = 0;
2346  for (DataIterator dit = a_data.dataIterator(); dit.ok(); ++dit)
2347  {
2348  unsigned globalIdx = a_data.boxLayout().index(dit());
2349  {
2350  const long long allIntvSize =
2351  (m_offsets[0][globalIdx + 1] - m_offsets[0][globalIdx]);
2352  CH_assert(allIntvSize >= 0);
2353  if (allIntvSize > bufferSize)
2354  {
2355  bufferSize = allIntvSize;
2356  }
2357  }
2358  }
2359  // Get buffer size in bytes, and adjusted for write interval
2360  bufferSize = ((bufferSize/m_allIntvFile.size())*a_intvMem.size())*
2361  H5Tget_size(m_types[0]);
2362  linearBuffer = mallocMT(bufferSize);
2363  if (linearBuffer == NULL)
2364  {
2365  pout() << " bufferCapacity = " << (int)bufferSize << endl;
2366  MayDay::Error("Memory error in buffer allocation write");
2367  }
2368  }
2369 
2370  // Step 2: linearize the data and then write to file
2371  DataIterator dit = a_data.dataIterator();
2372  for (int iBox = 0; iBox != m_maxBoxPerProc; ++iBox)
2373  {
2374  hid_t memDataSpace;
2375  if (dit.ok())
2376  {
2377  const T& data = a_data[dit];
2378  unsigned globalIdx = a_data.boxLayout().index(dit());
2379  Box rbox = a_data.box(dit()); // Read box
2380  rbox.grow(m_outputGhost);
2381  {
2382  CH_TIMELEAF("linearize");
2383  data.linearOut(linearBuffer, rbox, a_intvMem);
2384  }
2385  const hsize_t dpnts = rbox.numPts(); // Points per comp
2386  // Create the dataspace in memory
2387  memCount[0] = dpnts*a_intvMem.size();
2388  memDataSpace = H5Screate_simple(1, memCount, NULL);
2389  // Create the hyperslab in the file dataspace
2390  fileOffset[0] = m_offsets[0][globalIdx];
2391  fileCount[0] = m_offsets[0][globalIdx + 1] - fileOffset[0];
2392  if (fileCount[0] > 0) // Else catches more processes than boxes
2393  {
2394  // Revise offsets based on selection of interval
2395  fileOffset[0] += dpnts*a_intvFile.begin();
2396  fileCount[0] = dpnts*a_intvFile.size();
2397  CH_assert(fileOffset[0] + fileCount[0] <=
2398  m_offsets[0][globalIdx+1]);
2399  err = H5Sselect_hyperslab(fileDataSpace, H5S_SELECT_SET,
2400  fileOffset, NULL, fileCount, NULL);
2401  CH_assert(err >= 0);
2402  }
2403  else // More processes than boxes
2404  {
2405  H5Sselect_none(memDataSpace);
2406  H5Sselect_none(fileDataSpace);
2407  }
2408  ++dit;
2409  }
2410  else // Does not have a box to participate in collective :(
2411  {
2412  memCount[0] = 0;
2413  memDataSpace = H5Screate_simple(SpaceDim+1, memCount, NULL);
2414  H5Sselect_none(memDataSpace);
2415  H5Sselect_none(fileDataSpace);
2416  }
2417  // Collective write
2418  err = H5Dwrite(m_dataSet, m_types[0], memDataSpace, fileDataSpace,
2419  DXPL, linearBuffer);
2420  CH_assert(err >= 0);
2421  H5Sclose(memDataSpace);
2422  if (err < 0)
2423  {
2424  ret = err;
2425  freeMT(linearBuffer);
2426  goto cleanup;
2427  }
2428  }
2429  freeMT(linearBuffer);
2430  }
2431 
2432  cleanup: ;
2433  H5Pclose(DXPL);
2434  H5Sclose(fileDataSpace);
2435  // m_dataSet closed in destructor
2436  return ret;
2437 }
2438 
2439 #include "NamespaceFooter.H"
2440 
2441 #else // CH_USE_HDF5
2442 
2443 // this is the only thing needed when HDF is not used
2444 #define HOFFSET(S,M) (offsetof(S,M))
2445 
2446 #endif // CH_USE_HDF5 not defined
2447 #endif // CH_HDF5_H
std::ostream & pout()
Use this in place of std::cout for program output.
int open(const std::string &a_filename, mode a_mode, const char *a_globalGroupName="Chombo_global")
{ File functions}
#define CH_TIMERS(name)
Definition: CH_Timer.H:133
map< std::string, int > m_int
Definition: CH_HDF5.H:548
virtual void define(const BoxLayout &boxes, int comps, const DataFactory< T > &factory=DefaultDataFactory< T >())
Definition: BoxLayoutDataI.H:87
static herr_t attributeScan(hid_t loc_id, const char *name, void *opdata)
IntVect & min(const IntVect &p)
Definition: IntVect.H:1134
#define D_TERM6(a, b, c, d, e, f)
Definition: CHArray.H:40
#define freeMT(a_a)
Definition: memtrack.H:160
#define CH_SPACEDIM
Definition: SPACE.H:51
#define CH_assert(cond)
Definition: CHArray.H:37
int popGroup()
data to be added to HDF5 files.
Definition: CH_HDF5.H:517
const hid_t & fileID() const
mode m_mode
Definition: CH_HDF5.H:462
int readBoxes(HDF5Handle &a_handle, Vector< Vector< Box > > &boxes)
reads the set of Boxes out from the level_* groups of a Chombo HDF5 AMR file
mode
Definition: CH_HDF5.H:314
Definition: CH_HDF5.H:62
void writeDataset(hid_t a_dataset, hid_t a_dataspace, const void *start, ch_offset_t off, hsize_t count)
A not-necessarily-disjoint collective of boxes.
Definition: BoxLayout.H:145
one dimensional dynamic array
Definition: Vector.H:53
int nComp() const
Definition: BoxLayoutData.H:306
map< std::string, Real > m_real
Definition: CH_HDF5.H:545
const IntVect m_outputGhost
Number of ghost cells written.
Definition: CH_HDF5.H:665
int readFArrayBox(HDF5Handle &a_handle, FArrayBox &a_fab, int a_level, int a_boxNumber, const Interval &a_components, const std::string &a_dataName="data")
FArrayBox-at-a-time read function. FArrayBox gets redefined in the function. Reads data field named b...
int writeToLocation(hid_t loc_id) const
std::string m_group
Definition: CH_HDF5.H:465
#define mallocMT(a_a)
Definition: memtrack.H:159
LayoutIterator layoutIterator() const
Iterator that processes through ALL the boxes in a BoxLayout.
#define CH_START(tpointer)
Definition: CH_Timer.H:145
void getOffsets(Vector< Vector< long long > > &offsets, const BoxLayoutData< T > &a_data, int types, const Interval &comps, const IntVect &outputGhost)
Definition: CH_HDF5.H:929
int size() const
Definition: DataIterator.H:218
void readDataset(hid_t a_dataset, hid_t a_dataspace, void *start, ch_offset_t off, hsize_t count)
IntVect size() const
size functions
Definition: Box.H:1819
void close()
void setGroupToLevel(int a_level)
{ Group functions}
Definition: DataIterator.H:190
void linearOut(void *const a_outBuf, const CH_XDIR::OffsetBuffer &a_inputT)
int readLevel(HDF5Handle &a_handle, const int &a_level, LevelData< T > &a_data, Real &a_dx, Real &a_dt, Real &a_time, Box &a_domain, int &a_refRatio, const Interval &a_comps=Interval(), bool setGhost=false)
Definition: CH_HDF5.H:1879
Definition: CH_HDF5.H:319
unsigned int size() const
Returns the total number of boxes in the BoxLayout.
std::ostream & operator<<(std::ostream &os, const HDF5HeaderData &data)
herr_t HDF5HeaderDataattributeScan(hid_t loc_id, const char *name, const H5A_info_t *info, void *opdata)
int LoadBalance(Vector< Vector< int > > &a_procAssignments, Real &a_effRatio, const Vector< Vector< Box > > &a_Grids, const Vector< Vector< long > > &a_ComputeLoads, const Vector< int > &a_RefRatios, int a_nProc=numProc())
int writeData(const BoxLayoutData< T > &a_data, const Interval &a_intvMem, const Interval &a_intvFile)
Write an interval of LevelData to the dataset.
Definition: CH_HDF5.H:2211
unsigned int numProc()
number of parallel processes
int uniqueProc(const SerialTask::task &a_task)
An Iterator based on a BoxLayout object.
Definition: LayoutIterator.H:35
void createDataset(hid_t &a_dataset, hid_t &a_dataspace, HDF5Handle &handle, const std::string &name, const T *dummy, hsize_t size)
Definition: CH_HDF5.H:737
int size() const
Definition: Interval.H:75
hid_t m_dataSet
Dataset for level data.
Definition: CH_HDF5.H:669
void dataTypes(Vector< hid_t > &a_types, const BaseFab< int > &dummy)
Definition: CH_HDF5.H:799
void dataSize(const BaseFab< int > &item, Vector< int > &a_sizes, const Box &box, const Interval &comps)
Definition: CH_HDF5.H:830
const int SpaceDim
Definition: SPACE.H:38
Definition: CH_HDF5.H:316
map< std::string, RealVect > m_realvect
Definition: CH_HDF5.H:560
#define CH_TIMER(name, tpointer)
Definition: CH_Timer.H:63
static const RealVect Unit
Definition: RealVect.H:427
Definition: EBInterface.H:45
static hid_t realvect_id
Definition: CH_HDF5.H:452
void resize(unsigned int isize)
Definition: Vector.H:346
static void initialize()
Definition: CH_HDF5.H:318
virtual void close()
map< std::string, IntVect > m_intvect
Definition: CH_HDF5.H:554
Interval m_allIntvFile
Interval for components in file.
Definition: CH_HDF5.H:670
A FArrayBox-like container for face-centered fluxes.
Definition: FluxBox.H:22
void gather(Vector< T > &a_outVec, const T &a_input, int a_dest)
Definition: SPMDI.H:197
void linearIn(CH_XDIR::OffsetBuffer &a_outputT, const void *const a_inBuf)
IOPolicy
Definition: CH_HDF5.H:59
char m_dataname[128]
Name for level data dataset.
Definition: CH_HDF5.H:668
Methods for writing multiple LevelData to an HDF5 file.
Definition: CH_HDF5.H:629
void push_back(const T &in)
Definition: Vector.H:295
Vector< hid_t > m_types
Type of data written.
Definition: CH_HDF5.H:678
static bool initialized
Definition: CH_HDF5.H:469
#define CH_TIME(name)
Definition: CH_Timer.H:82
Structure for passing component ranges in code.
Definition: Interval.H:23
Definition: CH_HDF5.H:64
static const IntVect Unit
Definition: IntVect.H:663
const char * name(const FArrayBox &a_dummySpecializationArg)
Definition: CH_HDF5.H:864
int linearSize(const CH_XDIR::OffsetBuffer &a_input)
new code
Definition: BoxLayoutData.H:170
long m_maxBoxPerProc
Maximum boxes written by any proc.
Definition: CH_HDF5.H:671
std::string m_filename
Definition: CH_HDF5.H:464
HDF5Handle::mode openMode() const
Definition: CH_HDF5.H:447
Interval interval() const
Definition: BoxLayoutData.H:312
Data on a BoxLayout.
Definition: BoxLayoutData.H:97
double Real
Definition: REAL.H:33
Box surroundingNodes(const Box &b, int dir)
Definition: Box.H:2161
unsigned int index(const LayoutIndex &index) const
Definition: BoxLayout.H:724
Definition: SPMD.H:281
Definition: CH_HDF5.H:762
virtual void define(const DisjointBoxLayout &dp, int comps, const IntVect &ghost=IntVect::Zero, const DataFactory< T > &a_factory=DefaultDataFactory< T >())
Definition: LevelDataI.H:80
void operator=(const OffsetBuffer &rhs)
int pushGroup(const std::string &grp)
bool isOpen() const
Vector< int > index
Definition: CH_HDF5.H:764
Definition: CH_HDF5.H:317
A BoxLayout that has a concept of disjointedness.
Definition: DisjointBoxLayout.H:30
void createData(hid_t &a_dataset, hid_t &a_dataspace, HDF5Handle &handle, const std::string &name, hid_t type, hsize_t size)
static void Error(const char *const a_msg=m_nullString, int m_exitCode=CH_DEFAULT_ERROR_CODE)
Print out message to cerr and exit with the specified exit code.
int write(HDF5Handle &a_handle, const BoxLayout &a_layout, const std::string &name="boxes")
writes BoxLayout to HDF5 file.
int begin() const
Definition: Interval.H:99
const BoxLayout & boxLayout() const
Definition: LayoutData.H:107
const int m_policyFlags
Policies.
Definition: CH_HDF5.H:664
static const IntVect Zero
Definition: IntVect.H:658
static hid_t box_id
Definition: CH_HDF5.H:450
const hid_t & groupID() const
Vector< Vector< int > > offsets
Definition: CH_HDF5.H:765
hssize_t ch_offset_t
Definition: CH_HDF5.H:693
static void Warning(const char *const a_msg=m_nullString)
Print out message to cerr and continue.
A Rectangular Domain on an Integer Lattice.
Definition: Box.H:469
A Real vector in SpaceDim-dimensional space.
Definition: RealVect.H:41
void define(int a_firstComp, int a_lastComp)
Definition: Interval.H:52
size_t numPts() const
hid_t m_fileID
Definition: CH_HDF5.H:460
int writeToFile(HDF5Handle &file) const
WriteMultiData & operator=(const WriteMultiData)
Vector< Vector< long long > > m_offsets
Definition: CH_HDF5.H:675
hid_t H5Type(const T *dummy)
Definition: CH_HDF5.H:722
map< std::string, Box > m_box
Definition: CH_HDF5.H:557
#define CH_STOP(tpointer)
Definition: CH_Timer.H:150
Handle to a particular group in an HDF file.
Definition: CH_HDF5.H:292
void dump() const
useful for debugging. dumps contents to std::cout
An integer Vector in SpaceDim-dimensional space.
Definition: CHArray.H:42
Definition: FArrayBox.H:45
const std::string & getGroup() const
DataIterator dataIterator() const
Definition: LayoutDataI.H:78
size_t size() const
Definition: Vector.H:192
HDF5Handle()
{ constructor}
WriteMultiData(HDF5Handle &a_handle, const BoxLayout &a_layout, const int a_numIntv, const string &a_name, const int a_policyFlags=CH_HDF5::IOPolicyDefault, const IntVect &a_outputGhost=IntVect::Zero, const bool a_newForm=false)
Constructor writes boxes and allocates dataset for LevelData.
Definition: CH_HDF5.H:2044
Definition: CH_HDF5.H:61
#define CH_TIMELEAF(name)
Definition: CH_Timer.H:100
Box & grow(int i)
grow functions
Definition: Box.H:2263
virtual bool ok() const
return true if this iterator is still in its Layout
Definition: LayoutIterator.H:117
int setGroup(const std::string &groupAbsPath)
void broadcast(T &a_inAndOut, int a_src)
broadcast to every process
Definition: SPMDI.H:207
int writeLevel(HDF5Handle &a_handle, const int &a_level, const T &a_data, const Real &a_dx, const Real &a_dt, const Real &a_time, const Box &a_domain, const int &a_refRatio, const IntVect &outputGhost=IntVect::Zero, const Interval &comps=Interval())
user-friendly function to write out data on a AMR level
Definition: CH_HDF5.H:1799
int readFromLocation(hid_t loc_id)
int end() const
Definition: Interval.H:104
map< std::string, std::string > m_string
Definition: CH_HDF5.H:551
static hid_t intvect_id
Definition: CH_HDF5.H:451
const IntVect & smallEnd() const
{ Accessors}
Definition: Box.H:1770
int procID()
local process ID
Box box(const DataIndex &a_index) const
Definition: LayoutDataI.H:66
int readFromFile(HDF5Handle &file)
HDF5Handle & operator=(const HDF5Handle &)
bool m_isOpen
Definition: CH_HDF5.H:463
static map< std::string, std::string > groups
Definition: CH_HDF5.H:453
void read(T &item, Vector< Vector< char > > &a_allocatedBuffers, const Box &box, const Interval &comps)
Definition: CH_HDF5.H:49
hid_t m_currentGroupID
Definition: CH_HDF5.H:461
const bool m_newForm
?
Definition: CH_HDF5.H:666
#define H5T_NATIVE_REAL
Definition: REAL.H:35
int m_level
Definition: CH_HDF5.H:466
const IntVect & ghostVect() const
Definition: LevelData.H:170
~WriteMultiData()
Destructor.
Definition: CH_HDF5.H:644