diff --git a/test/H5Fed/Makefile.am b/test/H5Fed/Makefile.am index 4f69a66..e74fe22 100644 --- a/test/H5Fed/Makefile.am +++ b/test/H5Fed/Makefile.am @@ -37,11 +37,11 @@ INC = $(HDFINC) $(MPIINC) $(H5INC) # What to build... make install will place these files in the $(prefix)/bin directory. bin_PROGRAMS = \ - write_tetmesh + write_tetmesh \ + read_tetmesh \ + write_trianglemesh \ + read_trianglemesh -# read_tetmesh \ -# write_trianglemesh \ -# read_trianglemesh \ # map_tet2globalid \ # map_triangle2globalid diff --git a/test/H5Fed/read_tetmesh.c b/test/H5Fed/read_tetmesh.c index 82dc590..cd83113 100644 --- a/test/H5Fed/read_tetmesh.c +++ b/test/H5Fed/read_tetmesh.c @@ -25,22 +25,27 @@ typedef struct tet tet_t; h5_err_t read_vertices ( - h5_file * f + h5_file_t * f ) { h5_id_t id, local_id; h5_float64_t P[3]; h5_size_t real_num = 0; h5_size_t num = H5FedGetNumVerticesTotal ( f ); - printf ( " Number of vertices on level: %d\n", num ); + printf ( " Number of vertices on level: %lld\n", num ); + + h5_err_t h5err = H5FedStartTraverseVertices ( f ); + if ( h5err < 0 ) return h5err; while ( (real_num < num) && ((local_id = H5FedTraverseVertices ( f, &id, P )) >= 0) ) { - printf ( " Vertex[%d]: local id: %d, coords: %f %f %f \n", + printf ( " Vertex[%lld]: local id: %lld, coords: %f %f %f \n", id, local_id, P[0], P[1], P[2] ); real_num++; } + H5FedEndTraverseVertices ( f ); + if ( real_num != num ) { - fprintf ( stderr, "!!! Got %d vertices, but expected %d.\n", + fprintf ( stderr, "!!! Got %lld vertices, but expected %lld.\n", real_num, num ); return -1; } @@ -49,21 +54,29 @@ read_vertices ( h5_err_t read_tets ( - h5_file * f + h5_file_t * f ) { h5_id_t id, local_id, parent_id, vids[4]; h5_size_t real_num = 0; - h5_size_t num = H5FedGetNumTetrahedraTotal ( f ); - printf ( " Number of tetrahedra on level: %d\n", num ); + h5_size_t num = H5FedGetNumElementsTotal ( f ); + printf ( " Number of tetrahedra on level: %lld\n", num ); + + h5_err_t h5err = H5FedStartTraverseElements ( f ); + if ( h5err < 0 ) return h5err; + while ( (real_num < num) && - ((local_id = H5FedTraverseTetrahedra ( f, &id, &parent_id, vids )) >= 0) ) { - printf ( " Tet[%d]: local id: %d, parent id: %d, vids: %d %d %d %d\n", - id, local_id, parent_id, vids[0], vids[1], vids[2], vids[3] ); + ((local_id = H5FedTraverseElements ( + f, &id, &parent_id, vids )) >= 0) ) { + printf ( " Tet[%lld]: local id: %lld, parent id: %lld," + " vids: %lld %lld %lld %lld\n", + id, local_id, parent_id, + vids[0], vids[1], vids[2], vids[3] ); real_num++; } + H5FedEndTraverseElements ( f ); if ( real_num != num ) { - fprintf ( stderr, "!!! Got %d tets, but expected %d.\n", + fprintf ( stderr, "!!! Got %lld tets, but expected %lld.\n", real_num, num ); return -1; } @@ -73,7 +86,7 @@ read_tets ( h5_err_t read_level ( - h5_file * f + h5_file_t * f ) { h5_err_t h5err = read_vertices ( f ); if ( h5err < 0 ) { @@ -88,19 +101,18 @@ read_level ( return H5_SUCCESS; } - h5_err_t read_mesh ( - h5_file * f + h5_file_t * f ) { h5_id_t level_id; h5_size_t num_levels = H5FedGetNumLevels ( f ); - printf ( " Number of levels in mesh: %d\n", num_levels ); + printf ( " Number of levels in mesh: %lld\n", num_levels ); for ( level_id = 0; level_id < num_levels; level_id++ ) { h5_err_t h5err = H5FedSetLevel ( f, level_id ); if ( h5err < 0 ) { - fprintf ( stderr, "!!! Can't set level %d.\n", level_id ); + fprintf ( stderr, "!!! Can't set level %lld.\n", level_id ); return -1; } h5err = read_level ( f ); @@ -120,22 +132,22 @@ main ( char *argv[] ) { - H5PartSetVerbosityLevel ( 4 ); + H5SetVerbosityLevel ( 4 ); - h5_file *f = H5OpenFile ( "simple_tet.h5", 0 ); + h5_file_t *f = H5OpenFile ( "simple_tet.h5", H5_O_RDONLY ); if ( f == NULL ) { fprintf ( stderr, "!!! Can't open file.\n" ); return -1; } - h5_size_t num_meshes = H5FedGetNumMeshes ( f, TETRAHEDRAL_MESH ); - printf ( " Number of meshes: %d\n", num_meshes ); + h5_size_t num_meshes = H5FedGetNumMeshes ( f, H5_TETRAHEDRAL_MESH ); + printf ( " Number of meshes: %lld\n", num_meshes ); h5_id_t mesh_id; for ( mesh_id = 0; mesh_id < num_meshes; mesh_id++ ) { - h5_err_t h5err = H5FedOpenMesh ( f, mesh_id, TETRAHEDRAL_MESH ); + h5_err_t h5err = H5FedOpenMesh ( f, mesh_id, H5_TETRAHEDRAL_MESH ); if ( h5err < 0 ) { - fprintf ( stderr, "!!! Can't open mesh %d\n", mesh_id ); + fprintf ( stderr, "!!! Can't open mesh %lld\n", mesh_id ); return -1; } h5err = read_mesh ( f ); diff --git a/test/H5Fed/read_trianglemesh.c b/test/H5Fed/read_trianglemesh.c index 916f319..382c33a 100644 --- a/test/H5Fed/read_trianglemesh.c +++ b/test/H5Fed/read_trianglemesh.c @@ -25,7 +25,7 @@ typedef struct entity entity_t; h5_err_t read_vertices ( - h5_file * f + h5_file_t * f ) { h5_id_t id, local_id; h5_float64_t P[3]; @@ -33,18 +33,17 @@ read_vertices ( h5_id_t level_id = H5FedGetLevel ( f ); h5_size_t num = H5FedGetNumVertices ( f ); - - printf ( " Number of vertices on level %d: %d\n", level_id, num ); + printf ( " Number of vertices on level %lld: %lld\n", level_id, num ); h5_err_t h5err = H5FedStartTraverseVertices ( f ); if ( h5err < 0 ) return h5err; while ( (real_num < num) && ((local_id = H5FedTraverseVertices ( f, &id, P )) >= 0) ) { - printf ( " Vertex[%d]: local id: %d, coords: %f %f %f \n", + printf ( " Vertex[%lld]: local id: %lld, coords: %f %f %f \n", id, local_id, P[0], P[1], P[2] ); real_num++; } if ( real_num != num ) { - fprintf ( stderr, "!!! Got %d vertices, but expected %d.\n", + fprintf ( stderr, "!!! Got %lld vertices, but expected %lld.\n", real_num, num ); return -1; } @@ -53,26 +52,28 @@ read_vertices ( h5_err_t read_entities ( - h5_file * f + h5_file_t * f ) { h5_id_t id, local_id, parent_id, vids[4]; h5_size_t real_num = 0; h5_id_t level_id = H5FedGetLevel ( f ); - h5_size_t num = H5FedGetNumTriangles ( f ); - printf ( " Number of triangles on level %d: %d\n", level_id, num ); + h5_size_t num = H5FedGetNumElements ( f ); + printf ( " Number of triangles on level %lld: %lld\n", level_id, num ); - h5_err_t h5err = H5FedStartTraverseTriangles ( f ); + h5_err_t h5err = H5FedStartTraverseElements ( f ); if ( h5err < 0 ) return h5err; while ( (real_num < num) && - ((local_id = H5FedTraverseTriangles ( f, &id, &parent_id, vids )) >= 0) ) { + ((local_id = H5FedTraverseElements ( + f, &id, &parent_id, vids )) >= 0) ) { printf ( - " entitiy[%d]: local id: %d, parent id: %d, vids: %d %d %d\n", + " entitiy[%lld]: local id: %lld, parent id: %lld, " + "vids: %lld %lld %lld\n", id, local_id, parent_id, vids[0], vids[1], vids[2] ); real_num++; } if ( real_num != num ) { - fprintf ( stderr, "!!! Got %d tets, but expected %d.\n", + fprintf ( stderr, "!!! Got %lld tets, but expected %lld.\n", real_num, num ); return -1; } @@ -82,7 +83,7 @@ read_entities ( h5_err_t read_level ( - h5_file * f + h5_file_t * f ) { h5_err_t h5err = read_vertices ( f ); if ( h5err < 0 ) { @@ -97,10 +98,9 @@ read_level ( return H5_SUCCESS; } - h5_err_t read_mesh ( - h5_file * f + h5_file_t * f ) { h5_id_t level_id; @@ -109,12 +109,12 @@ read_mesh ( fprintf ( stderr, "!!! Oops ...\n" ); return -1; } - printf ( " Number of levels in mesh: %d\n", num_levels ); + printf ( " Number of levels in mesh: %lld\n", num_levels ); for ( level_id = 0; level_id < num_levels; level_id++ ) { - printf ( " Going to level %d\n", level_id ); + printf ( " Going to level %lld\n", level_id ); h5_err_t h5err = H5FedSetLevel ( f, level_id ); if ( h5err < 0 ) { - fprintf ( stderr, "!!! Can't set level %d.\n", level_id ); + fprintf ( stderr, "!!! Can't set level %lld.\n", level_id ); return -1; } h5err = read_level ( f ); @@ -134,22 +134,22 @@ main ( char *argv[] ) { - H5PartSetVerbosityLevel ( 4 ); + H5SetVerbosityLevel ( 4 ); - h5_file *f = H5OpenFile ( "simple_triangle.h5", 0 ); + h5_file_t *f = H5OpenFile ( "simple_triangle.h5", H5_O_RDONLY ); if ( f == NULL ) { fprintf ( stderr, "!!! Can't open file.\n" ); return -1; } - h5_size_t num_meshes = H5FedGetNumMeshes ( f, TRIANGLE_MESH ); - printf ( " Number of meshes: %d\n", num_meshes ); + h5_size_t num_meshes = H5FedGetNumMeshes ( f, H5_TRIANGLE_MESH ); + printf ( " Number of meshes: %lld\n", num_meshes ); h5_id_t mesh_id; for ( mesh_id = 0; mesh_id < num_meshes; mesh_id++ ) { - h5_err_t h5err = H5FedOpenMesh ( f, mesh_id, TRIANGLE_MESH ); + h5_err_t h5err = H5FedOpenMesh ( f, mesh_id, H5_TRIANGLE_MESH ); if ( h5err < 0 ) { - fprintf ( stderr, "!!! Can't open mesh %d\n", mesh_id ); + fprintf ( stderr, "!!! Can't open mesh %lld\n", mesh_id ); return -1; } h5err = read_mesh ( f ); diff --git a/test/H5Fed/write_tetmesh.c b/test/H5Fed/write_tetmesh.c index 0dbb4b5..ab613b0 100644 --- a/test/H5Fed/write_tetmesh.c +++ b/test/H5Fed/write_tetmesh.c @@ -48,15 +48,15 @@ main ( int argc, char *argv[] ) { - H5SetVerbosityLevel ( 5 ); + H5SetVerbosityLevel ( 4 ); - h5_file_t *f = H5OpenFile ( "simple_tet.h5", 0 ); + h5_file_t *f = H5OpenFile ( "simple_tet.h5", H5_O_WRONLY ); if ( f == NULL ) { fprintf ( stderr, "!!! Can't open file.\n" ); return -1; } - h5_err_t h5err = H5FedAddTetMesh ( f, 2 ); + h5_err_t h5err = H5FedAddMesh ( f, 2, H5_TETRAHEDRAL_MESH ); if ( h5err < 0 ) { fprintf ( stderr, "!!! Can't add mesh.\n" ); return -1; @@ -84,7 +84,7 @@ main ( } } - h5_id_t level_id = H5FedAddLevel( f, 8 ); + h5_id_t level_id = H5FedAddLevel( f, 1 ); if ( level_id < 0 ) { fprintf ( stderr, "!!! Can't add level.\n" ); return -1; diff --git a/test/H5Fed/write_trianglemesh.c b/test/H5Fed/write_trianglemesh.c index 1471a14..517fedd 100644 --- a/test/H5Fed/write_trianglemesh.c +++ b/test/H5Fed/write_trianglemesh.c @@ -25,105 +25,72 @@ struct entity { typedef struct entity entity_t; -vertex_t V0[4] = { +vertex_t V[4] = { { 0, {-1.0, 0.0, 0.0} }, { 1, { 1.0, 0.0, 0.0} }, { 2, { 0.0, 1.0, 0.0} }, { 3, { 0.0, -1.0, 0.0} } }; -vertex_t V1[1] = { - { 4, {0.0, 0.0, 0.0 } } -}; - -entity_t T0[2] = { +entity_t T[2] = { { 1, -1, { 0, 1, 2 } }, { 0, -1, { 0, 1, 3 } } }; -entity_t T1[4] = { - { 2, 1, { 1, 2, 4 } }, - { 3, 1, { 0, 2, 4 } }, - { 4, 0, { 0, 3, 4 } }, - { 5, 0, { 1, 3, 4 } } -}; - -h5_err_t -add_level ( - h5_file *f, - vertex_t V[], - int num_verts, - entity_t T[], - int num_entities - ) { - - h5_err_t h5err = H5FedAddLevel ( f ); - if ( h5err < 0 ) { - fprintf ( stderr, "!!! Can't add level.\n" ); - return -1; - } - h5err = H5FedAddNumVertices ( f, num_verts ); - if ( h5err < 0 ) { - fprintf ( stderr, "!!! Can't set number of vertices.\n" ); - return -1; - } - - int i; - for ( i = 0; i