1 /**
2  * MOAB, a Mesh-Oriented datABase, is a software component for creating,
3  * storing and accessing finite element mesh data.
4  *
5  * Copyright 2004 Sandia Corporation.  Under the terms of Contract
6  * DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government
7  * retains certain rights in this software.
8  *
9  * This library is free software; you can redistribute it and/or
10  * modify it under the terms of the GNU Lesser General Public
11  * License as published by the Free Software Foundation; either
12  * version 2.1 of the License, or (at your option) any later version.
13  *
14  */
15 
16 
17 #ifndef WRITE_HDF5_HPP
18 #define WRITE_HDF5_HPP
19 
20 #include <list>
21 #include "moab/MOABConfig.h"
22 #ifdef MOAB_HAVE_MPI // include this before HDF5 headers to avoid conflicts
23 #  include "moab_mpi.h"
24 #endif
25 #include "moab_mpe.h"
26 #include "mhdf.h"
27 #include "moab/Forward.hpp"
28 #include "moab/Range.hpp"
29 #include "moab/WriterIface.hpp"
30 #include "moab/RangeMap.hpp"
31 #include "moab/WriteUtilIface.hpp"
32 #include "DebugOutput.hpp"
33 #include "HDF5Common.hpp"
34 
35 namespace moab {
36 
37 class IODebugTrack;
38 
39 /* If this define is not set, node->entity adjacencies will not be written */
40 #undef MB_H5M_WRITE_NODE_ADJACENCIES
41 
42 /**
43  * \brief  Write mesh database to MOAB's native HDF5-based file format.
44  * \author Jason Kraftcheck
45  * \date   01 April 2004
46  */
47 class WriteHDF5 : public WriterIface
48 {
49 
50 public:
51 
52   static WriterIface* factory( Interface* );
53 
54   WriteHDF5( Interface* iface );
55 
56   virtual ~WriteHDF5();
57 
58   /** Export specified meshsets to file
59    * \param filename     The filename to export.
60    * \param export_sets  Array of handles to sets to export, or NULL to export all.
61    * \param export_set_count Length of <code>export_sets</code> array.
62    */
63   ErrorCode write_file( const char* filename,
64                           const bool overwrite,
65                           const FileOptions& opts,
66                           const EntityHandle* export_sets,
67                           const int export_set_count,
68                           const std::vector<std::string>& qa_records,
69                           const Tag* tag_list = NULL,
70                           int num_tags = 0,
71                           int user_dimension = 3 );
72 
73   /** The type to use for entity IDs w/in the file.
74    *
75    * NOTE:  If this is changed, the value of id_type
76    *        MUST be changed accordingly.
77    */
78   typedef EntityHandle wid_t; // change the name,
79   //  to avoid conflicts to /usr/include/x86_64-linux-gnu/bits/types.h : id_t , which is unsigned int
80 
81   /** HDF5 type corresponding to type of wid_t */
82   static const hid_t id_type;
83 
84   struct ExportType
85   {
86       //! The type of the entities in the range
87     EntityType type;
88       //! The number of nodes per entity - not used for nodes and sets
89     int num_nodes;
90 
~ExportTypemoab::WriteHDF5::ExportType91     virtual ~ExportType()
92     { }
93 
operator ==moab::WriteHDF5::ExportType94     bool operator==(ExportType t) const
95       { return t.type == type && t.num_nodes == num_nodes; }
operator !=moab::WriteHDF5::ExportType96     bool operator!=(ExportType t) const
97       { return t.type != type || t.num_nodes != num_nodes; }
operator <moab::WriteHDF5::ExportType98     bool operator<(ExportType t) const
99       { return type < t.type || (type == t.type && num_nodes < t.num_nodes); }
100   };
101 
102   //! Range of entities, grouped by type, to export
103   struct ExportSet : public ExportType
104   {
105     //! The range of entities.
106     Range range;
107     //! The first Id allocated by the mhdf library.  Entities in range have sequential IDs.
108     wid_t first_id;
109     //! The offset at which to begin writing this processor's data.
110     //! Always zero except for parallel IO.
111     long offset;
112     //! Offset for adjacency data.  Always zero except for parallel IO
113     long adj_offset;
114     //! If doing parallel IO, largest number of entities to write
115     //! for any processor (needed to do collective IO).  Zero if unused.
116     long max_num_ents, max_num_adjs;
117     //! The total number of entities that will be written to the file
118     //! for this group.  For serial IO, this should always be range.size().
119     //! For parallel IO, it will be the sum of range size over all processors.
120     //! For parallel IO, this value is undefined except for on the root
121     //! processor.
122     long total_num_ents;
123 
operator <moab::WriteHDF5::ExportSet124     bool operator<( const ExportType& other ) const
125       { return type < other.type ||
126                (type == other.type && num_nodes < other.num_nodes); }
127 
operator <moab::WriteHDF5::ExportSet128     bool operator<( std::pair<int,int> other ) const
129       { return type < other.first ||
130                (type == other.first && num_nodes < other.second); }
131 
operator ==moab::WriteHDF5::ExportSet132     bool operator==( const ExportType& other ) const
133       { return (type == other.type && num_nodes == other.num_nodes); }
134 
operator ==moab::WriteHDF5::ExportSet135     bool operator==( std::pair<int,int> other ) const
136       { return (type == other.first && num_nodes == other.second); }
137 
138     const char* name() const;
139   };
140 
141   //! Tag to write to file.
142   struct TagDesc
143   {
144     //! The tag handle
145     Tag tag_id;
146     //! The offset at which to begin writting this processor's data.
147     //! Always zero except for parallel IO.
148     wid_t sparse_offset;
149     //! For variable-length tags, a second offset for the tag data table,
150     //! separate from the offset used for the ID and Index tables.
151     //! Always zero except for parallel IO.
152     wid_t var_data_offset;
153     //! Write sparse tag data (for serial, is always equal to !range.empty())
154     bool write_sparse;
155     //! If doing parallel IO, largest number, over all processes, of entities
156     //! for which to write tag data.  Zero if unused.
157     unsigned long max_num_ents;
158     //! For variable-length tags during parallel IO: the largest number
159     //! of tag values to be written on by any process, used to calculate
160     //! the total number of collective writes that all processes must do.
161     //! Zero for fixed-length tags or if not doing parallel IO.
162     unsigned long max_num_vals;
163 
164     //! List of entity groups for which to write tag data in
165     //! dense format
166     std::vector<ExportType> dense_list;
167 
have_densemoab::WriteHDF5::TagDesc168     bool have_dense( const ExportType& type ) const
169       { return std::find(dense_list.begin(), dense_list.end(), type) != dense_list.end(); }
170 
171     bool operator<(const TagDesc&) const;
172   };
173 
174   /** Create attributes holding the HDF5 type handle for the
175    *  type of a bunch of the default tags.
176    */
177   //static ErrorCode register_known_tag_types( Interface* );
178 
179   //! Store old HDF5 error handling function
180   struct HDF5ErrorHandler {
181     HDF5_Error_Func_Type func;
182     void* data;
183   };
184 
file_ptr()185   mhdf_FileHandle file_ptr() { return filePtr; }
186 
write_util()187   WriteUtilIface* write_util() { return writeUtil; }
188 
189 protected:
190 
191   //! Store old HDF5 error handling function
192   HDF5ErrorHandler errorHandler;
193 
194   /** Function to create the file.  Virtual to allow override
195    *  for parallel version.
196    */
197   virtual ErrorCode parallel_create_file( const char* filename,
198                                             bool overwrite,
199                                             const std::vector<std::string>& qa_records,
200                                             const FileOptions& opts,
201                                             const Tag* tag_list,
202                                             int num_tags,
203                                             int dimension = 3,
204                                             double* times = 0 );
205   virtual ErrorCode write_finished();
206   virtual void debug_barrier_line(int lineno);
207 
208   //! Gather tags
209   ErrorCode gather_tags( const Tag* user_tag_list, int user_tag_list_length );
210 
211   /** Check if tag values for a given ExportSet should be written in dense format
212    *
213    *\param ents        ExportSet to consider
214    *\param all_tagged  Range containing all the entities in ents.range for
215    *                   which an explicit tag value is stored.  Range may
216    *                   also contain entities not in ents.range, but may
217    *                   not contain entities in ents.range for which no tag
218    *                   value is stored.
219    *\param prefer_dense If true, will return true if at least 2/3 of the
220    *                   entities are tagged.  This should not be passed as
221    *                   true if the tag does not have a default value, as
222    *                   tag values must be stored for all entities in the
223    *                   ExportSet for dense-formatted data.
224    */
225   bool check_dense_format_tag( const ExportSet& ents,
226                                const Range& all_tagged,
227                                bool prefer_dense );
228 
229   /** Helper function for create-file
230    *
231    * Calculate the sum of the number of non-set adjacencies
232    * of all entities in the passed range.
233    */
234   ErrorCode count_adjacencies( const Range& elements, wid_t& result );
235 
236 public: // make these public so helper classes in WriteHDF5Parallel can use them
237 
238   /** Helper function for create-file
239    *
240    * Create zero-ed tables where element connectivity and
241    * adjacency data will be stored.
242    */
243   ErrorCode create_elem_table( const ExportSet& block, long num_ents, long& first_id_out );
244 
245   /** Helper function for create-file
246    *
247    * Create zero-ed table where set descriptions will be written
248    */
249   ErrorCode create_set_meta( long num_sets, long& first_id_out );
250 
251 protected:
252 
253   /** Helper function for create-file
254    *
255    * Calculate total length of set contents and child tables.
256    */
257   ErrorCode count_set_size( const Range& sets,
258                               long& contents_length_out,
259                               long& children_length_out,
260                               long& parents_length_out );
261 
262   //! Get information about a meshset
263   ErrorCode get_set_info( EntityHandle set,
264                             long& num_entities,
265                             long& num_children,
266                             long& num_parents,
267                             unsigned long& flags );
268 
269   /** Helper function for create-file
270    *
271    * Create zero-ed tables where set data will be written.
272    */
273   ErrorCode create_set_tables( long contents_length,
274                                long children_length,
275                                long parents_length );
276 
277   //! Write exodus-type QA info
278   ErrorCode write_qa( const std::vector<std::string>& list );
279 
280   //!\brief Get tagged entities for which to write tag values
281   ErrorCode get_num_sparse_tagged_entities( const TagDesc& tag, size_t& count );
282   //!\brief Get tagged entities for which to write tag values
283   ErrorCode get_sparse_tagged_entities( const TagDesc& tag, Range& range );
284   //!\brief Get entities that will be written to file
285   void get_write_entities( Range& range );
286 
287   //! The size of the data buffer (<code>dataBuffer</code>).
288   size_t bufferSize;
289   //! A memory buffer to use for all I/O operations.
290   char* dataBuffer;
291 
292   //! Interface pointer passed to constructor
293   Interface* iFace;
294   //! Cached pointer to writeUtil interface.
295   WriteUtilIface* writeUtil;
296 
297   //! The file handle from the mhdf library
298   mhdf_FileHandle filePtr;
299 
300   //! Map from entity handles to file IDs
301   RangeMap<EntityHandle,wid_t> idMap;
302 
303   //! The list elements to export.
304   std::list<ExportSet> exportList;
305   //! The list of nodes to export
306   ExportSet nodeSet;
307   //! The list of sets to export
308   ExportSet setSet;
309 
find(ExportType type) const310   const ExportSet* find( ExportType type ) const {
311     if (type.type == MBVERTEX)
312       return &nodeSet;
313     else if (type.type == MBENTITYSET)
314       return &setSet;
315     else {
316       std::list<ExportSet>::const_iterator it;
317       it = std::find( exportList.begin(), exportList.end(), type );
318       return it == exportList.end() ? 0 : &*it;
319     }
320   }
321 
322   //! Offset into set contents table (zero except for parallel)
323   unsigned long setContentsOffset;
324   //! Offset into set children table (zero except for parallel)
325   unsigned long setChildrenOffset, setParentsOffset;
326   //! The largest number of values to write
327   //! for any processor (needed to do collective IO).
328   long maxNumSetContents, maxNumSetChildren, maxNumSetParents;
329   //! Flags idicating if set data should be written.
330   //! For the normal (non-parallel) case, these values
331   //! will depend only on whether or not there is any
332   //! data to be written.  For parallel-meshes, opening
333   //! the data table is collective so the values must
334   //! depend on whether or not any processor has meshsets
335   //! to be written.
336   bool writeSets, writeSetContents, writeSetChildren, writeSetParents;
337 
338   //! Struct describing a set for which the contained and linked entity
339   //! lists are something other than the local values.  Used to store
340   //! data for shared sets owned by this process when writing in parallel.
341   struct SpecialSetData {
342     EntityHandle setHandle;
343     unsigned setFlags;
344     std::vector<wid_t> contentIds;
345     std::vector<wid_t> childIds;
346     std::vector<wid_t> parentIds;
347   };
348   struct SpecSetLess {
operator ()moab::WriteHDF5::SpecSetLess349     bool operator() (const SpecialSetData& a, SpecialSetData b) const
350       { return a.setHandle < b.setHandle; }
351   };
352 
353 
354   //! Array of special/shared sets, in order of handle value.
355   std::vector<SpecialSetData> specialSets;
find_set_data(EntityHandle h) const356   const SpecialSetData* find_set_data( EntityHandle h ) const
357     { return const_cast<WriteHDF5*>(this)->find_set_data(h); }
358   SpecialSetData* find_set_data( EntityHandle h );
359 
360   //! The list of tags to export
361   std::list<TagDesc> tagList;
362 
363   //! True if doing parallel write
364   bool parallelWrite;
365   //! True if using collective IO calls for parallel write
366   bool collectiveIO;
367   //! True if writing dense-formatted tag data
368   bool writeTagDense;
369 
370   //! Property set to pass to H5Dwrite calls.
371   //! For serial, should be H5P_DEFAULTS.
372   //! For parallel, may request collective IO.
373   hid_t writeProp;
374 
375   //! Utility to log debug output
376   DebugOutput dbgOut;
377 
378   static MPEState topState;
379   static MPEState subState;
380 
381   //! Look for overlapping and/or missing writes
382   bool debugTrack;
383 
384   void print_id_map() const;
385   void print_id_map( std::ostream& str, const char* prefix = "" ) const;
386 
387   /** Helper function for create-file
388    *
389    * Write tag meta-info and create zero-ed table where
390    * tag values will be written.
391    *\param num_entities  Number of entities for which to write tag data.
392    *\param var_len_total For variable-length tags, the total number of values
393    *                     in the data table.
394    */
395   ErrorCode create_tag( const TagDesc& tag_data,
396                         unsigned long num_entities,
397                         unsigned long var_len_total );
398 
399   /**\brief add entities to idMap */
400   ErrorCode assign_ids( const Range& entities, wid_t first_id );
401 
402   /** Get possibly compacted list of IDs for passed entities
403    *
404    * For the passed range of entities, determine if IDs
405    * can be compacted and write IDs to passed list.
406    *
407    * If the IDs are not compacted, the output list will contain
408    * a simple ordered list of IDs.
409    *
410    * If IDs are compacted, the output list will contain
411    * {start,count} pairs.
412    *
413    * If the ID list is compacted, ranged_list will be 'true'.
414    * Otherwise it will be 'false'.
415    */
416   ErrorCode range_to_blocked_list( const Range& input_range,
417                                    std::vector<wid_t>& output_id_list ,
418                                    bool& ranged_list );
419 
420   /** Get possibly compacted list of IDs for passed entities
421    *
422    * For the passed range of entities, determine if IDs
423    * can be compacted and write IDs to passed list.
424    *
425    * If the IDs are not compacted, the output list will contain
426    * a simple ordered list of IDs.
427    *
428    * If IDs are compacted, the output list will contain
429    * {start,count} pairs.
430    *
431    * If the ID list is compacted, ranged_list will be 'true'.
432    * Otherwise it will be 'false'.
433    */
434   ErrorCode range_to_blocked_list( const EntityHandle* input_ranges,
435                                    size_t num_input_ranges,
436                                    std::vector<wid_t>& output_id_list ,
437                                    bool& ranged_list );
438 
439 
440   ErrorCode range_to_id_list( const Range& input_range,
441                                 wid_t* array );
442   //! Get IDs for entities
443   ErrorCode vector_to_id_list( const std::vector<EntityHandle>& input,
444                                std::vector<wid_t>& output,
445                                bool remove_non_written = false );
446   //! Get IDs for entities
447   ErrorCode vector_to_id_list( const EntityHandle* input,
448                                wid_t* output,
449                                size_t num_entities );
450   //! Get IDs for entities
451   ErrorCode vector_to_id_list( const EntityHandle* input,
452                                size_t input_len,
453                                wid_t* output,
454                                size_t& output_len,
455                                bool remove_non_written );
456 
457   /** When writing tags containing EntityHandles to file, need to convert tag
458    *  data from EntityHandles to file IDs.  This function does that.
459    *
460    * If the handle is not valid or does not correspond to an entity that will
461    * be written to the file, the file ID is set to zero.
462    *\param data  The data buffer.  As input, an array of EntityHandles.  As
463    *             output an array of file IDS, where the size of each integral
464    *             file ID is the same as the size of EntityHandle.
465    *\param count The number of handles in the buffer.
466    *\return true if at least one of the handles is valid and will be written to
467    *             the file or at least one of the handles is NULL (zero). false
468    *             otherwise
469    */
470   bool convert_handle_tag( EntityHandle* data, size_t count ) const;
471   bool convert_handle_tag( const EntityHandle* source,
472                            EntityHandle* dest,
473                            size_t count ) const;
474 
475   /** Get IDs of adjacent entities.
476    *
477    * For all entities adjacent to the passed entity, if the
478    * adjacent entity is to be exported (ID is not zero), append
479    * the ID to the passed list.
480    */
481   ErrorCode get_adjacencies( EntityHandle entity, std::vector<wid_t>& adj );
482 
483   //! get sum of lengths of tag values (as number of type) for
484   //! variable length tag data.
485   ErrorCode get_tag_data_length( const TagDesc& tag_info,
486                                  const Range& range,
487                                  unsigned long& result );
488 
489 private:
490 
491   //! Do the actual work of write_file.  Separated from write_file
492   //! for easier resource cleanup.
493   ErrorCode write_file_impl( const char* filename,
494                                const bool overwrite,
495                                const FileOptions& opts,
496                                const EntityHandle* export_sets,
497                                const int export_set_count,
498                                const std::vector<std::string>& qa_records,
499                                const Tag* tag_list,
500                                int num_tags,
501                                int user_dimension = 3 );
502 
503   ErrorCode init();
504 
505   ErrorCode serial_create_file( const char* filename,
506                                   bool overwrite,
507                                   const std::vector<std::string>& qa_records,
508                                   const Tag* tag_list,
509                                   int num_tags,
510                                   int dimension = 3 );
511 
512   /** Get all mesh to export from given list of sets.
513    *
514    * Populate exportSets, nodeSet and setSet with lists of
515    * entities to write.
516    *
517    * \param export_sets  The list of meshsets to export
518    */
519   ErrorCode gather_mesh_info( const std::vector<EntityHandle>& export_sets );
520 
521   //! Same as gather_mesh_info, except for entire mesh
522   ErrorCode gather_all_mesh( );
523 
524   //! Initialize internal data structures from gathered mesh
525   ErrorCode initialize_mesh( const Range entities_by_dim[5] );
526 
527   /** Write out the nodes.
528    *
529    * Note: Assigns IDs to nodes.
530    */
531   ErrorCode write_nodes( );
532 
533   /** Write out element connectivity.
534    *
535    * Write connectivity for passed set of elements.
536    *
537    * Note: Assigns element IDs.
538    * Note: Must do write_nodes first so node IDs get assigned.
539    */
540   ErrorCode write_elems( ExportSet& elemset );
541 
542   /** Write out meshsets
543    *
544    * Write passed set of meshsets, including parent/child relations.
545    *
546    * Note: Must have written nodes and element connectivity
547    *       so entities have assigned IDs.
548    */
549   ErrorCode write_sets( double* times );
550 
551   /** Write set contents/parents/children lists
552    *
553    *\param which_data Which set data to write (contents, parents, or children)
554    *\param handle     HDF5 handle for data set in which to write data
555    *\param track      Debugging tool
556    *\param ranged     Will be populated with handles of sets for which
557    *                  contents were written in a range-compacted format.
558    *                  (mhdf_SET_RANGE_BIT).  Should be null for parents/children.
559    *\param null_stripped Will be populated with handles of sets for which
560    *                  invalid or null handles were stripped from the contents
561    *                  list.  This is only done for unordered sets.  This argument
562    *                  should be null if writing parents/children because those
563    *                  lists are always ordered.
564    *\param set_sizes  Will be populated with the length of the data written
565    *                  for those sets for which the handles were added to
566    *                  either \c ranged or \c null_stripped.  Values are
567    *                  in handle order.
568    */
569   ErrorCode write_set_data( const WriteUtilIface::EntityListType which_data,
570                             const hid_t handle,
571                             IODebugTrack& track,
572                             Range* ranged = 0,
573                             Range* null_stripped = 0,
574                             std::vector<long>* set_sizes = 0);
575 
576   /** Write adjacency info for passed set of elements
577    *
578    * Note: Must have written element connectivity so elements
579    *       have IDs assigned.
580    */
581   ErrorCode write_adjacencies( const ExportSet& export_set );
582 
583   /** Write tag information and data.
584    *
585    * Note: Must have already written nodes, elem connectivity and
586    *       sets so that entities have IDs assigned.
587    */
588 
589   //! Write tag for all entities.
590   ErrorCode write_tag( const TagDesc& tag_data,
591                        double* times );
592 
593   //! Get element connectivity
594   ErrorCode get_connectivity( Range::const_iterator begin,
595                               Range::const_iterator end,
596                               int nodes_per_element,
597                               wid_t* id_data_out );
598 
599   //! Get size data for tag
600   //!\param tag       MOAB tag ID
601   //!\param moab_type Output: DataType for tag
602   //!\param num_bytes Output: MOAB tag size (bits for bit tags).
603   //!                         MB_VARIABLE_LENGTH for variable-length tags.
604   //!\param elem_size Output: Size of of the base data type of the
605   //!                         tag data (e.g. sizeof(double) if
606   //!                         moab_type == MB_TYPE_DOUBLE).
607   //!                         One for bit and opaque tags.
608   //!\param array_size Output: The number of valeus of size elem_size
609   //!                          for each tag.  Always 1 for opaque data.
610   //!                          Nubmer of bits for bit tags.
611   //!\param file_type Output: mhdf type enumeration
612   //!\param hdf_type  Output: Handle to HDF5 type object.  Caller is
613   //!                         responsible for releasing this object
614   //!                         (calling H5Tclose).
615   ErrorCode get_tag_size( Tag tag,
616                           DataType& moab_type,
617                           int& num_bytes,
618                           int& elem_size,
619                           int& file_size,
620                           mhdf_TagDataType& file_type,
621                           hid_t& hdf_type );
622 
623   //! Write ID table for sparse tag
624   ErrorCode write_sparse_ids( const TagDesc& tag_data,
625                               const Range& range,
626                               hid_t table_handle,
627                               size_t table_size,
628                               const char* name = 0 );
629 
630   //! Write fixed-length tag data in sparse format
631   ErrorCode write_sparse_tag( const TagDesc& tag_data,
632                               const std::string& tag_name,
633                               DataType tag_data_type,
634                               hid_t hdf5_data_type,
635                               int hdf5_type_size );
636 
637   //! Write end index data_set for a variable-length tag
638   ErrorCode write_var_len_indices( const TagDesc& tag_data,
639                                    const Range& range,
640                                    hid_t idx_table,
641                                    size_t table_size,
642                                    int type_size,
643                                    const char* name = 0 );
644 
645   //! Write tag value data_set for a variable-length tag
646   ErrorCode write_var_len_data( const TagDesc& tag_data,
647                                 const Range& range,
648                                 hid_t table,
649                                 size_t table_size,
650                                 bool handle_tag,
651                                 hid_t hdf_type,
652                                 int type_size,
653                                 const char* name = 0 );
654 
655   //! Write varialbe-length tag data
656   ErrorCode write_var_len_tag( const TagDesc& tag_info,
657                                const std::string& tag_name,
658                                DataType tag_data_type,
659                                hid_t hdf5_type,
660                                int hdf5_type_size );
661 
662   //! Write dense-formatted tag data
663   ErrorCode write_dense_tag( const TagDesc& tag_data,
664                              const ExportSet& elem_data,
665                              const std::string& tag_name,
666                              DataType tag_data_type,
667                              hid_t hdf5_data_type,
668                              int hdf5_type_size );
669 
670   //! Write data for fixed-size tag
671   ErrorCode write_tag_values( Tag tag_id,
672                               hid_t data_table,
673                               unsigned long data_offset,
674                               const Range& range,
675                               DataType tag_data_type,
676                               hid_t hdf5_data_type,
677                               int hdf5_type_size,
678                               unsigned long max_num_ents,
679                               IODebugTrack& debug_track );
680 
681 protected:
682 
683   enum TimingValues {
684        TOTAL_TIME = 0,
685          GATHER_TIME,
686          CREATE_TIME,
687            CREATE_NODE_TIME,
688            NEGOTIATE_TYPES_TIME,
689            CREATE_ELEM_TIME,
690            FILEID_EXCHANGE_TIME,
691            CREATE_ADJ_TIME,
692            CREATE_SET_TIME,
693              SHARED_SET_IDS,
694              SHARED_SET_CONTENTS,
695              SET_OFFSET_TIME,
696            CREATE_TAG_TIME,
697          COORD_TIME,
698          CONN_TIME,
699          SET_TIME,
700            SET_META,
701            SET_CONTENT,
702            SET_PARENT,
703            SET_CHILD,
704          ADJ_TIME,
705          TAG_TIME,
706            DENSE_TAG_TIME,
707            SPARSE_TAG_TIME,
708            VARLEN_TAG_TIME,
709          NUM_TIMES };
710 
711 
712   virtual void print_times( const double times[NUM_TIMES] ) const;
713 };
714 
715 } // namespace moab
716 
717 #endif
718