diff --git a/src/mesh/GmshReader.cpp b/src/mesh/GmshReader.cpp index 71284c13ed850c7ecd0da8c3bcab938ef93c5bf9..9866c9ad8389c509565ba26c8e94915c3f597931 100644 --- a/src/mesh/GmshReader.cpp +++ b/src/mesh/GmshReader.cpp @@ -237,6 +237,7 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData& descriptor.cell_type_vector[j] = CellType::Tetrahedron; descriptor.cell_number_vector[j] = gmsh_data.__tetrahedra_number[j]; } + const size_t nb_hexahedra = gmsh_data.__hexahedra.size(); for (size_t j = 0; j < nb_hexahedra; ++j) { const size_t jh = nb_tetrahedra + j; @@ -248,6 +249,17 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData& descriptor.cell_number_vector[jh] = gmsh_data.__hexahedra_number[j]; } + const size_t nb_pyramids = gmsh_data.__pyramids.size(); + for (size_t j = 0; j < nb_pyramids; ++j) { + const size_t jh = nb_tetrahedra + nb_hexahedra + j; + descriptor.cell_to_node_vector[jh].resize(5); + for (int r = 0; r < 5; ++r) { + descriptor.cell_to_node_vector[jh][r] = gmsh_data.__pyramids[j][r]; + } + descriptor.cell_type_vector[jh] = CellType::Pyramid; + descriptor.cell_number_vector[jh] = gmsh_data.__pyramids_number[j]; + } + std::map<unsigned int, std::vector<unsigned int>> ref_cells_map; for (unsigned int r = 0; r < gmsh_data.__tetrahedra_ref.size(); ++r) { const unsigned int elem_number = gmsh_data.__tetrahedra_ref[r]; @@ -261,6 +273,12 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData& ref_cells_map[ref].push_back(elem_number); } + for (unsigned int j = 0; j < gmsh_data.__pyramids_ref.size(); ++j) { + const size_t elem_number = nb_tetrahedra + nb_hexahedra + j; + const unsigned int& ref = gmsh_data.__pyramids_ref[j]; + ref_cells_map[ref].push_back(elem_number); + } + for (const auto& ref_cell_list : ref_cells_map) { Array<CellId> cell_list(ref_cell_list.second.size()); for (size_t j = 0; j < ref_cell_list.second.size(); ++j) { @@ -309,8 +327,8 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData& const unsigned int& ref = gmsh_data.__triangles_ref[f]; ref_faces_map[ref].push_back(face_id); - if (descriptor.face_number_vector[face_id] != gmsh_data.__quadrangles_number[f]) { - if (auto i_face = face_number_id_map.find(gmsh_data.__quadrangles_number[f]); + if (descriptor.face_number_vector[face_id] != gmsh_data.__triangles_number[f]) { + if (auto i_face = face_number_id_map.find(gmsh_data.__triangles_number[f]); i_face != face_number_id_map.end()) { const int other_face_id = i_face->second; std::swap(descriptor.face_number_vector[face_id], descriptor.face_number_vector[other_face_id]); @@ -322,7 +340,7 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData& face_number_id_map[descriptor.face_number_vector[other_face_id]] = other_face_id; } else { face_number_id_map.erase(descriptor.face_number_vector[face_id]); - descriptor.face_number_vector[face_id] = gmsh_data.__quadrangles_number[f]; + descriptor.face_number_vector[face_id] = gmsh_data.__triangles_number[f]; face_number_id_map[descriptor.face_number_vector[face_id]] = face_id; } } @@ -517,7 +535,7 @@ GmshReader::GmshReader(const std::string& filename) : m_filename(filename) __primitivesNames[5] = "prisms"; __supportedPrimitives[5] = false; __primitivesNames[6] = "pyramids"; - __supportedPrimitives[6] = false; + __supportedPrimitives[6] = true; __primitivesNames[7] = "second order edges"; __supportedPrimitives[7] = false; __primitivesNames[8] = "second order triangles"; @@ -1009,6 +1027,12 @@ GmshReader::__proceedData() m_mesh_data.__hexahedra_number.resize(elementNumber[i]); break; } + case 6: { // pyramid + m_mesh_data.__pyramids = Array<GmshData::Pyramid>(elementNumber[i]); + m_mesh_data.__pyramids_ref.resize(elementNumber[i]); + m_mesh_data.__pyramids_number.resize(elementNumber[i]); + break; + } case 14: { // point m_mesh_data.__points = Array<GmshData::Point>(elementNumber[i]); m_mesh_data.__points_ref.resize(elementNumber[i]); @@ -1017,7 +1041,6 @@ GmshReader::__proceedData() } // Unsupported entities case 5: // prism - case 6: // pyramid case 7: // second order edge case 8: // second order triangle case 9: // second order quadrangle @@ -1146,6 +1169,24 @@ GmshReader::__proceedData() m_mesh_data.__hexahedra_number[hexahedronNumber] = m_mesh_data.__elementNumber[i]; hexahedronNumber++; // one more hexahedron break; + } + case 6: { // pyramid + int& pyramid_number = elementNumber[6]; + + const int a = m_mesh_data.__verticesCorrepondance[elementVertices[0]]; + const int b = m_mesh_data.__verticesCorrepondance[elementVertices[1]]; + const int c = m_mesh_data.__verticesCorrepondance[elementVertices[2]]; + const int d = m_mesh_data.__verticesCorrepondance[elementVertices[3]]; + const int e = m_mesh_data.__verticesCorrepondance[elementVertices[4]]; + if ((a < 0) or (b < 0) or (c < 0) or (d < 0) or (e < 0)) { + throw NormalError("reading file '" + m_filename + "': error reading element " + std::to_string(i) + + " [bad vertices definition]"); + } + m_mesh_data.__pyramids[pyramid_number] = GmshData::Pyramid(d, c, b, a, e); + m_mesh_data.__pyramids_ref[pyramid_number] = m_mesh_data.__references[i]; + m_mesh_data.__pyramids_number[pyramid_number] = m_mesh_data.__elementNumber[i]; + pyramid_number++; // one more hexahedron + break; } // Unsupported entities case 14: { // point @@ -1158,7 +1199,6 @@ GmshReader::__proceedData() break; } case 5: // prism - case 6: // pyramid case 7: // second order edge case 8: // second order triangle case 9: // second order quadrangle diff --git a/src/mesh/GmshReader.hpp b/src/mesh/GmshReader.hpp index ae18d30ca3acd989f7a48d610afb458de8c4342d..36e21508f10eda492ec1128816215d5a4d8c624d 100644 --- a/src/mesh/GmshReader.hpp +++ b/src/mesh/GmshReader.hpp @@ -104,6 +104,11 @@ class GmshReader : public MeshBuilderBase std::vector<int> __tetrahedra_ref; std::vector<int> __tetrahedra_number; + using Pyramid = TinyVector<5, unsigned int>; + Array<Pyramid> __pyramids; + std::vector<int> __pyramids_ref; + std::vector<int> __pyramids_number; + using Hexahedron = TinyVector<8, unsigned int>; Array<Hexahedron> __hexahedra; std::vector<int> __hexahedra_ref;