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 #ifndef READ_HDF5_HPP
17 #define READ_HDF5_HPP
18 
19 #include <stdlib.h>
20 #include <list>
21 #include "mhdf.h"
22 #include "moab/Forward.hpp"
23 #include "moab/ReadUtilIface.hpp"
24 #include "moab/Range.hpp"
25 #include "moab/ReaderIface.hpp"
26 #include "moab/RangeMap.hpp"
27 #include "DebugOutput.hpp"
28 #include "HDF5Common.hpp"
29 
30 #ifdef MOAB_HAVE_MPI
31 # include <moab_mpi.h>
32 #endif
33 
34 namespace moab {
35 
36 class ParallelComm;
37 class ReadHDF5Dataset;
38 class CpuTimer;
39 
40 /**
41  * \brief  Read mesh from MOAB HDF5 (.h5m) file.
42  * \author Jason Kraftcheck
43  * \date   18 April 2004
44  */
45 class ReadHDF5 : public ReaderIface
46 {
47 public:
48 
49 #ifdef MOAB_HAVE_MPI
50   typedef MPI_Comm Comm;
51 #else
52   typedef int Comm;
53 #endif
54 
55   static ReaderIface* factory(Interface*);
56 
57   ReadHDF5(Interface* iface);
58 
59   virtual ~ReadHDF5();
60 
61   /** Export specified meshsets to file
62    * \param filename     The filename to export. Must end in <em>.mhdf</em>
63    * \param export_sets  Array of handles to sets to export, or NULL to export all.
64    * \param export_set_count Length of <code>export_sets</code> array.
65    */
66   ErrorCode load_file(const char* file_name,
67                       const EntityHandle* file_set,
68                       const FileOptions& opts,
69                       const SubsetList* subset_list = 0,
70                       const Tag* file_id_tag = 0);
71 
72   ErrorCode read_tag_values(const char* file_name,
73                             const char* tag_name,
74                             const FileOptions& opts,
75                             std::vector<int>& tag_values_out,
76                             const SubsetList* subset_list = 0);
moab() const77   Interface* moab() const
78     { return iFace; }
79 
80   //! Store old HDF5 error handling function
81   struct HDF5ErrorHandler {
82     HDF5_Error_Func_Type func;
83     void* data;
84   };
85 
86 protected:
87 
88   ErrorCode load_file_impl(const FileOptions& opts);
89 
90   ErrorCode load_file_partial(const ReaderIface::IDTag* subset_list,
91                               int subset_list_length,
92                               int num_parts,
93                               int part_number,
94                               const FileOptions& opts);
95 
96   ErrorCode read_tag_values_all(int tag_index, std::vector<int>& results);
97   ErrorCode read_tag_values_partial(int tag_index, const Range& file_ids,
98                                     std::vector<int>& results);
99 
100   enum ReadTimingValues {
101     TOTAL_TIME = 0,
102     SET_META_TIME,
103     SUBSET_IDS_TIME,
104     GET_PARTITION_TIME,
105     GET_SET_IDS_TIME,
106     GET_SET_CONTENTS_TIME,
107     GET_POLYHEDRA_TIME,
108     GET_ELEMENTS_TIME,
109     GET_NODES_TIME,
110     GET_NODEADJ_TIME,
111     GET_SIDEELEM_TIME,
112     UPDATECONN_TIME,
113     ADJACENCY_TIME,
114     DELETE_NON_SIDEELEM_TIME,
115     READ_SET_IDS_RECURS_TIME,
116     FIND_SETS_CONTAINING_TIME,
117     READ_SETS_TIME,
118     READ_TAGS_TIME,
119     STORE_FILE_IDS_TIME,
120     READ_QA_TIME,
121     NUM_TIMES
122   };
123 
124 
125   void print_times( );
126 
127 private:
128   ErrorCode init();
129 
is_error(mhdf_Status & status)130   inline int is_error(mhdf_Status& status) {
131     int i;
132     if ((i = mhdf_isError(&status)))
133       MB_SET_ERR_CONT(mhdf_message(&status));
134     return i;
135   }
136 
137   //! The size of the data buffer (<code>dataBuffer</code>).
138   int bufferSize;
139   //! A memory buffer to use for all I/O operations.
140   char* dataBuffer;
141 
142   //! Interface pointer passed to constructor
143   Interface* iFace;
144 
145   //! The file handle from the mhdf library
146   mhdf_FileHandle filePtr;
147 
148   //! File summary
149   mhdf_FileDesc* fileInfo;
150 
151   //! Map from File ID to MOAB handle
152   typedef RangeMap<long, EntityHandle> IDMap;
153   IDMap idMap;
154 
155   //! Cache pointer to read util
156   ReadUtilIface* readUtil;
157 
158   //! The type of an EntityHandle
159   hid_t handleType;
160 
161   //! List of connectivity arrays for which conversion from file ID
162   //! to handle was deferred until later.
163   struct IDConnectivity {
164     EntityHandle handle; // Start_handle
165     size_t count; // Num entities
166     int nodes_per_elem; // Per-element connectivity length
167     EntityHandle* array; // Connectivity array
168   };
169   //! List of connectivity arrays for which conversion from file ID
170   //! to handle was deferred until later.
171   std::vector<IDConnectivity> idConnectivityList;
172 
173   //! Read/write property handle
174   //! indepIO -> independent IO during true parallel read
175   //! collIO  -> collective IO during true parallel read
176   //! Both are H5P_DEFAULT for serial IO and collective
177   //! when reading the entire file on all processors.
178   hid_t indepIO, collIO;
179 
180   ParallelComm* myPcomm;
181 
182   //! Use IODebugTrack instances to verify reads.
183   //! Enable with the DEBUG_OVERLAPS option.
184   bool debugTrack;
185   //! Debug output. Verbosity controlled with DEBUG_FORMAT option.
186   DebugOutput dbgOut;
187   //! Doing true parallel read (PARALLEL=READ_PART)
188   bool nativeParallel;
189   //! MPI_Comm value (unused if \c !nativeParallel)
190   Comm* mpiComm;
191 
192   //! Flags for some behavior that can be changed through
193   //! reader options
194   bool blockedCoordinateIO;
195   bool bcastSummary;
196   bool bcastDuplicateReads;
197 
198   //! Store old HDF5 error handling function
199   HDF5ErrorHandler errorHandler;
200 
201   long (*setMeta)[4];
202 
203   double _times[NUM_TIMES];
204   CpuTimer * timer;
205   bool cputime;
206   ErrorCode read_all_set_meta();
207 
208   ErrorCode set_up_read(const char* file_name, const FileOptions& opts);
209   ErrorCode clean_up_read(const FileOptions& opts);
210 
211   //! Given a list of tags and values, get the file ids for the
212   //! corresponding entities in the file.
213   ErrorCode get_subset_ids(const ReaderIface::IDTag* subset_list,
214                            int subset_list_length,
215                            Range& file_ids_out);
216 
217   /**\brief Remove all but the specified fraction of sets from the passed range
218    *
219    * Select a subset of the gathered set of file ids to read in based on
220    * communicator size and rank.
221    *\param tmp_file_ids As input: sets to be read on all procs.
222    *                    As output: sets to read on this proc.
223    *\param num_parts    communicator size
224    *\param part_number  communicator rank
225    */
226   ErrorCode get_partition(Range& tmp_file_ids, int num_parts, int part_number);
227 
228   ErrorCode read_nodes(const Range& node_file_ids);
229 
230   // Read elements in fileInfo->elems[index]
231   ErrorCode read_elems(int index);
232 
233   // Read subset of elements in fileInfo->elems[index]
234   ErrorCode read_elems(int index, const Range& file_ids, Range* node_ids = 0);
235 
236   //! Read element connectivity.
237   //!
238   //!\param node_ids  If this is non-null, the union of the connectivity list
239   //!                 for all elements is passed back as FILE IDS in this
240   //!                 range AND the connectivity list is left in the form
241   //!                 of file IDs (NOT NODE HANDLES).
242   ErrorCode read_elems(const mhdf_ElemDesc& elems, const Range& file_ids, Range* node_ids = 0);
243 
244   //! Update connectivity data for all element groups for which read_elems
245   //! was called with a non-null \c node_ids argument.
246   ErrorCode update_connectivity();
247 
248   // Read connectivity data for a list of element file ids.
249   // passing back the file IDs for the element connectivity
250   // w/out actually creating any elements in MOAB.
251   ErrorCode read_elems(int index, const Range& element_file_ids, Range& node_file_ids);
252 
253   // Scan all elements in group. For each element for which all vertices
254   // are contained in idMap (class member), read the element. All new
255   // elements are added to idMap.
256   //
257   // NOTE: Collective IO calls in parallel.
258   ErrorCode read_node_adj_elems(const mhdf_ElemDesc& group,
259                                 Range* read_entities = 0);
260 
261   // Scan all elements in specified file table. For each element for
262   // which all vertices are contained in idMap (class member), read the
263   // element. All new elements are added to idMap.
264   //
265   // NOTE: Collective IO calls in parallel.
266   ErrorCode read_node_adj_elems(const mhdf_ElemDesc& group,
267                                 hid_t connectivity_handle,
268                                 Range* read_entities = 0);
269 
270   //! Read poly(gons|hedra)
271   ErrorCode read_poly(const mhdf_ElemDesc& elems, const Range& file_ids);
272 
273   //! Clean up elements that were a) read because we had read all of the
274   //! nodes and b) weren't actually sides of the top-dimension elements
275   //! we were trying to read.
276   ErrorCode delete_non_side_elements(const Range& side_ents);
277 
278   //! Read sets
279   ErrorCode read_sets(const Range& set_file_ids);
280 
281   ErrorCode read_adjacencies(hid_t adjacency_table,
282                              long table_length);
283 
284   //! Create tag and read all data.
285   ErrorCode read_tag(int index);
286 
287   //! Create new tag or verify type matches existing tag
288   ErrorCode create_tag(const mhdf_TagDesc& info, Tag& handle, hid_t& type);
289 
290   //! Read dense tag for all entities
291   ErrorCode read_dense_tag(Tag tag_handle,
292                            const char* ent_name,
293                            hid_t hdf_read_type,
294                            hid_t data_table,
295                            long start_id,
296                            long count);
297 
298   //! Read sparse tag for all entities.
299   ErrorCode read_sparse_tag(Tag tag_handle,
300                             hid_t hdf_read_type,
301                             hid_t ent_table,
302                             hid_t val_table,
303                             long num_entities);
304 
305   //! Read variable-length tag for all entities.
306   ErrorCode read_var_len_tag(Tag tag_handle,
307                              hid_t hdf_read_type,
308                              hid_t ent_table,
309                              hid_t val_table,
310                              hid_t off_table,
311                              long num_entities,
312                              long num_values);
313 
314   /**\brief Read index table for sparse tag.
315    *
316    * Read ID table for a sparse or variable-length tag, returning
317    * the handles and offsets within the table for each file ID
318    * that corresponds to an entity we've read from the file (an
319    * entity that is in \c idMap).
320    *
321    * \param id_table     The MOAB handle for the tag
322    * \param start_offset Some non-zero value because ranges (in this case
323    *                     the offset_range) cannot contain zeros.
324    * \param offset_range Output: The offsets in the id table for which IDs
325    *                     that occur in \c idMap were found. All values
326    *                     are increased by \c start_offset to avoid
327    *                     putting zeros in the range.
328    * \param handle_range Output: For each valid ID read from the table,
329    *                     the corresponding entity handle. Note: if
330    *                     the IDs did not occur in handle order, then
331    *                     this will be empty. Use \c handle_vect instead.
332    * \param handle_vect  Output: For each valid ID read from the table,
333    *                     the corresponding entity handle. Note: if
334    *                     the IDs occurred in handle order, then
335    *                     this will be empty. Use \c handle_range instead.
336    */
337   ErrorCode read_sparse_tag_indices(const char* name,
338                                     hid_t id_table,
339                                     EntityHandle start_offset,
340                                     Range& offset_range,
341                                     Range& handle_range,
342                                     std::vector<EntityHandle>& handle_vect);
343 
344   ErrorCode read_qa(EntityHandle file_set);
345 
346 public:
347   ErrorCode convert_id_to_handle(EntityHandle* in_out_array,
348                                  size_t array_length);
349 
convert_id_to_handle(EntityHandle * in_out_array,size_t array_length,size_t & array_length_out) const350   void convert_id_to_handle(EntityHandle* in_out_array,
351                             size_t array_length,
352                             size_t& array_length_out) const
353     { return convert_id_to_handle(in_out_array, array_length, array_length_out, idMap); }
354 
355   ErrorCode convert_range_to_handle(const EntityHandle* ranges,
356                                     size_t num_ranges,
357                                     Range& merge);
358 
359   static
360   void convert_id_to_handle(EntityHandle* in_out_array,
361                             size_t array_length,
362                             const RangeMap<long, EntityHandle>& id_map);
363 
364   static
365   void convert_id_to_handle(EntityHandle* in_out_array,
366                             size_t array_length,
367                             size_t& array_length_out,
368                             const RangeMap<long, EntityHandle>& id_map);
369 
370   static
371   void convert_range_to_handle(const EntityHandle* ranges,
372                                size_t num_ranges,
373                                const RangeMap<long, EntityHandle>& id_map,
374                                Range& merge);
375 
376   ErrorCode insert_in_id_map(const Range& file_ids,
377                              EntityHandle start_id);
378 
379   ErrorCode insert_in_id_map(long file_id, EntityHandle handle);
380 
381 private:
382 
383   /**\brief Search for entities with specified tag values
384    *
385    *\NOTE For parallel reads, this function does collective IO.
386    *
387    *\param tag_index  Index into info->tags specifying which tag to search.
388    *\param sorted_values  List of tag values to check for, in ascending sorted
389    *                  order.
390    *\param file_ids_out  File IDs for entities with specified tag values.
391    */
392   ErrorCode search_tag_values(int tag_index,
393                               const std::vector<int>& sorted_values,
394                               Range& file_ids_out,
395                               bool sets_only = false);
396 
397   /**\brief Search for entities with specified tag
398    *
399    *\NOTE For parallel reads, this function does collective IO.
400    *
401    *\param tag_index  Index into info->tags specifying which tag to search.
402    *\param file_ids_out  File IDs for entities with specified tag values.
403    */
404   ErrorCode get_tagged_entities(int tag_index, Range& file_ids_out);
405 
406   /**\brief Search a table of tag data for a specified set of values.
407    *
408    * Search a table of tag values, returning the indices into the table
409    * at which matches were found.
410    *\NOTE For parallel reads, this function does collective IO.
411    *
412    *\param info       Summary of data contained in file.
413    *\param tag_table     HDF5/mhdf handle for tag values
414    *\param table_size    Number of values in table
415    *\param sorted_values Sorted list of values to search for.
416    *\param value_indices Output: Offsets into the table of data at which
417    *                       matching values were found.
418    */
419   ErrorCode search_tag_values(hid_t tag_table,
420                               unsigned long table_size,
421                               const std::vector<int>& sorted_values,
422                               std::vector<EntityHandle>& value_indices);
423 
424   /**\brief Get the file IDs for nodes and elements contained in sets.
425    *
426    * Read the contents for the specified sets and return the file IDs
427    * of all nodes and elements contained within those sets.
428    *\param sets       Container of file IDs designating entity sets.
429    *\param file_ids   Output: File IDs of entities contained in sets.
430    */
431   ErrorCode get_set_contents(const Range& sets, Range& file_ids);
432 
433   /** Given a list of file IDs for entity sets, find all contained
434    *  or child sets (at any depth) and append them to the Range
435    *  of file IDs.
436    */
437   ErrorCode read_set_ids_recursive(Range& sets_in_out,
438                                    bool containted_sets,
439                                    bool child_sets);
440 
441   /** Find all sets containing one or more entities read from the file
442    *  and added to idMap
443    */
444   ErrorCode find_sets_containing(Range& sets_out, bool read_set_containing_parents);
445 
446   /**\brief Read sets from file into MOAB for partial read of file.
447    *
448    * Given the file IDs for entity sets (sets_in) and elements and
449    * nodes (id_map), read in all sets containing any of the elements
450    * or nodes and all sets that are (recursively) children of any
451    * other set to be read (those in sets_in or those containing any
452    * already-read element or node.)
453    *\param sets_in    File IDs for sets to read (unconditionally)
454    */
455   ErrorCode read_sets_partial(const Range& sets_in);
456 
457   /** Find file IDs of sets containing any entities in the passed id_map */
458   ErrorCode find_sets_containing(hid_t content_handle,
459                                  hid_t content_type,
460                                  long content_len,
461                                  bool read_set_containing_parents,
462                                  Range& file_ids);
463 
464 public:
465   enum SetMode {CONTENT = 0, CHILD = 1, PARENT = 2};
466 
467 private:
468   // Read the set data specified by by which_data.
469   // Update set data in MOAB according to which_data,
470   // unless file_ids_out is non-null. If file_ids_out is
471   // non-null, then return the file IDs in the passed
472   // range and do not make any changes to MOAB data
473   // structures. If file_ids_out is NULL, then set_start_handle
474   // is ignored.
475   ErrorCode read_set_data(const Range& set_file_ids,
476                           EntityHandle set_start_handle,
477                           ReadHDF5Dataset& set_data_set,
478                           SetMode which_data,
479                           Range* file_ids_out = 0);
480 
481   /**\brief Store file IDS in tag values
482    *
483    * Copy file ID from IDMap for each entity read from file
484    * into a tag value on the entity.
485    */
486   ErrorCode store_file_ids(Tag tag);
487 
488   /**\brief Store sets file IDS in a custom tag
489    *
490    * Copy sets file ID from IDMap for each set read from file
491    *  into a custom tag on the sets
492    *  needed now for the VisIt plugin, to identify sets read
493    */
494   ErrorCode store_sets_file_ids();
495 
496   /**\brief Find index in mhdf_FileDesc* fileInfo for specified tag name
497    *
498    * Given a tag name, find its index in fileInfo and verify that
499    * each tag value is a single integer.
500    */
501   ErrorCode find_int_tag(const char* name, int& index_out);
502 
503   void debug_barrier_line(int lineno);
504 };
505 
506 } // namespace moab
507 
508 #endif
509