From a43a7c91bfa8a3798b1f1dbaf529971a8cd31e9d Mon Sep 17 00:00:00 2001 From: Mikael Vaillant Date: Thu, 5 Sep 2024 13:40:34 -0400 Subject: [PATCH] Modified read_msh function to take into account comments and other msh file structures --- source/grid/grid_in.cc | 1156 +++++++++-------- tests/grid/grid_in_msh_version_2.cc | 2 + tests/grid/grid_in_msh_version_2.output | 568 ++++++++ .../hole81_commented.msh | 217 ++++ 4 files changed, 1369 insertions(+), 574 deletions(-) create mode 100644 tests/grid/grid_in_msh_version_2/hole81_commented.msh diff --git a/source/grid/grid_in.cc b/source/grid/grid_in.cc index b793d88d00..bfe44ae0ea 100644 --- a/source/grid/grid_in.cc +++ b/source/grid/grid_in.cc @@ -2009,23 +2009,54 @@ GridIn::read_msh(std::istream &in) unsigned int n_vertices; unsigned int n_cells; + int n_entity_blocks = 1; unsigned int dummy; std::string line; + + // set up vector of vertices. This vector is filled using information from the $Elements or $ELM section + std::vector> vertices; + + // set up mapping between numbering in msh-file (nod) and in the vertices vector + std::map vertex_indices; + + // set up array of cells and subcells (faces). In 1d, there is currently no + // standard way in deal.II to pass boundary indicators attached to + // individual vertices, so do this by hand via the boundary_ids_1d array + std::vector> cells; + SubCellData subcelldata; + std::map boundary_ids_1d; + + // Track the number of times each vertex is used in 1D. This determines + // whether or not we can assign a boundary id to a vertex. This is necessary + // because sometimes gmsh saves internal vertices in the $ELEM list in codim + // 1 or codim 2. + std::map vertex_counts; + // This array stores maps from the 'entities' to the 'physical tags' for // points, curves, surfaces and volumes. We use this information later to // assign boundary ids. std::array, 4> tag_maps; + + // Section markers according to mesh format + static const std::string begin_nodes_marker[] = {"$NOD", "$Nodes"}; + static const std::string begin_elements_marker[] = {"$ELM", "$Elements"}; in >> line; - // first determine file format + // first determine file format and skip all preliminary information such as comments unsigned int gmsh_file_format = 0; - if (line == "$NOD") - gmsh_file_format = 10; - else if (line == "$MeshFormat") - gmsh_file_format = 20; - else - AssertThrow(false, ExcInvalidGMSHInput(line)); + const std::string first_line = line; + while (true){ + if (line == "$NOD") + {gmsh_file_format = 10; + break;} + else if (line == "$MeshFormat") + {gmsh_file_format = 20; + break;} + else if (in.eof()) + AssertThrow(false, ExcInvalidGMSHInput(first_line)); + in >> line; + } // if file format is 2.0 or greater then we also have to read the rest of // the header @@ -2046,622 +2077,599 @@ GridIn::read_msh(std::istream &in) // description to synch ourselves with the format 1 handling above in >> line; AssertThrow(line == "$EndMeshFormat", ExcInvalidGMSHInput(line)); + } + + // Loop through every line of the file. Processes information in $Entities, $Nodes/$NOD and $Elements/$ELLM section of the .msh file to build the triangulation + while (!in.eof()){ + + // if the next block is of kind $Entities, parse it + if (line == "$Entities"){ + unsigned long n_points, n_curves, n_surfaces, n_volumes; - in >> line; - // if the next block is of kind $PhysicalNames, ignore it - if (line == "$PhysicalNames") - { - do - { - in >> line; - } - while (line != "$EndPhysicalNames"); - in >> line; - } + in >> n_points >> n_curves >> n_surfaces >> n_volumes; + for (unsigned int i = 0; i < n_points; ++i) + { + // parse point ids + int tag; + unsigned int n_physicals; + double box_min_x, box_min_y, box_min_z, box_max_x, box_max_y, + box_max_z; - // if the next block is of kind $Entities, parse it - if (line == "$Entities") - { - unsigned long n_points, n_curves, n_surfaces, n_volumes; + // we only care for 'tag' as key for tag_maps[0] + if (gmsh_file_format > 40) + { + in >> tag >> box_min_x >> box_min_y >> box_min_z >> + n_physicals; + box_max_x = box_min_x; + box_max_y = box_min_y; + box_max_z = box_min_z; + } + else + { + in >> tag >> box_min_x >> box_min_y >> box_min_z >> + box_max_x >> box_max_y >> box_max_z >> n_physicals; + } + // if there is a physical tag, we will use it as boundary id + // below + AssertThrow(n_physicals < 2, + ExcMessage("More than one tag is not supported!")); + // if there is no physical tag, use 0 as default + int physical_tag = 0; + for (unsigned int j = 0; j < n_physicals; ++j) + in >> physical_tag; + tag_maps[0][tag] = physical_tag; + } + for (unsigned int i = 0; i < n_curves; ++i) + { + // parse curve ids + int tag; + unsigned int n_physicals; + double box_min_x, box_min_y, box_min_z, box_max_x, box_max_y, + box_max_z; + + // we only care for 'tag' as key for tag_maps[1] + in >> tag >> box_min_x >> box_min_y >> box_min_z >> box_max_x >> + box_max_y >> box_max_z >> n_physicals; + // if there is a physical tag, we will use it as boundary id + // below + AssertThrow(n_physicals < 2, + ExcMessage("More than one tag is not supported!")); + // if there is no physical tag, use 0 as default + int physical_tag = 0; + for (unsigned int j = 0; j < n_physicals; ++j) + in >> physical_tag; + tag_maps[1][tag] = physical_tag; + // we don't care about the points associated to a curve, but + // have to parse them anyway because their format is + // unstructured + in >> n_points; + for (unsigned int j = 0; j < n_points; ++j) + in >> tag; + } - in >> n_points >> n_curves >> n_surfaces >> n_volumes; - for (unsigned int i = 0; i < n_points; ++i) + for (unsigned int i = 0; i < n_surfaces; ++i) + { + // parse surface ids + int tag; + unsigned int n_physicals; + double box_min_x, box_min_y, box_min_z, box_max_x, box_max_y, + box_max_z; + + // we only care for 'tag' as key for tag_maps[2] + in >> tag >> box_min_x >> box_min_y >> box_min_z >> box_max_x >> + box_max_y >> box_max_z >> n_physicals; + // if there is a physical tag, we will use it as boundary id + // below + AssertThrow(n_physicals < 2, + ExcMessage("More than one tag is not supported!")); + // if there is no physical tag, use 0 as default + int physical_tag = 0; + for (unsigned int j = 0; j < n_physicals; ++j) + in >> physical_tag; + tag_maps[2][tag] = physical_tag; + // we don't care about the curves associated to a surface, but + // have to parse them anyway because their format is + // unstructured + in >> n_curves; + for (unsigned int j = 0; j < n_curves; ++j) + in >> tag; + } + for (unsigned int i = 0; i < n_volumes; ++i) + { + // parse volume ids + int tag; + unsigned int n_physicals; + double box_min_x, box_min_y, box_min_z, box_max_x, box_max_y, + box_max_z; + + // we only care for 'tag' as key for tag_maps[3] + in >> tag >> box_min_x >> box_min_y >> box_min_z >> box_max_x >> + box_max_y >> box_max_z >> n_physicals; + // if there is a physical tag, we will use it as boundary id + // below + AssertThrow(n_physicals < 2, + ExcMessage("More than one tag is not supported!")); + // if there is no physical tag, use 0 as default + int physical_tag = 0; + for (unsigned int j = 0; j < n_physicals; ++j) + in >> physical_tag; + tag_maps[3][tag] = physical_tag; + // we don't care about the surfaces associated to a volume, but + // have to parse them anyway because their format is + // unstructured + in >> n_surfaces; + for (unsigned int j = 0; j < n_surfaces; ++j) + in >> tag; + } + in >> line; + AssertThrow(line == "$EndEntities", ExcInvalidGMSHInput(line)); + in >> line; + } + + // if the next block is of kind $Nodes or $NOD, parse it + if (line == begin_nodes_marker[gmsh_file_format == 10 ? 0 : 1]){ + + // now read the nodes list + int n_entity_blocks = 1; + if (gmsh_file_format > 40) + { + int min_node_tag; + int max_node_tag; + in >> n_entity_blocks >> n_vertices >> min_node_tag >> max_node_tag; + } + else if (gmsh_file_format == 40) + { + in >> n_entity_blocks >> n_vertices; + } + else + in >> n_vertices; + + vertices.resize(n_vertices); + + { + unsigned int global_vertex = 0; + for (int entity_block = 0; entity_block < n_entity_blocks; ++entity_block) { - // parse point ids - int tag; - unsigned int n_physicals; - double box_min_x, box_min_y, box_min_z, box_max_x, box_max_y, - box_max_z; + int parametric; + unsigned long numNodes; - // we only care for 'tag' as key for tag_maps[0] - if (gmsh_file_format > 40) + if (gmsh_file_format < 40) { - in >> tag >> box_min_x >> box_min_y >> box_min_z >> - n_physicals; - box_max_x = box_min_x; - box_max_y = box_min_y; - box_max_z = box_min_z; + numNodes = n_vertices; + parametric = 0; } else { - in >> tag >> box_min_x >> box_min_y >> box_min_z >> - box_max_x >> box_max_y >> box_max_z >> n_physicals; + // for gmsh_file_format 4.1 the order of tag and dim is reversed, + // but we are ignoring both anyway. + int tagEntity, dimEntity; + in >> tagEntity >> dimEntity >> parametric >> numNodes; } - // if there is a physical tag, we will use it as boundary id - // below - AssertThrow(n_physicals < 2, - ExcMessage("More than one tag is not supported!")); - // if there is no physical tag, use 0 as default - int physical_tag = 0; - for (unsigned int j = 0; j < n_physicals; ++j) - in >> physical_tag; - tag_maps[0][tag] = physical_tag; - } - for (unsigned int i = 0; i < n_curves; ++i) - { - // parse curve ids - int tag; - unsigned int n_physicals; - double box_min_x, box_min_y, box_min_z, box_max_x, box_max_y, - box_max_z; - - // we only care for 'tag' as key for tag_maps[1] - in >> tag >> box_min_x >> box_min_y >> box_min_z >> box_max_x >> - box_max_y >> box_max_z >> n_physicals; - // if there is a physical tag, we will use it as boundary id - // below - AssertThrow(n_physicals < 2, - ExcMessage("More than one tag is not supported!")); - // if there is no physical tag, use 0 as default - int physical_tag = 0; - for (unsigned int j = 0; j < n_physicals; ++j) - in >> physical_tag; - tag_maps[1][tag] = physical_tag; - // we don't care about the points associated to a curve, but - // have to parse them anyway because their format is - // unstructured - in >> n_points; - for (unsigned int j = 0; j < n_points; ++j) - in >> tag; - } - for (unsigned int i = 0; i < n_surfaces; ++i) - { - // parse surface ids - int tag; - unsigned int n_physicals; - double box_min_x, box_min_y, box_min_z, box_max_x, box_max_y, - box_max_z; - - // we only care for 'tag' as key for tag_maps[2] - in >> tag >> box_min_x >> box_min_y >> box_min_z >> box_max_x >> - box_max_y >> box_max_z >> n_physicals; - // if there is a physical tag, we will use it as boundary id - // below - AssertThrow(n_physicals < 2, - ExcMessage("More than one tag is not supported!")); - // if there is no physical tag, use 0 as default - int physical_tag = 0; - for (unsigned int j = 0; j < n_physicals; ++j) - in >> physical_tag; - tag_maps[2][tag] = physical_tag; - // we don't care about the curves associated to a surface, but - // have to parse them anyway because their format is - // unstructured - in >> n_curves; - for (unsigned int j = 0; j < n_curves; ++j) - in >> tag; - } - for (unsigned int i = 0; i < n_volumes; ++i) - { - // parse volume ids - int tag; - unsigned int n_physicals; - double box_min_x, box_min_y, box_min_z, box_max_x, box_max_y, - box_max_z; - - // we only care for 'tag' as key for tag_maps[3] - in >> tag >> box_min_x >> box_min_y >> box_min_z >> box_max_x >> - box_max_y >> box_max_z >> n_physicals; - // if there is a physical tag, we will use it as boundary id - // below - AssertThrow(n_physicals < 2, - ExcMessage("More than one tag is not supported!")); - // if there is no physical tag, use 0 as default - int physical_tag = 0; - for (unsigned int j = 0; j < n_physicals; ++j) - in >> physical_tag; - tag_maps[3][tag] = physical_tag; - // we don't care about the surfaces associated to a volume, but - // have to parse them anyway because their format is - // unstructured - in >> n_surfaces; - for (unsigned int j = 0; j < n_surfaces; ++j) - in >> tag; - } - in >> line; - AssertThrow(line == "$EndEntities", ExcInvalidGMSHInput(line)); - in >> line; - } + std::vector vertex_numbers; + int vertex_number; + if (gmsh_file_format > 40) + for (unsigned long vertex_per_entity = 0; + vertex_per_entity < numNodes; + ++vertex_per_entity) + { + in >> vertex_number; + vertex_numbers.push_back(vertex_number); + } - // if the next block is of kind $PartitionedEntities, ignore it - if (line == "$PartitionedEntities") - { - do - { - in >> line; - } - while (line != "$EndPartitionedEntities"); - in >> line; - } + for (unsigned long vertex_per_entity = 0; vertex_per_entity < numNodes; + ++vertex_per_entity, ++global_vertex) + { + int vertex_number; + double x[3]; - // but the next thing should, - // in any case, be the list of - // nodes: - AssertThrow(line == "$Nodes", ExcInvalidGMSHInput(line)); - } + // read vertex + if (gmsh_file_format > 40) + { + vertex_number = vertex_numbers[vertex_per_entity]; + in >> x[0] >> x[1] >> x[2]; + } + else + in >> vertex_number >> x[0] >> x[1] >> x[2]; - // now read the nodes list - int n_entity_blocks = 1; - if (gmsh_file_format > 40) - { - int min_node_tag; - int max_node_tag; - in >> n_entity_blocks >> n_vertices >> min_node_tag >> max_node_tag; - } - else if (gmsh_file_format == 40) - { - in >> n_entity_blocks >> n_vertices; - } - else - in >> n_vertices; - std::vector> vertices(n_vertices); - // set up mapping between numbering - // in msh-file (nod) and in the - // vertices vector - std::map vertex_indices; + for (unsigned int d = 0; d < spacedim; ++d) + vertices[global_vertex][d] = x[d]; + // store mapping + vertex_indices[vertex_number] = global_vertex; - { - unsigned int global_vertex = 0; - for (int entity_block = 0; entity_block < n_entity_blocks; ++entity_block) - { - int parametric; - unsigned long numNodes; + // ignore parametric coordinates + if (parametric != 0) + { + double u = 0.; + double v = 0.; + in >> u >> v; + (void)u; + (void)v; + } + } + } + AssertDimension(global_vertex, n_vertices); + } - if (gmsh_file_format < 40) + // Assert we reached the end of the block + in >> line; + static const std::string end_nodes_marker[] = {"$ENDNOD", "$EndNodes"}; + AssertThrow(line == end_nodes_marker[gmsh_file_format == 10 ? 0 : 1], + ExcInvalidGMSHInput(line)); + + } + + // if the next block is of kind $Elements or $ELM, parse it + if (line == begin_elements_marker[gmsh_file_format == 10 ? 0 : 1]){ + // now read the cell list + if (gmsh_file_format > 40) { - numNodes = n_vertices; - parametric = 0; + int min_node_tag; + int max_node_tag; + in >> n_entity_blocks >> n_cells >> min_node_tag >> max_node_tag; } - else + else if (gmsh_file_format == 40) { - // for gmsh_file_format 4.1 the order of tag and dim is reversed, - // but we are ignoring both anyway. - int tagEntity, dimEntity; - in >> tagEntity >> dimEntity >> parametric >> numNodes; + in >> n_entity_blocks >> n_cells; } - - std::vector vertex_numbers; - int vertex_number; - if (gmsh_file_format > 40) - for (unsigned long vertex_per_entity = 0; - vertex_per_entity < numNodes; - ++vertex_per_entity) - { - in >> vertex_number; - vertex_numbers.push_back(vertex_number); - } - - for (unsigned long vertex_per_entity = 0; vertex_per_entity < numNodes; - ++vertex_per_entity, ++global_vertex) + else { - int vertex_number; - double x[3]; - - // read vertex - if (gmsh_file_format > 40) - { - vertex_number = vertex_numbers[vertex_per_entity]; - in >> x[0] >> x[1] >> x[2]; - } - else - in >> vertex_number >> x[0] >> x[1] >> x[2]; + n_entity_blocks = 1; + in >> n_cells; + } - for (unsigned int d = 0; d < spacedim; ++d) - vertices[global_vertex][d] = x[d]; - // store mapping - vertex_indices[vertex_number] = global_vertex; + // // set up array of cells and subcells (faces). In 1d, there is currently no + // // standard way in deal.II to pass boundary indicators attached to + // // individual vertices, so do this by hand via the boundary_ids_1d array + // std::vector> cells; + // SubCellData subcelldata; + // std::map boundary_ids_1d; - // ignore parametric coordinates - if (parametric != 0) - { - double u = 0.; - double v = 0.; - in >> u >> v; - (void)u; - (void)v; - } - } - } - AssertDimension(global_vertex, n_vertices); - } + // // Track the number of times each vertex is used in 1D. This determines + // // whether or not we can assign a boundary id to a vertex. This is necessary + // // because sometimes gmsh saves internal vertices in the $ELEM list in codim + // // 1 or codim 2. + // std::map vertex_counts; - // Assert we reached the end of the block - in >> line; - static const std::string end_nodes_marker[] = {"$ENDNOD", "$EndNodes"}; - AssertThrow(line == end_nodes_marker[gmsh_file_format == 10 ? 0 : 1], - ExcInvalidGMSHInput(line)); + { + static constexpr std::array local_vertex_numbering = { + {0, 1, 5, 4, 2, 3, 7, 6}}; + unsigned int global_cell = 0; + for (int entity_block = 0; entity_block < n_entity_blocks; ++entity_block) + { + unsigned int material_id; + unsigned long numElements; + int cell_type; - // Now read in next bit - in >> line; - static const std::string begin_elements_marker[] = {"$ELM", "$Elements"}; - AssertThrow(line == begin_elements_marker[gmsh_file_format == 10 ? 0 : 1], - ExcInvalidGMSHInput(line)); + if (gmsh_file_format < 40) + { + material_id = 0; + cell_type = 0; + numElements = n_cells; + } + else if (gmsh_file_format == 40) + { + int tagEntity, dimEntity; + in >> tagEntity >> dimEntity >> cell_type >> numElements; + material_id = tag_maps[dimEntity][tagEntity]; + } + else + { + // for gmsh_file_format 4.1 the order of tag and dim is reversed, + int tagEntity, dimEntity; + in >> dimEntity >> tagEntity >> cell_type >> numElements; + material_id = tag_maps[dimEntity][tagEntity]; + } - // now read the cell list - if (gmsh_file_format > 40) - { - int min_node_tag; - int max_node_tag; - in >> n_entity_blocks >> n_cells >> min_node_tag >> max_node_tag; - } - else if (gmsh_file_format == 40) - { - in >> n_entity_blocks >> n_cells; - } - else - { - n_entity_blocks = 1; - in >> n_cells; - } + for (unsigned int cell_per_entity = 0; cell_per_entity < numElements; + ++cell_per_entity, ++global_cell) + { + // note that since in the input + // file we found the number of + // cells at the top, there + // should still be input here, + // so check this: + AssertThrow(in.fail() == false, ExcIO()); - // set up array of cells and subcells (faces). In 1d, there is currently no - // standard way in deal.II to pass boundary indicators attached to - // individual vertices, so do this by hand via the boundary_ids_1d array - std::vector> cells; - SubCellData subcelldata; - std::map boundary_ids_1d; + unsigned int nod_num; - // Track the number of times each vertex is used in 1D. This determines - // whether or not we can assign a boundary id to a vertex. This is necessary - // because sometimes gmsh saves internal vertices in the $ELEM list in codim - // 1 or codim 2. - std::map vertex_counts; + /* + For file format version 1, the format of each cell is as + follows: elm-number elm-type reg-phys reg-elem number-of-nodes + node-number-list - { - static constexpr std::array local_vertex_numbering = { - {0, 1, 5, 4, 2, 3, 7, 6}}; - unsigned int global_cell = 0; - for (int entity_block = 0; entity_block < n_entity_blocks; ++entity_block) - { - unsigned int material_id; - unsigned long numElements; - int cell_type; + However, for version 2, the format reads like this: + elm-number elm-type number-of-tags < tag > ... + node-number-list - if (gmsh_file_format < 40) - { - material_id = 0; - cell_type = 0; - numElements = n_cells; - } - else if (gmsh_file_format == 40) - { - int tagEntity, dimEntity; - in >> tagEntity >> dimEntity >> cell_type >> numElements; - material_id = tag_maps[dimEntity][tagEntity]; - } - else - { - // for gmsh_file_format 4.1 the order of tag and dim is reversed, - int tagEntity, dimEntity; - in >> dimEntity >> tagEntity >> cell_type >> numElements; - material_id = tag_maps[dimEntity][tagEntity]; - } + For version 4, we have: + tag(int) numVert(int) ... - for (unsigned int cell_per_entity = 0; cell_per_entity < numElements; - ++cell_per_entity, ++global_cell) - { - // note that since in the input - // file we found the number of - // cells at the top, there - // should still be input here, - // so check this: - AssertThrow(in.fail() == false, ExcIO()); - - unsigned int nod_num; - - /* - For file format version 1, the format of each cell is as - follows: elm-number elm-type reg-phys reg-elem number-of-nodes - node-number-list - - However, for version 2, the format reads like this: - elm-number elm-type number-of-tags < tag > ... - node-number-list - - For version 4, we have: - tag(int) numVert(int) ... - - In the following, we will ignore the element number (we simply - enumerate them in the order in which we read them, and we will - take reg-phys (version 1) or the first tag (version 2, if any - tag is given at all) as material id. For version 4, we already - read the material and the cell type in above. - */ - - unsigned int elm_number = 0; - if (gmsh_file_format < 40) - { - in >> elm_number // ELM-NUMBER - >> cell_type; // ELM-TYPE - } + In the following, we will ignore the element number (we simply + enumerate them in the order in which we read them, and we will + take reg-phys (version 1) or the first tag (version 2, if any + tag is given at all) as material id. For version 4, we already + read the material and the cell type in above. + */ - if (gmsh_file_format < 20) - { - in >> material_id // REG-PHYS - >> dummy // reg_elm - >> nod_num; - } - else if (gmsh_file_format < 40) - { - // read the tags; ignore all but the first one which we will - // interpret as the material_id (for cells) or boundary_id - // (for faces) - unsigned int n_tags; - in >> n_tags; - if (n_tags > 0) - in >> material_id; - else - material_id = 0; + unsigned int elm_number = 0; + if (gmsh_file_format < 40) + { + in >> elm_number // ELM-NUMBER + >> cell_type; // ELM-TYPE + } - for (unsigned int i = 1; i < n_tags; ++i) - in >> dummy; - - if (cell_type == 1) // line - nod_num = 2; - else if (cell_type == 2) // tri - nod_num = 3; - else if (cell_type == 3) // quad - nod_num = 4; - else if (cell_type == 4) // tet - nod_num = 4; - else if (cell_type == 5) // hex - nod_num = 8; - } - else // file format version 4.0 and later - { - // ignore tag - int tag; - in >> tag; - - if (cell_type == 1) // line - nod_num = 2; - else if (cell_type == 2) // tri - nod_num = 3; - else if (cell_type == 3) // quad - nod_num = 4; - else if (cell_type == 4) // tet - nod_num = 4; - else if (cell_type == 5) // hex - nod_num = 8; - } + if (gmsh_file_format < 20) + { + in >> material_id // REG-PHYS + >> dummy // reg_elm + >> nod_num; + } + else if (gmsh_file_format < 40) + { + // read the tags; ignore all but the first one which we will + // interpret as the material_id (for cells) or boundary_id + // (for faces) + unsigned int n_tags; + in >> n_tags; + if (n_tags > 0) + in >> material_id; + else + material_id = 0; + + for (unsigned int i = 1; i < n_tags; ++i) + in >> dummy; + + if (cell_type == 1) // line + nod_num = 2; + else if (cell_type == 2) // tri + nod_num = 3; + else if (cell_type == 3) // quad + nod_num = 4; + else if (cell_type == 4) // tet + nod_num = 4; + else if (cell_type == 5) // hex + nod_num = 8; + } + else // file format version 4.0 and later + { + // ignore tag + int tag; + in >> tag; + + if (cell_type == 1) // line + nod_num = 2; + else if (cell_type == 2) // tri + nod_num = 3; + else if (cell_type == 3) // quad + nod_num = 4; + else if (cell_type == 4) // tet + nod_num = 4; + else if (cell_type == 5) // hex + nod_num = 8; + } - /* `ELM-TYPE' - defines the geometrical type of the N-th element: - `1' - Line (2 nodes, 1 edge). + /* `ELM-TYPE' + defines the geometrical type of the N-th element: + `1' + Line (2 nodes, 1 edge). - `2' - Triangle (3 nodes, 3 edges). + `2' + Triangle (3 nodes, 3 edges). - `3' - Quadrangle (4 nodes, 4 edges). + `3' + Quadrangle (4 nodes, 4 edges). - `4' - Tetrahedron (4 nodes, 6 edges, 6 faces). + `4' + Tetrahedron (4 nodes, 6 edges, 6 faces). - `5' - Hexahedron (8 nodes, 12 edges, 6 faces). + `5' + Hexahedron (8 nodes, 12 edges, 6 faces). - `15' - Point (1 node). - */ + `15' + Point (1 node). + */ - if (((cell_type == 1) && (dim == 1)) || // a line in 1d - ((cell_type == 2) && (dim == 2)) || // a triangle in 2d - ((cell_type == 3) && (dim == 2)) || // a quadrilateral in 2d - ((cell_type == 4) && (dim == 3)) || // a tet in 3d - ((cell_type == 5) && (dim == 3))) // a hex in 3d - // found a cell - { - unsigned int vertices_per_cell = 0; - if (cell_type == 1) // line - vertices_per_cell = 2; - else if (cell_type == 2) // tri - vertices_per_cell = 3; - else if (cell_type == 3) // quad - vertices_per_cell = 4; - else if (cell_type == 4) // tet - vertices_per_cell = 4; - else if (cell_type == 5) // hex - vertices_per_cell = 8; - - AssertThrow(nod_num == vertices_per_cell, - ExcMessage( - "Number of nodes does not coincide with the " - "number required for this object")); - - // allocate and read indices - cells.emplace_back(); - CellData &cell = cells.back(); - cell.vertices.resize(vertices_per_cell); - for (unsigned int i = 0; i < vertices_per_cell; ++i) - { - // hypercube cells need to be reordered - if (vertices_per_cell == - GeometryInfo::vertices_per_cell) - in >> cell.vertices[dim == 3 ? - local_vertex_numbering[i] : - GeometryInfo::ucd_to_deal[i]]; - else - in >> cell.vertices[i]; - } + if (((cell_type == 1) && (dim == 1)) || // a line in 1d + ((cell_type == 2) && (dim == 2)) || // a triangle in 2d + ((cell_type == 3) && (dim == 2)) || // a quadrilateral in 2d + ((cell_type == 4) && (dim == 3)) || // a tet in 3d + ((cell_type == 5) && (dim == 3))) // a hex in 3d + // found a cell + { + unsigned int vertices_per_cell = 0; + if (cell_type == 1) // line + vertices_per_cell = 2; + else if (cell_type == 2) // tri + vertices_per_cell = 3; + else if (cell_type == 3) // quad + vertices_per_cell = 4; + else if (cell_type == 4) // tet + vertices_per_cell = 4; + else if (cell_type == 5) // hex + vertices_per_cell = 8; + + AssertThrow(nod_num == vertices_per_cell, + ExcMessage( + "Number of nodes does not coincide with the " + "number required for this object")); + + // allocate and read indices + cells.emplace_back(); + CellData &cell = cells.back(); + cell.vertices.resize(vertices_per_cell); + for (unsigned int i = 0; i < vertices_per_cell; ++i) + { + // hypercube cells need to be reordered + if (vertices_per_cell == + GeometryInfo::vertices_per_cell) + in >> cell.vertices[dim == 3 ? + local_vertex_numbering[i] : + GeometryInfo::ucd_to_deal[i]]; + else + in >> cell.vertices[i]; + } - // to make sure that the cast won't fail - Assert(material_id <= - std::numeric_limits::max(), - ExcIndexRange( - material_id, - 0, - std::numeric_limits::max())); - // we use only material_ids in the range from 0 to - // numbers::invalid_material_id-1 - AssertIndexRange(material_id, numbers::invalid_material_id); - - cell.material_id = material_id; - - // transform from gmsh to consecutive numbering - for (unsigned int i = 0; i < vertices_per_cell; ++i) - { - AssertThrow(vertex_indices.find(cell.vertices[i]) != - vertex_indices.end(), - ExcInvalidVertexIndexGmsh(cell_per_entity, - elm_number, - cell.vertices[i])); - - const auto vertex = vertex_indices[cell.vertices[i]]; - if (dim == 1) - vertex_counts[vertex] += 1u; - cell.vertices[i] = vertex; - } - } - else if ((cell_type == 1) && - ((dim == 2) || (dim == 3))) // a line in 2d or 3d - // boundary info - { - subcelldata.boundary_lines.emplace_back(); - in >> subcelldata.boundary_lines.back().vertices[0] >> - subcelldata.boundary_lines.back().vertices[1]; - - // to make sure that the cast won't fail - Assert(material_id <= - std::numeric_limits::max(), - ExcIndexRange( - material_id, - 0, - std::numeric_limits::max())); - // we use only boundary_ids in the range from 0 to - // numbers::internal_face_boundary_id-1 - AssertIndexRange(material_id, - numbers::internal_face_boundary_id); - - subcelldata.boundary_lines.back().boundary_id = - static_cast(material_id); - - // transform from ucd to - // consecutive numbering - for (unsigned int &vertex : - subcelldata.boundary_lines.back().vertices) - if (vertex_indices.find(vertex) != vertex_indices.end()) - // vertex with this index exists - vertex = vertex_indices[vertex]; - else + // to make sure that the cast won't fail + Assert(material_id <= + std::numeric_limits::max(), + ExcIndexRange( + material_id, + 0, + std::numeric_limits::max())); + // we use only material_ids in the range from 0 to + // numbers::invalid_material_id-1 + AssertIndexRange(material_id, numbers::invalid_material_id); + + cell.material_id = material_id; + + // transform from gmsh to consecutive numbering + for (unsigned int i = 0; i < vertices_per_cell; ++i) + { + AssertThrow(vertex_indices.find(cell.vertices[i]) != + vertex_indices.end(), + ExcInvalidVertexIndexGmsh(cell_per_entity, + elm_number, + cell.vertices[i])); + + const auto vertex = vertex_indices[cell.vertices[i]]; + if (dim == 1) + vertex_counts[vertex] += 1u; + cell.vertices[i] = vertex; + } + } + else if ((cell_type == 1) && + ((dim == 2) || (dim == 3))) // a line in 2d or 3d + // boundary info { - // no such vertex index - AssertThrow(false, - ExcInvalidVertexIndex(cell_per_entity, - vertex)); - vertex = numbers::invalid_unsigned_int; + subcelldata.boundary_lines.emplace_back(); + in >> subcelldata.boundary_lines.back().vertices[0] >> + subcelldata.boundary_lines.back().vertices[1]; + + // to make sure that the cast won't fail + Assert(material_id <= + std::numeric_limits::max(), + ExcIndexRange( + material_id, + 0, + std::numeric_limits::max())); + // we use only boundary_ids in the range from 0 to + // numbers::internal_face_boundary_id-1 + AssertIndexRange(material_id, + numbers::internal_face_boundary_id); + + subcelldata.boundary_lines.back().boundary_id = + static_cast(material_id); + + // transform from ucd to + // consecutive numbering + for (unsigned int &vertex : + subcelldata.boundary_lines.back().vertices) + if (vertex_indices.find(vertex) != vertex_indices.end()) + // vertex with this index exists + vertex = vertex_indices[vertex]; + else + { + // no such vertex index + AssertThrow(false, + ExcInvalidVertexIndex(cell_per_entity, + vertex)); + vertex = numbers::invalid_unsigned_int; + } } - } - else if ((cell_type == 2 || cell_type == 3) && - (dim == 3)) // a triangle or a quad in 3d - // boundary info - { - unsigned int vertices_per_cell = 0; - // check cell type - if (cell_type == 2) // tri - vertices_per_cell = 3; - else if (cell_type == 3) // quad - vertices_per_cell = 4; - - subcelldata.boundary_quads.emplace_back(); - - // resize vertices - subcelldata.boundary_quads.back().vertices.resize( - vertices_per_cell); - // for loop - for (unsigned int i = 0; i < vertices_per_cell; ++i) - in >> subcelldata.boundary_quads.back().vertices[i]; - - // to make sure that the cast won't fail - Assert(material_id <= - std::numeric_limits::max(), - ExcIndexRange( - material_id, - 0, - std::numeric_limits::max())); - // we use only boundary_ids in the range from 0 to - // numbers::internal_face_boundary_id-1 - AssertIndexRange(material_id, - numbers::internal_face_boundary_id); - - subcelldata.boundary_quads.back().boundary_id = - static_cast(material_id); - - // transform from gmsh to - // consecutive numbering - for (unsigned int &vertex : - subcelldata.boundary_quads.back().vertices) - if (vertex_indices.find(vertex) != vertex_indices.end()) - // vertex with this index exists - vertex = vertex_indices[vertex]; - else + else if ((cell_type == 2 || cell_type == 3) && + (dim == 3)) // a triangle or a quad in 3d + // boundary info { - // no such vertex index - Assert(false, - ExcInvalidVertexIndex(cell_per_entity, vertex)); - vertex = numbers::invalid_unsigned_int; + unsigned int vertices_per_cell = 0; + // check cell type + if (cell_type == 2) // tri + vertices_per_cell = 3; + else if (cell_type == 3) // quad + vertices_per_cell = 4; + + subcelldata.boundary_quads.emplace_back(); + + // resize vertices + subcelldata.boundary_quads.back().vertices.resize( + vertices_per_cell); + // for loop + for (unsigned int i = 0; i < vertices_per_cell; ++i) + in >> subcelldata.boundary_quads.back().vertices[i]; + + // to make sure that the cast won't fail + Assert(material_id <= + std::numeric_limits::max(), + ExcIndexRange( + material_id, + 0, + std::numeric_limits::max())); + // we use only boundary_ids in the range from 0 to + // numbers::internal_face_boundary_id-1 + AssertIndexRange(material_id, + numbers::internal_face_boundary_id); + + subcelldata.boundary_quads.back().boundary_id = + static_cast(material_id); + + // transform from gmsh to + // consecutive numbering + for (unsigned int &vertex : + subcelldata.boundary_quads.back().vertices) + if (vertex_indices.find(vertex) != vertex_indices.end()) + // vertex with this index exists + vertex = vertex_indices[vertex]; + else + { + // no such vertex index + Assert(false, + ExcInvalidVertexIndex(cell_per_entity, vertex)); + vertex = numbers::invalid_unsigned_int; + } } - } - else if (cell_type == 15) - { - // read the indices of nodes given - unsigned int node_index = 0; - if (gmsh_file_format < 20) - { - // For points (cell_type==15), we can only ever - // list one node index. - AssertThrow(nod_num == 1, ExcInternalError()); - in >> node_index; - } - else - { - in >> node_index; - } + else if (cell_type == 15) + { + // read the indices of nodes given + unsigned int node_index = 0; + if (gmsh_file_format < 20) + { + // For points (cell_type==15), we can only ever + // list one node index. + AssertThrow(nod_num == 1, ExcInternalError()); + in >> node_index; + } + else + { + in >> node_index; + } - // we only care about boundary indicators assigned to - // individual vertices in 1d (because otherwise the vertices - // are not faces) - if (dim == 1) - boundary_ids_1d[vertex_indices[node_index]] = material_id; - } - else - { - AssertThrow(false, ExcGmshUnsupportedGeometry(cell_type)); - } - } + // we only care about boundary indicators assigned to + // individual vertices in 1d (because otherwise the vertices + // are not faces) + if (dim == 1) + boundary_ids_1d[vertex_indices[node_index]] = material_id; + } + else + { + AssertThrow(false, ExcGmshUnsupportedGeometry(cell_type)); + } + } + } + AssertDimension(global_cell, n_cells); + } + // Assert that we reached the end of the block + in >> line; + static const std::string end_elements_marker[] = {"$ENDELM", "$EndElements"}; + AssertThrow(line == end_elements_marker[gmsh_file_format == 10 ? 0 : 1], + ExcInvalidGMSHInput(line)); + + } - AssertDimension(global_cell, n_cells); + // skip the line if not recognized + in >> line; } - // Assert that we reached the end of the block - in >> line; - static const std::string end_elements_marker[] = {"$ENDELM", "$EndElements"}; - AssertThrow(line == end_elements_marker[gmsh_file_format == 10 ? 0 : 1], - ExcInvalidGMSHInput(line)); - AssertThrow(in.fail() == false, ExcIO()); - + // check that we actually read some cells. AssertThrow(cells.size() > 0, ExcGmshNoCellInformation(subcelldata.boundary_lines.size(), diff --git a/tests/grid/grid_in_msh_version_2.cc b/tests/grid/grid_in_msh_version_2.cc index 59c362190f..74e7b868e1 100644 --- a/tests/grid/grid_in_msh_version_2.cc +++ b/tests/grid/grid_in_msh_version_2.cc @@ -51,6 +51,8 @@ filename_resolution() GridIn<2>::msh); check_file<2>(std::string(SOURCE_DIR "/grid_in_msh_version_2/hole8170.msh"), GridIn<2>::msh); + check_file<2>(std::string(SOURCE_DIR "/grid_in_msh_version_2/hole81_commented.msh"), + GridIn<2>::msh); } diff --git a/tests/grid/grid_in_msh_version_2.output b/tests/grid/grid_in_msh_version_2.output index e9b06639de..ff350da1de 100644 --- a/tests/grid/grid_in_msh_version_2.output +++ b/tests/grid/grid_in_msh_version_2.output @@ -57758,3 +57758,571 @@ DEAL:: 8449 8170 0.25961 0.84911 0 1 +DEAL:: 109 81 +0.20547 0.33579 0 1 +0.26169 0.25359 0 1 +0.32322 0.32322 0 1 +0.26903 0.40433 0 1 +0.20547 0.33579 0 1 + + +0.16890 0.42183 0 1 +0.20547 0.33579 0 1 +0.26903 0.40433 0 1 +0.25000 0.50000 0 1 +0.16890 0.42183 0 1 + + +0.35191 0.20190 0 1 +0.40433 0.26903 0 1 +0.32322 0.32322 0 1 +0.26169 0.25359 0 1 +0.35191 0.20190 0 1 + + +0.26745 0.11592 0 1 +0.35191 0.20190 0 1 +0.26169 0.25359 0 1 +0.15106 0.14697 0 1 +0.26745 0.11592 0 1 + + +0.11560 0.25755 0 1 +0.15106 0.14697 0 1 +0.26169 0.25359 0 1 +0.20547 0.33579 0 1 +0.11560 0.25755 0 1 + + +0.092879 0.35590 0 1 +0.11560 0.25755 0 1 +0.20547 0.33579 0 1 +0.16890 0.42183 0 1 +0.092879 0.35590 0 1 + + +0.47567 0.17307 0 1 +0.50000 0.25000 0 1 +0.40433 0.26903 0 1 +0.35191 0.20190 0 1 +0.47567 0.17307 0 1 + + +0.0000 0.20000 0 1 +0.0000 0.10000 0 1 +0.15106 0.14697 0 1 +0.11560 0.25755 0 1 +0.0000 0.20000 0 1 + + +0.0000 0.30000 0 1 +0.0000 0.20000 0 1 +0.11560 0.25755 0 1 +0.092879 0.35590 0 1 +0.0000 0.30000 0 1 + + +0.0000 0.40000 0 1 +0.0000 0.30000 0 1 +0.092879 0.35590 0 1 +0.070625 0.44267 0 1 +0.0000 0.40000 0 1 + + +0.070625 0.44267 0 1 +0.092879 0.35590 0 1 +0.16890 0.42183 0 1 +0.12425 0.50205 0 1 +0.070625 0.44267 0 1 + + +0.12425 0.50205 0 1 +0.16890 0.42183 0 1 +0.25000 0.50000 0 1 +0.17004 0.57976 0 1 +0.12425 0.50205 0 1 + + +0.21001 0.66461 0 1 +0.17004 0.57976 0 1 +0.25000 0.50000 0 1 +0.26903 0.59567 0 1 +0.21001 0.66461 0 1 + + +0.27365 0.72750 0 1 +0.21001 0.66461 0 1 +0.26903 0.59567 0 1 +0.32322 0.67678 0 1 +0.27365 0.72750 0 1 + + +0.33293 0.76189 0 1 +0.27365 0.72750 0 1 +0.32322 0.67678 0 1 +0.40433 0.73097 0 1 +0.33293 0.76189 0 1 + + +0.50692 0.82456 0 1 +0.40960 0.82316 0 1 +0.40433 0.73097 0 1 +0.50000 0.75000 0 1 +0.50692 0.82456 0 1 + + +0.60494 0.81159 0 1 +0.50692 0.82456 0 1 +0.50000 0.75000 0 1 +0.59567 0.73097 0 1 +0.60494 0.81159 0 1 + + +0.70756 0.77928 0 1 +0.60494 0.81159 0 1 +0.59567 0.73097 0 1 +0.67678 0.67678 0 1 +0.70756 0.77928 0 1 + + +0.0000 0.50000 0 1 +0.0000 0.40000 0 1 +0.070625 0.44267 0 1 +0.047080 0.50155 0 1 +0.0000 0.50000 0 1 + + +0.047080 0.50155 0 1 +0.070625 0.44267 0 1 +0.12425 0.50205 0 1 +0.070925 0.56122 0 1 +0.047080 0.50155 0 1 + + +0.070925 0.56122 0 1 +0.12425 0.50205 0 1 +0.17004 0.57976 0 1 +0.096134 0.64928 0 1 +0.070925 0.56122 0 1 + + +0.12988 0.76254 0 1 +0.096134 0.64928 0 1 +0.17004 0.57976 0 1 +0.21001 0.66461 0 1 +0.12988 0.76254 0 1 + + +0.23514 0.80526 0 1 +0.12988 0.76254 0 1 +0.21001 0.66461 0 1 +0.27365 0.72750 0 1 +0.23514 0.80526 0 1 + + +0.32220 0.82686 0 1 +0.23514 0.80526 0 1 +0.27365 0.72750 0 1 +0.33293 0.76189 0 1 +0.32220 0.82686 0 1 + + +0.30000 1.0000 0 1 +0.30978 0.91020 0 1 +0.40726 0.91083 0 1 +0.40000 1.0000 0 1 +0.30000 1.0000 0 1 + + +0.50489 0.91016 0 1 +0.40726 0.91083 0 1 +0.40960 0.82316 0 1 +0.50692 0.82456 0 1 +0.50489 0.91016 0 1 + + +0.60598 0.90353 0 1 +0.50489 0.91016 0 1 +0.50692 0.82456 0 1 +0.60494 0.81159 0 1 +0.60598 0.90353 0 1 + + +0.71419 0.88840 0 1 +0.60598 0.90353 0 1 +0.60494 0.81159 0 1 +0.70756 0.77928 0 1 +0.71419 0.88840 0 1 + + +0.50000 1.0000 0 1 +0.40000 1.0000 0 1 +0.40726 0.91083 0 1 +0.50489 0.91016 0 1 +0.50000 1.0000 0 1 + + +0.30000 1.0000 0 1 +0.20000 1.0000 0 1 +0.21346 0.90127 0 1 +0.30978 0.91020 0 1 +0.30000 1.0000 0 1 + + +0.20000 1.0000 0 1 +0.10000 1.0000 0 1 +0.10985 0.89198 0 1 +0.21346 0.90127 0 1 +0.20000 1.0000 0 1 + + +0.0000 1.0000 0 1 +0.0000 0.90000 0 1 +0.10985 0.89198 0 1 +0.10000 1.0000 0 1 +0.0000 1.0000 0 1 + + +0.60000 1.0000 0 1 +0.50000 1.0000 0 1 +0.50489 0.91016 0 1 +0.60598 0.90353 0 1 +0.60000 1.0000 0 1 + + +0.0000 0.90000 0 1 +0.0000 0.80000 0 1 +0.12988 0.76254 0 1 +0.10985 0.89198 0 1 +0.0000 0.90000 0 1 + + +0.0000 0.80000 0 1 +0.0000 0.70000 0 1 +0.096134 0.64928 0 1 +0.12988 0.76254 0 1 +0.0000 0.80000 0 1 + + +0.0000 0.60000 0 1 +0.070925 0.56122 0 1 +0.096134 0.64928 0 1 +0.0000 0.70000 0 1 +0.0000 0.60000 0 1 + + +0.0000 0.50000 0 1 +0.047080 0.50155 0 1 +0.070925 0.56122 0 1 +0.0000 0.60000 0 1 +0.0000 0.50000 0 1 + + +0.21346 0.90127 0 1 +0.10985 0.89198 0 1 +0.12988 0.76254 0 1 +0.23514 0.80526 0 1 +0.21346 0.90127 0 1 + + +0.30978 0.91020 0 1 +0.21346 0.90127 0 1 +0.23514 0.80526 0 1 +0.32220 0.82686 0 1 +0.30978 0.91020 0 1 + + +0.82051 0.73439 0 1 +0.70756 0.77928 0 1 +0.67678 0.67678 0 1 +0.73097 0.59567 0 1 +0.82051 0.73439 0 1 + + +1.0000 0.70000 0 1 +0.82051 0.73439 0 1 +0.73097 0.59567 0 1 +1.0000 0.60000 0 1 +1.0000 0.70000 0 1 + + +0.85111 0.48662 0 1 +1.0000 0.60000 0 1 +0.73097 0.59567 0 1 +0.75000 0.50000 0 1 +0.85111 0.48662 0 1 + + +0.81603 0.38904 0 1 +0.85111 0.48662 0 1 +0.75000 0.50000 0 1 +0.73097 0.40433 0 1 +0.81603 0.38904 0 1 + + +0.84281 0.86350 0 1 +0.71419 0.88840 0 1 +0.70756 0.77928 0 1 +0.82051 0.73439 0 1 +0.84281 0.86350 0 1 + + +0.90941 0.38287 0 1 +0.92017 0.45650 0 1 +0.85111 0.48662 0 1 +0.81603 0.38904 0 1 +0.90941 0.38287 0 1 + + +1.0000 0.90000 0 1 +0.91427 0.92117 0 1 +0.84281 0.86350 0 1 +1.0000 0.80000 0 1 +1.0000 0.90000 0 1 + + +1.0000 1.0000 0 1 +0.90000 1.0000 0 1 +0.91427 0.92117 0 1 +1.0000 0.90000 0 1 +1.0000 1.0000 0 1 + + +0.90941 0.38287 0 1 +1.0000 0.40000 0 1 +1.0000 0.50000 0 1 +0.92017 0.45650 0 1 +0.90941 0.38287 0 1 + + +0.90000 1.0000 0 1 +0.80000 1.0000 0 1 +0.84281 0.86350 0 1 +0.91427 0.92117 0 1 +0.90000 1.0000 0 1 + + +0.80000 1.0000 0 1 +0.70000 1.0000 0 1 +0.71419 0.88840 0 1 +0.84281 0.86350 0 1 +0.80000 1.0000 0 1 + + +0.70000 1.0000 0 1 +0.60000 1.0000 0 1 +0.60598 0.90353 0 1 +0.71419 0.88840 0 1 +0.70000 1.0000 0 1 + + +0.89263 0.28721 0 1 +1.0000 0.30000 0 1 +1.0000 0.40000 0 1 +0.90941 0.38287 0 1 +0.89263 0.28721 0 1 + + +0.62245 0.20943 0 1 +0.75883 0.27502 0 1 +0.67678 0.32322 0 1 +0.59567 0.26903 0 1 +0.62245 0.20943 0 1 + + +0.47567 0.17307 0 1 +0.62245 0.20943 0 1 +0.59567 0.26903 0 1 +0.50000 0.25000 0 1 +0.47567 0.17307 0 1 + + +0.57957 0.061246 0 1 +0.60000 0.0000 0 1 +0.70000 0.0000 0 1 +0.67900 0.074028 0 1 +0.57957 0.061246 0 1 + + +0.67900 0.074028 0 1 +0.70000 0.0000 0 1 +0.80000 0.0000 0 1 +0.78694 0.086571 0 1 +0.67900 0.074028 0 1 + + +0.89357 0.094341 0 1 +0.78694 0.086571 0 1 +0.80000 0.0000 0 1 +0.90000 0.0000 0 1 +0.89357 0.094341 0 1 + + +0.89357 0.094341 0 1 +0.90000 0.0000 0 1 +1.0000 0.0000 0 1 +1.0000 0.10000 0 1 +0.89357 0.094341 0 1 + + +0.88980 0.18952 0 1 +0.89357 0.094341 0 1 +1.0000 0.10000 0 1 +1.0000 0.20000 0 1 +0.88980 0.18952 0 1 + + +0.89263 0.28721 0 1 +0.88980 0.18952 0 1 +1.0000 0.20000 0 1 +1.0000 0.30000 0 1 +0.89263 0.28721 0 1 + + +0.88980 0.18952 0 1 +0.77200 0.17509 0 1 +0.78694 0.086571 0 1 +0.89357 0.094341 0 1 +0.88980 0.18952 0 1 + + +1.0000 0.80000 0 1 +0.84281 0.86350 0 1 +0.82051 0.73439 0 1 +1.0000 0.70000 0 1 +1.0000 0.80000 0 1 + + +0.92017 0.45650 0 1 +1.0000 0.50000 0 1 +1.0000 0.60000 0 1 +0.85111 0.48662 0 1 +0.92017 0.45650 0 1 + + +0.75883 0.27502 0 1 +0.81603 0.38904 0 1 +0.73097 0.40433 0 1 +0.67678 0.32322 0 1 +0.75883 0.27502 0 1 + + +0.65340 0.14414 0 1 +0.67900 0.074028 0 1 +0.78694 0.086571 0 1 +0.77200 0.17509 0 1 +0.65340 0.14414 0 1 + + +0.62245 0.20943 0 1 +0.65340 0.14414 0 1 +0.77200 0.17509 0 1 +0.75883 0.27502 0 1 +0.62245 0.20943 0 1 + + +0.47567 0.17307 0 1 +0.54011 0.11687 0 1 +0.65340 0.14414 0 1 +0.62245 0.20943 0 1 +0.47567 0.17307 0 1 + + +0.54011 0.11687 0 1 +0.57957 0.061246 0 1 +0.67900 0.074028 0 1 +0.65340 0.14414 0 1 +0.54011 0.11687 0 1 + + +0.89263 0.28721 0 1 +0.90941 0.38287 0 1 +0.81603 0.38904 0 1 +0.75883 0.27502 0 1 +0.89263 0.28721 0 1 + + +0.89263 0.28721 0 1 +0.75883 0.27502 0 1 +0.77200 0.17509 0 1 +0.88980 0.18952 0 1 +0.89263 0.28721 0 1 + + +0.37351 0.096370 0 1 +0.30000 0.0000 0 1 +0.40000 0.0000 0 1 +0.45966 0.072308 0 1 +0.37351 0.096370 0 1 + + +0.37351 0.096370 0 1 +0.26745 0.11592 0 1 +0.20000 0.0000 0 1 +0.30000 0.0000 0 1 +0.37351 0.096370 0 1 + + +0.51267 0.044159 0 1 +0.50000 0.0000 0 1 +0.60000 0.0000 0 1 +0.57957 0.061246 0 1 +0.51267 0.044159 0 1 + + +0.45966 0.072308 0 1 +0.40000 0.0000 0 1 +0.50000 0.0000 0 1 +0.51267 0.044159 0 1 +0.45966 0.072308 0 1 + + +0.45966 0.072308 0 1 +0.51267 0.044159 0 1 +0.57957 0.061246 0 1 +0.54011 0.11687 0 1 +0.45966 0.072308 0 1 + + +0.37351 0.096370 0 1 +0.47567 0.17307 0 1 +0.35191 0.20190 0 1 +0.26745 0.11592 0 1 +0.37351 0.096370 0 1 + + +0.37351 0.096370 0 1 +0.45966 0.072308 0 1 +0.54011 0.11687 0 1 +0.47567 0.17307 0 1 +0.37351 0.096370 0 1 + + +0.30978 0.91020 0 1 +0.32220 0.82686 0 1 +0.40960 0.82316 0 1 +0.40726 0.91083 0 1 +0.30978 0.91020 0 1 + + +0.32220 0.82686 0 1 +0.33293 0.76189 0 1 +0.40433 0.73097 0 1 +0.40960 0.82316 0 1 +0.32220 0.82686 0 1 + + +0.26745 0.11592 0 1 +0.15106 0.14697 0 1 +0.10000 0.0000 0 1 +0.20000 0.0000 0 1 +0.26745 0.11592 0 1 + + +0.0000 0.10000 0 1 +0.0000 0.0000 0 1 +0.10000 0.0000 0 1 +0.15106 0.14697 0 1 +0.0000 0.10000 0 1 + + diff --git a/tests/grid/grid_in_msh_version_2/hole81_commented.msh b/tests/grid/grid_in_msh_version_2/hole81_commented.msh new file mode 100644 index 0000000000..6484b8b305 --- /dev/null +++ b/tests/grid/grid_in_msh_version_2/hole81_commented.msh @@ -0,0 +1,217 @@ +$Comments +------------------------------------------------------- +File exported by XYZ +------------------------------------------------------- +$EndComments +$MeshFormat +2.0 0 8 +$EndMeshFormat +$Comments +------------------------------------------------------- +NumNodes +Id x y z +------------------------------------------------------- +$EndComments +$Nodes +109 +1 0.323223 0.323223 0 +2 0.26903 0.404329 0 +3 0.205471 0.335785 0 +4 0.261688 0.253588 0 +5 0.25 0.5 0 +6 0.168897 0.421825 0 +7 0.351915 0.201897 0 +8 0.404329 0.26903 0 +9 0.151061 0.146973 0 +10 0.267453 0.115922 0 +11 0.115602 0.257549 0 +12 0.0928787 0.355903 0 +13 0.475666 0.173068 0 +14 0.5 0.25 0 +15 0 0.2 0 +16 0 0.1 0 +17 0 0.3 0 +18 0.0706253 0.44267 0 +19 0 0.4 0 +20 0.12425 0.502054 0 +21 0.170037 0.579759 0 +22 0.26903 0.595671 0 +23 0.210011 0.664612 0 +24 0.323223 0.676777 0 +25 0.273645 0.727505 0 +26 0.404329 0.73097 0 +27 0.332934 0.761889 0 +28 0.5 0.75 0 +29 0.506916 0.824562 0 +30 0.409603 0.823159 0 +31 0.595671 0.73097 0 +32 0.604938 0.811586 0 +33 0.676777 0.676777 0 +34 0.707555 0.779276 0 +35 0.0470804 0.501553 0 +36 0 0.5 0 +37 0.0709247 0.561217 0 +38 0.0961342 0.649279 0 +39 0.129878 0.762537 0 +40 0.235142 0.805255 0 +41 0.322202 0.826855 0 +42 0.407261 0.910834 0 +43 0.4 1 0 +44 0.3 1 0 +45 0.309776 0.910203 0 +46 0.504887 0.910157 0 +47 0.605977 0.903532 0 +48 0.714192 0.888402 0 +49 0.5 1 0 +50 0.2 1 0 +51 0.213462 0.90127 0 +52 0.1 1 0 +53 0.109849 0.891977 0 +54 0 1 0 +55 0 0.9 0 +56 0.6 1 0 +57 0 0.8 0 +58 0 0.7 0 +59 0 0.6 0 +60 0.73097 0.595671 0 +61 0.820507 0.734392 0 +62 1 0.6 0 +63 1 0.7 0 +64 0.75 0.5 0 +65 0.851105 0.486625 0 +66 0.73097 0.404329 0 +67 0.816035 0.389041 0 +68 0.842812 0.863497 0 +69 0.909413 0.382874 0 +70 0.920173 0.4565 0 +71 1 0.8 0 +72 1 0.9 0 +73 0.914271 0.921166 0 +74 1 1 0 +75 0.9 1 0 +76 1 0.4 0 +77 1 0.5 0 +78 0.8 1 0 +79 0.7 1 0 +80 0.892629 0.287213 0 +81 1 0.3 0 +82 0.676777 0.323223 0 +83 0.595671 0.26903 0 +84 0.622453 0.20943 0 +85 0.758832 0.275022 0 +86 0.6 0 0 +87 0.7 0 0 +88 0.679002 0.0740277 0 +89 0.579569 0.0612463 0 +90 0.8 0 0 +91 0.786941 0.0865707 0 +92 0.9 0 0 +93 0.893574 0.0943409 0 +94 1 0 0 +95 1 0.1 0 +96 1 0.2 0 +97 0.889796 0.189523 0 +98 0.772001 0.17509 0 +99 0.653399 0.144139 0 +100 0.540113 0.116873 0 +101 0.4 0 0 +102 0.459658 0.0723085 0 +103 0.373509 0.0963702 0 +104 0.3 0 0 +105 0.2 0 0 +106 0.5 0 0 +107 0.512671 0.0441593 0 +108 0.1 0 0 +109 0 0 0 +$EndNodes +$Comments +------------------------------------------------------- +Types: bar(1) tri(2) quad(3) tet(4) hex(5) prism(6) pyramid(7) +NumElements +Id Type NumTags PhysGrp ElemGrp IndexList +------------------------------------------------------- +$EndComments +$Elements +81 +81 3 1 1 1 2 3 4 +80 3 1 1 2 5 6 3 +79 3 1 1 1 4 7 8 +78 3 1 1 4 9 10 7 +77 3 1 1 4 3 11 9 +76 3 1 1 3 6 12 11 +75 3 1 1 8 7 13 14 +74 3 1 1 9 11 15 16 +73 3 1 1 11 12 17 15 +72 3 1 1 17 12 18 19 +71 3 1 1 12 6 20 18 +70 3 1 1 6 5 21 20 +69 3 1 1 5 22 23 21 +68 3 1 1 22 24 25 23 +67 3 1 1 24 26 27 25 +66 3 1 1 26 28 29 30 +65 3 1 1 28 31 32 29 +64 3 1 1 31 33 34 32 +63 3 1 1 19 18 35 36 +62 3 1 1 18 20 37 35 +61 3 1 1 20 21 38 37 +60 3 1 1 21 23 39 38 +59 3 1 1 23 25 40 39 +58 3 1 1 25 27 41 40 +57 3 1 1 42 43 44 45 +56 3 1 1 30 29 46 42 +55 3 1 1 29 32 47 46 +54 3 1 1 32 34 48 47 +53 3 1 1 42 46 49 43 +52 3 1 1 44 50 51 45 +51 3 1 1 50 52 53 51 +50 3 1 1 52 54 55 53 +49 3 1 1 46 47 56 49 +48 3 1 1 55 57 39 53 +47 3 1 1 57 58 38 39 +46 3 1 1 58 59 37 38 +45 3 1 1 37 59 36 35 +44 3 1 1 51 53 39 40 +43 3 1 1 40 41 45 51 +42 3 1 1 33 60 61 34 +41 3 1 1 60 62 63 61 +40 3 1 1 60 64 65 62 +39 3 1 1 64 66 67 65 +38 3 1 1 34 61 68 48 +37 3 1 1 65 67 69 70 +36 3 1 1 71 72 73 68 +35 3 1 1 72 74 75 73 +34 3 1 1 70 69 76 77 +33 3 1 1 75 78 68 73 +32 3 1 1 78 79 48 68 +31 3 1 1 48 79 56 47 +30 3 1 1 76 69 80 81 +29 3 1 1 82 83 84 85 +28 3 1 1 83 14 13 84 +27 3 1 1 86 87 88 89 +26 3 1 1 87 90 91 88 +25 3 1 1 90 92 93 91 +24 3 1 1 92 94 95 93 +23 3 1 1 95 96 97 93 +22 3 1 1 96 81 80 97 +21 3 1 1 91 93 97 98 +20 3 1 1 68 61 63 71 +19 3 1 1 62 65 70 77 +18 3 1 1 67 66 82 85 +17 3 1 1 91 98 99 88 +16 3 1 1 98 85 84 99 +15 3 1 1 99 84 13 100 +14 3 1 1 88 99 100 89 +13 3 1 1 80 69 67 85 +12 3 1 1 80 85 98 97 +11 3 1 1 101 102 103 104 +10 3 1 1 104 103 10 105 +9 3 1 1 106 86 89 107 +8 3 1 1 101 106 107 102 +7 3 1 1 107 89 100 102 +6 3 1 1 10 103 13 7 +5 3 1 1 103 102 100 13 +4 3 1 1 45 41 30 42 +3 3 1 1 41 27 26 30 +2 3 1 1 108 105 10 9 +1 3 1 1 9 16 109 108 +$EndElements -- 2.39.5