MOAB: Mesh Oriented datABase  (version 5.5.0)
mhdf_parallel.c
Go to the documentation of this file.
1 /* Test mhdf library on top of Parallel HDF5.
2  *
3  * Call mhdf API similar to how WriteHDF5Parallel does,
4  * but avoid any of our own parallel communiation.
5  *
6  * Output will be composed of:
7  * a 1-D array of hexes, one per processor, where the first
8  * process writes all the nodes for its hex and every other
9  * process writes the right 4 nodes of its hex.
10  * a set containing all of the hexes
11  * a set containing all of the nodes
12  * a set per processor containing all of the entities written by that proc
13  * a tag specifying the process ID that wrote each entity.
14  */
15 
16 #include "mhdf.h"
17 
18 #include <time.h>
19 #include <stdlib.h>
20 #include <stdio.h>
21 #include <assert.h>
22 
23 #include <H5Ppublic.h>
24 #include <H5Tpublic.h>
25 #include <H5Epublic.h>
26 #include <H5FDmpi.h>
27 #include <H5FDmpio.h>
28 
29 #ifdef H5_HAVE_PARALLEL
30 
31 const char filename[] = "mhdf_ll.h5m";
32 const char proc_tag_name[] = "proc_id";
33 const char elem_handle[] = "Hex";
34 
35 int RANK;
36 int NUM_PROC;
37 
38 #define CHECK( A ) check( &( A ), __LINE__ )
39 void check( mhdf_Status* status, int line )
40 {
41  if( mhdf_isError( status ) )
42  {
43  fprintf( stderr, "%d: ERROR at line %d: \"%s\"\n", RANK, line, mhdf_message( status ) );
44  abort();
45  }
46 }
47 
48 typedef int handle_id_t;
49 #define id_type H5T_NATIVE_INT
50 
51 /* create file layout in serial */
52 void create_file()
53 {
54  const char* elem_types[1] = { elem_handle };
55  const int num_elem_types = sizeof( elem_types ) / sizeof( elem_types[0] );
56  const int num_nodes = 4 + 4 * NUM_PROC;
57  const int num_hexes = NUM_PROC;
58  const int num_sets = 2 + NUM_PROC;
59  const int set_data_size = num_hexes + num_nodes + num_hexes + num_nodes;
60  const int num_tag_vals = num_hexes + num_nodes + num_sets - 2;
61  const char* history[] = { __FILE__, NULL };
62  const int history_size = sizeof( history ) / sizeof( history[0] );
63  int default_tag_val = NUM_PROC;
64  hid_t data, tagdata[2];
65  long junk;
66 
67  mhdf_Status status;
68  mhdf_FileHandle file;
69 
70  time_t t;
71  t = time( NULL );
72  history[1] = ctime( &t );
73 
74  file = mhdf_createFile( filename, 1, elem_types, num_elem_types, id_type, &status );
75  CHECK( status );
76 
77  mhdf_writeHistory( file, history, history_size, &status );
78  CHECK( status );
79 
80  data = mhdf_createNodeCoords( file, 3, num_nodes, &junk, &status );
81  CHECK( status );
82  mhdf_closeData( file, data, &status );
83  CHECK( status );
84 
85  mhdf_addElement( file, elem_types[0], 0, &status );
86  CHECK( status );
87  data = mhdf_createConnectivity( file, elem_handle, 8, num_hexes, &junk, &status );
88  CHECK( status );
89  mhdf_closeData( file, data, &status );
90  CHECK( status );
91 
92  data = mhdf_createSetMeta( file, num_sets, &junk, &status );
93  CHECK( status );
94  mhdf_closeData( file, data, &status );
95  CHECK( status );
96 
97  data = mhdf_createSetData( file, set_data_size, &status );
98  CHECK( status );
99  mhdf_closeData( file, data, &status );
100  CHECK( status );
101 
102  mhdf_createTag( file, proc_tag_name, mhdf_INTEGER, 1, mhdf_DENSE_TYPE, &default_tag_val, &default_tag_val, 0, 0,
103  &status );
104  CHECK( status );
105 
106  mhdf_createSparseTagData( file, proc_tag_name, num_tag_vals, tagdata, &status );
107  CHECK( status );
108  mhdf_closeData( file, tagdata[0], &status );
109  CHECK( status );
110  mhdf_closeData( file, tagdata[1], &status );
111  CHECK( status );
112 
113  mhdf_closeFile( file, &status );
114  CHECK( status );
115 }
116 
117 void write_file_data()
118 {
119  const int total_num_nodes = 4 + 4 * NUM_PROC;
120  const int total_num_hexes = NUM_PROC;
121  long first_node, first_elem, first_set, count, ntag;
122  unsigned long ucount;
123  mhdf_index_t set_desc[4] = { 0, -1, -1, 0 };
124  hid_t handle, handles[2];
125  mhdf_Status status;
126  mhdf_FileHandle file;
127  int num_node, offset, dim;
128  double coords[3 * 8] = { 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 1.0, 1.0, 0.0, 0.0, 1.0,
129  0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 1.0, 1.0, 0.0, 0.0, 1.0 };
130  handle_id_t list[10];
131  int i, tagdata[10];
132  for( i = 0; i < 10; i++ )
133  tagdata[i] = RANK;
134 
135  /* open file */
136  handle = H5Pcreate( H5P_FILE_ACCESS );
137  H5Pset_fapl_mpio( handle, MPI_COMM_WORLD, MPI_INFO_NULL );
138  file = mhdf_openFileWithOpt( filename, 1, &ucount, id_type, handle, &status );
139  CHECK( status );
140  H5Pclose( handle );
141 
142  /* write node coordinates */
143  if( 0 == RANK )
144  {
145  num_node = 8;
146  offset = 0;
147  coords[12] = coords[15] = coords[18] = coords[21] = 1.0;
148  }
149  else
150  {
151  num_node = 4;
152  offset = 4 + 4 * RANK;
153  coords[0] = coords[3] = coords[6] = coords[9] = 1.0 + RANK;
154  }
155  handle = mhdf_openNodeCoords( file, &count, &dim, &first_node, &status );
156  CHECK( status );
157  assert( count == total_num_nodes );
158  assert( dim == 3 );
159  mhdf_writeNodeCoords( handle, offset, num_node, coords, &status );
160  CHECK( status );
161  mhdf_closeData( file, handle, &status );
162  CHECK( status );
163 
164  /* write hex connectivity */
165  for( i = 0; i < 8; ++i )
166  list[i] = 4 * RANK + i + first_node;
167  handle = mhdf_openConnectivity( file, elem_handle, &dim, &count, &first_elem, &status );
168  CHECK( status );
169  assert( count == total_num_hexes );
170  assert( 8 == dim );
171  mhdf_writeConnectivity( handle, RANK, 1, id_type, list, &status );
172  CHECK( status );
173  mhdf_closeData( file, handle, &status );
174  CHECK( status );
175 
176  /* write set descriptions */
177  handle = mhdf_openSetMeta( file, &count, &first_set, &status );
178  CHECK( status );
179  assert( count == 2 + NUM_PROC );
180 
181  /* write set descriptions for per-proc sets */
182  set_desc[0] = total_num_nodes + total_num_hexes + 9 + 5 * RANK - 1;
183  mhdf_writeSetMeta( handle, 2 + RANK, 1, MHDF_INDEX_TYPE, set_desc, &status );
184  CHECK( status );
185 
186  /* write set descriptions for multi-proc sets */
187  if( 0 == RANK )
188  {
189  set_desc[0] = total_num_nodes - 1;
190  mhdf_writeSetMeta( handle, 0, 1, MHDF_INDEX_TYPE, set_desc, &status );
191  CHECK( status );
192  set_desc[0] += total_num_hexes;
193  mhdf_writeSetMeta( handle, 1, 1, MHDF_INDEX_TYPE, set_desc, &status );
194  CHECK( status );
195  }
196 
197  mhdf_closeData( file, handle, &status );
198  CHECK( status );
199 
200  /* write set contents */
201 
202  handle = mhdf_openSetData( file, &count, &status );
203  CHECK( status );
204  assert( 2 * total_num_nodes + 2 * total_num_hexes == count );
205 
206  /* write per-proc sets */
207  if( 0 == RANK )
208  {
209  count = 9;
210  offset = total_num_nodes + total_num_hexes;
211  for( i = 0; i < 8; ++i )
212  list[i] = first_node + i;
213  list[8] = first_elem;
214  }
215  else
216  {
217  count = 5;
218  offset = total_num_nodes + total_num_hexes + 4 + 5 * RANK;
219  for( i = 0; i < 4; ++i )
220  list[i] = 4 + 4 * RANK + first_node + i;
221  list[4] = RANK + first_elem;
222  }
223  mhdf_writeSetData( handle, offset, count, id_type, list, &status );
224  CHECK( status );
225 
226  /* write multi-proc sets */
227  /* write nodes */
228  offset = ( 0 == RANK ) ? 0 : 4 + 4 * RANK;
229  mhdf_writeSetData( handle, offset, count - 1, id_type, list, &status );
230  CHECK( status );
231  /* write hex */
232  offset = total_num_nodes + RANK;
233  mhdf_writeSetData( handle, offset, 1, id_type, list + count - 1, &status );
234  CHECK( status );
235 
236  mhdf_closeData( file, handle, &status );
237  CHECK( status );
238 
239  /* write tag data */
240  offset = ( 0 == RANK ) ? 0 : 4 + 4 * RANK;
241  offset += 2 * RANK;
242  list[count++] = first_set + 2 + RANK;
243  mhdf_openSparseTagData( file, proc_tag_name, &ntag, &ntag, handles, &status );
244  CHECK( status );
245  assert( ntag == total_num_nodes + total_num_hexes + NUM_PROC );
246  mhdf_writeSparseTagEntities( handles[0], offset, count, id_type, list, &status );
247  CHECK( status );
248  mhdf_writeTagValues( handles[1], offset, count, H5T_NATIVE_INT, tagdata, &status );
249  CHECK( status );
250  mhdf_closeData( file, handles[0], &status );
251  CHECK( status );
252  mhdf_closeData( file, handles[1], &status );
253  CHECK( status );
254 
255  /* done */
256  mhdf_closeFile( file, &status );
257  CHECK( status );
258 }
259 
260 /* Define a dummy error handler to register with HDF5.
261  * This function doesn't do anything except pass the
262  * error on to the default handler that would have been
263  * called anyway. It's only purpose is to provide a
264  * spot to set a break point so we can figure out where
265  * (in our code) that we made an invalid call into HDF5
266  */
267 #if defined( H5E_auto_t_vers ) && H5E_auto_t_vers > 1
268 herr_t ( *default_handler )( hid_t, void* );
269 static herr_t handle_hdf5_error( hid_t stack, void* data )
270 #else
271 herr_t ( *default_handler )( void* );
272 static herr_t handle_hdf5_error( void* data )
273 #endif
274 {
275  herr_t result = 0;
276  if( default_handler )
277 #if defined( H5E_auto_t_vers ) && H5E_auto_t_vers > 1
278  result = (default_handler)( stack, data );
279 #else
280  result = (default_handler)( data );
281 #endif
282  assert( 0 );
283  return result;
284 }
285 
286 #endif /* #ifdef H5_HAVE_PARALLEL */
287 
288 int main( int argc, char* argv[] )
289 {
290 #ifdef H5_HAVE_PARALLEL
291  int rval;
292  void* data;
293  herr_t err;
294 
295 #if defined( H5Eget_auto_vers ) && H5Eget_auto_vers > 1
296  err = H5Eget_auto( H5E_DEFAULT, &default_handler, &data );
297 #else
298  err = H5Eget_auto( &default_handler, &data );
299 #endif
300  if( err >= 0 )
301  {
302 #if defined( H5Eset_auto_vers ) && H5Eset_auto_vers > 1
303  H5Eset_auto( H5E_DEFAULT, &handle_hdf5_error, data );
304 #else
305  H5Eset_auto( &handle_hdf5_error, data );
306 #endif
307  }
308 
309  rval = MPI_Init( &argc, &argv );
310  if( rval ) return rval;
311  rval = MPI_Comm_rank( MPI_COMM_WORLD, &RANK );
312  if( rval ) return rval;
313  rval = MPI_Comm_size( MPI_COMM_WORLD, &NUM_PROC );
314  if( rval ) return rval;
315 
316  if( RANK == 0 ) create_file();
317 
318  /* Wait for rank 0 to finish creating the file, otherwise rank 1 may find it to be invalid */
319  rval = MPI_Barrier( MPI_COMM_WORLD );
320  if( rval ) return rval;
321 
322  write_file_data();
323 
324  MPI_Finalize();
325 #endif
326  return 0;
327 }