diff --git a/src/h5core/h5t_openclose.c b/src/h5core/h5t_openclose.c index de7be72..3289e68 100644 --- a/src/h5core/h5t_openclose.c +++ b/src/h5core/h5t_openclose.c @@ -222,7 +222,6 @@ _init_fdata ( t->num_meshes = -1; t->cur_mesh = -1; t->num_levels = -1; - t->new_level = -1; t->cur_level = -1; t->last_stored_vid = -1; t->last_stored_eid = -1; @@ -544,14 +543,19 @@ h5t_close_mesh ( h5_err_t h5t_set_level ( h5_file_t* const f, - const h5_id_t id + const h5_id_t level_id ) { h5t_fdata_t* t = f->t; - if ((id < 0) || (id >= t->num_levels)) - return HANDLE_H5_OUT_OF_RANGE_ERR (f, "Level", id); - t->cur_level = id; + if ((level_id < 0) || (level_id >= t->num_levels)) + return HANDLE_H5_OUT_OF_RANGE_ERR (f, "Level", level_id); + h5_id_t prev_level = t->cur_level; + t->cur_level = level_id; + + if (level_id >= t->num_loaded_levels) { + TRY( (t->methods.adjacency->update_internal_structs)(f, prev_level+1) ); + } return H5_SUCCESS; } diff --git a/src/h5core/h5t_readwrite.c b/src/h5core/h5t_readwrite.c index 2ec9cca..a3ffde1 100644 --- a/src/h5core/h5t_readwrite.c +++ b/src/h5core/h5t_readwrite.c @@ -284,7 +284,8 @@ read_elems ( ) { h5t_fdata_t* t = f->t; - TRY( (*t->methods.store->alloc_elems)(f, 0, t->num_elems[t->num_levels-1]) ); + TRY( (*t->methods.store->alloc_elems)( + f, 0, t->num_elems[t->num_levels-1]) ); TRY( h5priv_read_dataset_by_name ( f, t->mesh_gid, @@ -295,9 +296,8 @@ read_elems ( TRY( h5tpriv_sort_elems (f) ); TRY( h5tpriv_rebuild_global_2_local_map_of_elems (f) ); - TRY( build_elems_ldta (f) ); - TRY( (*t->methods.adjacency->update_internal_structs)(f) ); + return H5_SUCCESS; } @@ -314,6 +314,7 @@ h5_err_t h5tpriv_read_mesh ( h5_file_t* const f ) { + h5t_fdata_t* t = f->t; if (f->t->mesh_gid < 0) { TRY( h5tpriv_open_mesh_group (f) ); } @@ -322,6 +323,8 @@ h5tpriv_read_mesh ( TRY( read_num_elems (f) ); TRY( read_vertices (f) ); TRY( read_elems (f) ); + TRY( (t->methods.adjacency->update_internal_structs)(f, 0) ); TRY( read_mtags (f) ); + t->num_loaded_levels = t->num_levels; return H5_SUCCESS; } diff --git a/src/h5core/h5t_storemesh.c b/src/h5core/h5t_storemesh.c index 96ea7aa..e939797 100644 --- a/src/h5core/h5t_storemesh.c +++ b/src/h5core/h5t_storemesh.c @@ -25,63 +25,59 @@ h5t_add_mesh ( /* - * Assign unique global id's to vertices. Vertices already have - unique id's assigned by the mesher but this id's may not be - consecutive numbered starting from 0. - * Set the global vertex id's in element definitions. + * Assign unique global indices to vertices. */ static h5_err_t -assign_global_vertex_ids ( +assign_global_vertex_indices ( h5_file_t* const f ) { h5t_fdata_t* const t = f->t; - h5_id_t local_id; if (t->cur_level < 0) return H5_SUCCESS; /* no level defined */ /* simple in serial runs: global_id = local_id */ - for (local_id = 0; - local_id < t->num_vertices[t->num_levels-1]; - local_id++) { - t->vertices[local_id].global_vid = local_id; + h5_id_t local_idx = (t->cur_level == 0) ? 0 : t->num_vertices[t->cur_level-1]; + for (local_idx = 0; + local_idx < t->num_vertices[t->num_levels-1]; + local_idx++) { + t->vertices[local_idx].global_vid = local_idx; } return H5_SUCCESS; } /*! - Assign unique global IDs to new elements. + Assign unique global indices to new elements. */ static h5_err_t -assign_global_elem_ids ( +assign_global_elem_indices ( h5_file_t* const f ) { h5t_fdata_t* const t = f->t; - h5_id_t local_eid; - if ( t->cur_level < 0 ) return H5_SUCCESS; /* no level defined */ + if (t->cur_level < 0) return H5_SUCCESS; /* no level defined */ /* - simple in serial runs: global_id = local_id + simple in serial runs: global index = local index */ - for (local_eid = (t->cur_level == 0) ? 0 : t->num_elems[t->cur_level-1]; - local_eid < t->num_elems[t->cur_level]; - local_eid++) { + h5_id_t local_idx = (t->cur_level == 0) ? 0 : t->num_elems[t->cur_level-1]; + + for (; local_idx < t->num_elems[t->cur_level]; local_idx++) { h5_elem_t *elem; - h5_elem_ldta_t *elem_ldta = &t->elems_ldta[local_eid]; + h5_elem_ldta_t *elem_ldta = &t->elems_ldta[local_idx]; switch ( t->mesh_type ) { case H5_OID_TETRAHEDRON: - elem = (h5_elem_t*)&t->elems.tets[local_eid]; + elem = (h5_elem_t*)&t->elems.tets[local_idx]; break; case H5_OID_TRIANGLE: - elem = (h5_elem_t*)&t->elems.tris[local_eid]; + elem = (h5_elem_t*)&t->elems.tris[local_idx]; break; default: return h5_error_internal (f,__FILE__,__func__,__LINE__); } - elem->global_eid = local_eid; + elem->global_eid = local_idx; elem->global_parent_eid = elem_ldta->local_parent_eid; elem->global_child_eid = elem_ldta->local_child_eid; int i; @@ -106,9 +102,10 @@ h5t_add_level ( /* t->num_levels will be set to zero on file creation(!) */ if ((t->cur_mesh < 0) || (t->num_levels == -1)) { - return h5tpriv_error_undef_mesh ( f ); + return h5tpriv_error_undef_mesh (f); } t->cur_level = t->num_levels++; + t->num_loaded_levels = t->num_levels; t->dsinfo_num_vertices.dims[0] = t->num_levels; t->dsinfo_num_elems.dims[0] = t->num_levels; t->dsinfo_num_elems_on_level.dims[0] = t->num_levels; @@ -123,7 +120,6 @@ h5t_add_level ( f, t->num_elems_on_level, num_bytes) ); t->num_elems_on_level[t->cur_level] = -1; - t->new_level = t->cur_level; if (t->cur_level == 0) { /* nothing stored yet */ t->last_stored_vid = -1; @@ -146,7 +142,6 @@ h5t_begin_store_vertices ( if (t->cur_level < 0) { return h5tpriv_error_undef_level(f); } - t->storing_data = 1; h5_size_t cur_num_vertices = (t->cur_level > 0 ? t->num_vertices[t->cur_level-1] : 0); t->num_vertices[t->cur_level] = cur_num_vertices+num; @@ -175,7 +170,6 @@ h5t_store_vertex ( if (t->cur_level < 0) return h5tpriv_error_undef_level(f); - t->level_changed = 1; h5_id_t local_id = ++t->last_stored_vid; h5_vertex_t *vertex = &t->vertices[local_id]; vertex->global_vid = global_vid; /* ID from mesher, replaced later!*/ @@ -188,10 +182,9 @@ h5t_end_store_vertices ( h5_file_t* const f ) { h5t_fdata_t* const t = f->t; - t->storing_data = 0; t->num_vertices[t->cur_level] = t->last_stored_vid+1; - TRY( assign_global_vertex_ids (f) ); + TRY( assign_global_vertex_indices (f) ); TRY( h5tpriv_sort_vertices (f) ); TRY( h5tpriv_rebuild_global_2_local_map_of_vertices (f) ); return H5_SUCCESS; @@ -210,7 +203,6 @@ h5t_begin_store_elems ( ) { h5t_fdata_t* const t = f->t; - t->storing_data = 1; size_t cur = t->cur_level > 0 ? t->num_elems[t->cur_level-1] : 0; size_t new = num + cur; t->num_elems[t->cur_level] = new; @@ -220,7 +212,7 @@ h5t_begin_store_elems ( num + t->num_elems_on_level[t->cur_level-1] : num; /* We allocate a hash table for a minimum of 2^21 edges to - prevent resizing. + avoid resizing. */ size_t nel = 2097152 > 5*new ? 2097152 : 5*new; TRY( h5tpriv_resize_te_htab (f, nel) ); @@ -243,7 +235,7 @@ h5t_store_elem ( const h5_id_t elem_idx_of_parent, const h5_id_t* vertices ) { - h5t_fdata_t *t = f->t; + h5t_fdata_t* t = f->t; /* level set? */ if (t->cur_level < 0) @@ -256,18 +248,16 @@ h5t_store_elem ( t->num_elems[t->cur_level] ); /* check parent id */ - if ( - (t->cur_level == 0 && elem_idx_of_parent != -1) || - (t->cur_level > 0 && elem_idx_of_parent < 0) || - (t->cur_level > 0 - && elem_idx_of_parent >= t->num_elems[t->cur_level-1]) + if ((t->cur_level == 0 && elem_idx_of_parent != -1) || + (t->cur_level > 0 && elem_idx_of_parent < 0) || + (t->cur_level > 0 + && elem_idx_of_parent >= t->num_elems[t->cur_level-1]) ) { return HANDLE_H5_PARENT_ID_ERR ( f, h5tpriv_map_oid2str (t->mesh_type), elem_idx_of_parent); } - t->level_changed = 1; h5_id_t elem_idx = ++t->last_stored_eid; h5_elem_ldta_t* elem_ldta = &t->elems_ldta[elem_idx]; @@ -297,13 +287,13 @@ h5t_end_store_elems ( h5_file_t* const f ) { h5t_fdata_t* const t = f->t; - t->storing_data = 0; t->num_elems[t->cur_level] = t->last_stored_eid+1; - TRY( assign_global_elem_ids (f) ); + TRY( assign_global_elem_indices (f) ); TRY( h5tpriv_sort_elems (f) ); TRY( h5tpriv_rebuild_global_2_local_map_of_elems (f) ); + TRY( (t->methods.adjacency->update_internal_structs)(f, t->cur_level) ); return H5_SUCCESS; } @@ -362,7 +352,6 @@ h5t_pre_refine ( default: return h5_error_internal (f, __FILE__, __func__, __LINE__); } - t->storing_data = 1; TRY( h5t_begin_store_vertices (f, num_vertices_to_add) ); TRY( h5t_begin_store_elems (f, num_elems_to_add) ); @@ -377,7 +366,6 @@ h5t_refine ( h5_file_t* const f ) { h5t_fdata_t* const t = f->t; - t->storing_data = 1; int i; for (i = 0; i < t->marked_entities.num_items; i++) { TRY( h5t_refine_elem (f, t->marked_entities.items[i]) ); @@ -390,6 +378,8 @@ h5t_post_refine ( h5_file_t* const f ) { h5t_fdata_t* const t = f->t; + TRY( h5t_end_store_vertices (f) ); + TRY( h5t_end_store_elems (f) ); return h5priv_free_idlist_items (f, &t->marked_entities); } diff --git a/src/h5core/h5t_tetm_adjacencies.c b/src/h5core/h5t_tetm_adjacencies.c index 5780bde..839c7de 100644 --- a/src/h5core/h5t_tetm_adjacencies.c +++ b/src/h5core/h5t_tetm_adjacencies.c @@ -20,21 +20,24 @@ compute T(V) from current level up to highest levels. */ static h5_err_t -compute_tets_of_vertices ( - h5_file_t* const f +compute_elems_of_vertices ( + h5_file_t* const f, + const h5_id_t from_lvl ) { h5t_fdata_t* t = f->t; - h5_id_t eid = (t->cur_level <= 0) ? 0 : t->num_elems[t->cur_level-1]; - h5_elem_ldta_t* tet = tet = &t->elems_ldta[eid]; + h5_id_t elem_idx = (from_lvl <= 0) ? 0 : t->num_elems[from_lvl-1]; + h5_elem_ldta_t *el = &t->elems_ldta[elem_idx]; h5_id_t num_elems = t->num_elems[t->num_levels-1]; - for (;eid < num_elems; eid++, tet++) { - int i; - for (i = 0; i < 4; i++) { - h5_id_t vid = tet->local_vids[i]; - TRY ( h5priv_append_to_idlist ( - f, - &t->vertices_data[vid].tv, - h5tpriv_build_vertex_id( i, eid)) ); + for (;elem_idx < num_elems; elem_idx++, el++) { + int face_idx; + int num_faces = t->ref_element->num_faces[0]; + for (face_idx = 0; face_idx < num_faces; face_idx++) { + h5_id_t vidx = el->local_vids[face_idx]; + TRY( h5priv_append_to_idlist ( + f, + &t->vertices_data[vidx].tv, + h5tpriv_build_vertex_id ( + face_idx, elem_idx)) ); } } return H5_SUCCESS; @@ -44,38 +47,44 @@ compute_tets_of_vertices ( Compute T(E) from current level up to highest levels. */ static h5_err_t -compute_tets_of_edges ( - h5_file_t* const f +compute_elems_of_edges ( + h5_file_t* const f, + const h5_id_t from_lvl ) { - h5t_fdata_t* t = f->t; - h5_id_t eid = (t->cur_level <= 0) ? 0 : t->num_elems[t->cur_level-1]; - h5_elem_ldta_t* tet = tet = &t->elems_ldta[eid]; + h5t_fdata_t *t = f->t; + h5_id_t elem_idx = (from_lvl <= 0) ? 0 : t->num_elems[from_lvl-1]; + h5_elem_ldta_t *el = &t->elems_ldta[elem_idx]; h5_id_t num_elems = t->num_elems[t->num_levels-1]; h5_idlist_t *retval = NULL; - TRY( h5tpriv_resize_te_htab (f, 4*(num_elems-eid)) ); - for (; eid < num_elems; eid++, tet++) { - h5_id_t face; - for (face = 0; face < 6; face++) { - TRY( h5tpriv_search_te2 (f, face, eid, &retval) ); + TRY( h5tpriv_resize_te_htab (f, 4*(num_elems-elem_idx)) ); + for (;elem_idx < num_elems; elem_idx++, el++) { + int face_idx; + int num_faces = t->ref_element->num_faces[1]; + for (face_idx = 0; face_idx < num_faces; face_idx++) { + TRY( h5tpriv_search_te2 ( + f, face_idx, elem_idx, &retval) ); } } return H5_SUCCESS; } static h5_err_t -compute_tets_of_triangles ( - h5_file_t* const f +compute_elems_of_triangles ( + h5_file_t* const f, + const h5_id_t from_lvl ) { h5t_fdata_t* t = f->t; - h5_id_t eid = (t->cur_level <= 0) ? 0 : t->num_elems[t->cur_level-1]; - h5_elem_ldta_t* tet = tet = &t->elems_ldta[eid]; + h5_id_t elem_idx = (from_lvl <= 0) ? 0 : t->num_elems[from_lvl-1]; + h5_elem_ldta_t *el = &t->elems_ldta[elem_idx]; h5_id_t num_elems = t->num_elems[t->num_levels-1]; - h5_idlist_t* retval = NULL; - TRY( h5tpriv_resize_td_htab (f, 4*(num_elems-eid)) ); - for (; eid < num_elems; eid++, tet++) { - h5_id_t face; - for (face = 0; face < 4; face++) { - TRY( h5tpriv_search_td2 (f, face, eid, &retval) ); + h5_idlist_t *retval = NULL; + TRY( h5tpriv_resize_te_htab (f, 4*(num_elems-elem_idx)) ); + for (;elem_idx < num_elems; elem_idx++, el++) { + int face_idx; + int num_faces = t->ref_element->num_faces[1]; + for (face_idx = 0; face_idx < num_faces; face_idx++) { + TRY( h5tpriv_search_td2 ( + f, face_idx, elem_idx, &retval) ); } } return H5_SUCCESS; @@ -274,7 +283,7 @@ compute_sections_of_triangle ( h5_id_t* tri = td->items; h5_id_t* end = td->items+td->num_items; int refined = 0; - for (; tri < end; tri++ { + for (; tri < end; tri++) { h5_id_t eid = h5tpriv_get_elem_idx (*tri); h5_id_t face_idx = h5tpriv_get_face_idx (*tri); h5_elem_ldta_t *tet = &t->elems_ldta[eid]; @@ -643,20 +652,21 @@ get_triangles_downadjacent_to_tet ( static h5_err_t update_internal_structs ( - h5_file_t* const f + h5_file_t* const f, + const h5_id_t from_lvl ) { clock_t t1 = clock(); - TRY( compute_tets_of_vertices (f) ); + TRY( compute_elems_of_vertices (f, from_lvl) ); clock_t t2 = clock(); fprintf (stderr, "compute_tets_of_vertices(): %f\n", (float)(t2-t1)/CLOCKS_PER_SEC); t1 = clock(); - TRY( compute_tets_of_edges (f) ); + TRY( compute_elems_of_edges (f, from_lvl) ); t2 = clock(); fprintf (stderr, "compute_tets_of_edge(): %f\n", (float)(t2-t1)/CLOCKS_PER_SEC); t1 = clock(); - TRY( compute_tets_of_triangles (f) ); + TRY( compute_elems_of_triangles (f, from_lvl) ); t2 = clock(); fprintf (stderr, "compute_tets_of_triangle(): %f\n", (float)(t2-t1)/CLOCKS_PER_SEC); diff --git a/src/h5core/h5t_tetm_store.c b/src/h5core/h5t_tetm_store.c index 4942199..61c7e16 100644 --- a/src/h5core/h5t_tetm_store.c +++ b/src/h5core/h5t_tetm_store.c @@ -75,11 +75,11 @@ get_direct_children_of_edge ( {2,3} // edge 5 }; int num_faces = f->t->ref_element->num_faces[1]; - if ( ( face_idx < 0 ) || ( face_idx >= num_faces ) ) { - return h5_error_internal ( f, __FILE__, __func__, __LINE__ ); + if ((face_idx < 0) || ( face_idx >= num_faces)) { + return h5_error_internal (f, __FILE__, __func__, __LINE__); } - children[0] = h5tpriv_build_edge_id ( face_idx, elem_idx+offs[face_idx][0] ); - children[1] = h5tpriv_build_edge_id ( face_idx, elem_idx+offs[face_idx][1] ); + children[0] = h5tpriv_build_edge_id (face_idx, elem_idx+offs[face_idx][0]); + children[1] = h5tpriv_build_edge_id (face_idx, elem_idx+offs[face_idx][1]); return H5_SUCCESS; } @@ -118,14 +118,14 @@ bisect_edge ( h5_id_t face_id = h5tpriv_get_face_idx ( retval->items[i] ); h5_id_t kids[2], edge0[2], edge1[2]; - TRY ( get_direct_children_of_edge ( - f, - face_id, - tet->local_child_eid, - kids ) ); - TRY ( h5t_get_vertex_indices_of_edge ( f, kids[0], edge0 ) ); - TRY ( h5t_get_vertex_indices_of_edge ( f, kids[1], edge1 ) ); - if ( (edge0[0] == edge1[0]) || (edge0[0] == edge1[1]) ) + TRY( get_direct_children_of_edge ( + f, + face_id, + tet->local_child_eid, + kids) ); + TRY( h5t_get_vertex_indices_of_edge (f, kids[0], edge0) ); + TRY( h5t_get_vertex_indices_of_edge (f, kids[1], edge1) ); + if ((edge0[0] == edge1[0]) || (edge0[0] == edge1[1])) return edge0[0]; else return edge0[1]; diff --git a/src/h5core/h5t_trim_adjacencies.c b/src/h5core/h5t_trim_adjacencies.c index 8cc827a..17598d7 100644 --- a/src/h5core/h5t_trim_adjacencies.c +++ b/src/h5core/h5t_trim_adjacencies.c @@ -21,11 +21,11 @@ */ static h5_err_t compute_elems_of_vertices ( - h5_file_t * const f + h5_file_t* const f, + const h5_id_t from_lvl ) { h5t_fdata_t *t = f->t; - h5_id_t elem_idx = (t->cur_level <= 0) ? - 0 : t->num_elems[t->cur_level-1]; + h5_id_t elem_idx = (from_lvl <= 0) ? 0 : t->num_elems[from_lvl-1]; h5_elem_ldta_t *el = &t->elems_ldta[elem_idx]; h5_id_t num_elems = t->num_elems[t->num_levels-1]; for (;elem_idx < num_elems; elem_idx++, el++) { @@ -48,11 +48,11 @@ compute_elems_of_vertices ( */ static h5_err_t compute_elems_of_edges ( - h5_file_t * const f + h5_file_t* const f, + const h5_id_t from_lvl ) { h5t_fdata_t *t = f->t; - h5_id_t elem_idx = (t->cur_level <= 0) ? - 0 : t->num_elems[t->cur_level-1]; + h5_id_t elem_idx = (from_lvl <= 0) ? 0 : t->num_elems[from_lvl-1]; h5_elem_ldta_t *el = &t->elems_ldta[elem_idx]; h5_id_t num_elems = t->num_elems[t->num_levels-1]; h5_idlist_t *retval = NULL; @@ -68,24 +68,6 @@ compute_elems_of_edges ( return H5_SUCCESS; } -static h5_err_t -update_internal_structs ( - h5_file_t * const f - ) { - clock_t t1 = clock(); - TRY ( compute_elems_of_vertices ( f ) ); - clock_t t2 = clock(); - fprintf ( stderr, "compute_elems_of_vertices(): %f\n", - (float)(t2-t1)/CLOCKS_PER_SEC ); - t1 = clock(); - TRY ( compute_elems_of_edges ( f ) ); - t2 = clock(); - fprintf ( stderr, "compute_elems_of_edge(): %f\n", - (float)(t2-t1)/CLOCKS_PER_SEC ); - - return H5_SUCCESS; -} - static h5_err_t compute_children_of_edge ( h5_file_t * const f, @@ -346,7 +328,7 @@ get_edges_downadjacent_to_triangle ( h5tpriv_build_edge_id (face_idx, elem_idx), children) ); } -- TRY ( h5priv_alloc_idlist (f, list, 8) ); + TRY ( h5priv_alloc_idlist (f, list, 8) ); h5_id_t *edge_idp = children->items; int i; for (i = 0; i < children->num_items; i++, edge_idp++) { @@ -356,10 +338,28 @@ get_edges_downadjacent_to_triangle ( return H5_SUCCESS; } +static h5_err_t +update_internal_structs ( + h5_file_t* const f, + const h5_id_t from_lvl + ) { + clock_t t1 = clock(); + TRY( compute_elems_of_vertices (f, from_lvl) ); + clock_t t2 = clock(); + fprintf (stderr, "compute_elems_of_vertices(): %f\n", + (float)(t2-t1)/CLOCKS_PER_SEC); + t1 = clock(); + TRY( compute_elems_of_edges (f, from_lvl ) ); + t2 = clock(); + fprintf (stderr, "compute_elems_of_edge(): %f\n", + (float)(t2-t1)/CLOCKS_PER_SEC); + + return H5_SUCCESS; +} static h5_err_t release_internal_structs ( - h5_file_t * const f + h5_file_t* const f ) { /* TO BE WRITTEN @@@ */ return H5_SUCCESS; diff --git a/src/h5core/h5t_types_private.h b/src/h5core/h5t_types_private.h index e07989f..84cca6f 100644 --- a/src/h5core/h5t_types_private.h +++ b/src/h5core/h5t_types_private.h @@ -12,8 +12,8 @@ struct h5t_retrieve_methods { }; struct h5t_adjacency_methods { - h5_err_t (*update_internal_structs)(h5_file_t * const); - h5_err_t (*release_internal_structs)(h5_file_t * const); + h5_err_t (*update_internal_structs)(h5_file_t* const, h5_id_t); + h5_err_t (*release_internal_structs)(h5_file_t* const); h5_err_t (*get_edges_upadjacent_to_vertex)( h5_file_t * const, const h5_id_t, h5_idlist_t**); h5_err_t (*get_triangles_upadjacent_to_vertex)( @@ -139,11 +139,10 @@ typedef struct h5t_fdata { h5_id_t mesh_changed; /* true if new or has been changed */ h5_id_t num_meshes; /* number of meshes */ h5_id_t cur_level; /* id of current level */ - h5_id_t new_level; /* idx of the first new level or -1 */ h5_size_t num_levels; /* number of levels */ - h5_id_t level_changed; - h5_id_t storing_data; + h5_id_t num_loaded_levels; + /*** HDF5 IDs ***/ hid_t topo_gid; /* grp id of mesh in current level */ hid_t meshes_gid; /* HDF5 id */