Skip to content
Snippets Groups Projects

Compare revisions

Changes are shown as if the source revision was being merged into the target revision. Learn more about comparing revisions.

Source

Select target project
No results found
Select Git revision
  • Nodal_diffusion
  • develop
  • feature/Navier-Stokes
  • feature/Nodal_diffusion
  • feature/composite-scheme
  • feature/composite-scheme-other-fluxes
  • feature/composite-scheme-sources
  • feature/coupling_module
  • feature/discontinuous-galerkin
  • feature/escobar-smoother
  • feature/explicit-gp-cfl
  • feature/gks
  • feature/gmsh-reader
  • feature/hypoelasticity
  • feature/hypoelasticity-clean
  • feature/implicit-solver
  • feature/implicit-solver-o2
  • feature/kinetic-schemes
  • feature/local-dt-fsi
  • feature/merge-local-dt-fsi
  • feature/polynomials
  • feature/reconstruction
  • feature/serraille
  • feature/variational-hydro
  • hyperplastic
  • master
  • navier-stokes
  • origin/stage/bouguettaia
  • save_clemence
  • test/voronoi1d
  • Kidder
  • v0
  • v0.0.1
  • v0.0.2
  • v0.0.3
  • v0.0.4
  • v0.1.0
  • v0.2.0
  • v0.3.0
  • v0.4.0
  • v0.4.1
  • v0.5.0
42 results

Target

Select target project
  • code/pugs
1 result
Select Git revision
  • Nodal_diffusion
  • develop
  • feature/Navier-Stokes
  • feature/Nodal_diffusion
  • feature/composite-scheme
  • feature/composite-scheme-other-fluxes
  • feature/composite-scheme-sources
  • feature/coupling_module
  • feature/discontinuous-galerkin
  • feature/escobar-smoother
  • feature/explicit-gp-cfl
  • feature/gks
  • feature/gmsh-reader
  • feature/hypoelasticity
  • feature/hypoelasticity-clean
  • feature/implicit-solver
  • feature/implicit-solver-o2
  • feature/kinetic-schemes
  • feature/local-dt-fsi
  • feature/merge-local-dt-fsi
  • feature/polynomials
  • feature/reconstruction
  • feature/serraille
  • feature/variational-hydro
  • hyperplastic
  • master
  • navier-stokes
  • origin/stage/bouguettaia
  • save_clemence
  • test/voronoi1d
  • Kidder
  • v0
  • v0.0.1
  • v0.0.2
  • v0.0.3
  • v0.0.4
  • v0.1.0
  • v0.2.0
  • v0.3.0
  • v0.4.0
  • v0.4.1
  • v0.5.0
42 results
Show changes

Commits on Source 2

Show whitespace changes
Inline Side-by-side
src/mesh/GmshReader.cpp
View file @ dbc6302d
......@@ -16,10 +16,9 @@
#include <pegtl.hpp>
#include <rang.hpp>
#include <fstream>
#include <charconv>
#include <iostream>
#include <map>
#include <regex>
#include <set>
#include <sstream>
#include <unordered_map>
......@@ -36,17 +35,18 @@ GmshConnectivityBuilder<1>::GmshConnectivityBuilder(const GmshReader::GmshData&
{
ConnectivityDescriptor descriptor;
descriptor.setNodeNumberVector(gmsh_data.__verticesNumbers);
descriptor.setNodeNumberVector(gmsh_data.node_number);
Array<CellType> cell_type_vector(nb_cells);
Array<int> cell_number_vector(nb_cells);
Array<unsigned int> cell_to_node_row(nb_cells + 1);
parallel_for(cell_to_node_row.size(), PUGS_LAMBDA(const CellId cell_id) { cell_to_node_row[cell_id] = 2 * cell_id; });
parallel_for(
cell_to_node_row.size(), PUGS_LAMBDA(const CellId cell_id) { cell_to_node_row[cell_id] = 2 * cell_id; });
Array<unsigned int> cell_to_node_list(cell_to_node_row[cell_to_node_row.size() - 1]);
for (CellId cell_id = 0; cell_id < nb_cells; ++cell_id) {
for (int i_node = 0; i_node < 2; ++i_node) {
cell_to_node_list[2 * cell_id + i_node] = gmsh_data.__edges[cell_id][i_node];
cell_to_node_list[2 * cell_id + i_node] = gmsh_data.edge_entity[cell_id][i_node];
}
}
descriptor.setCellToNodeMatrix(ConnectivityMatrix(cell_to_node_row, cell_to_node_list));
......@@ -55,16 +55,14 @@ GmshConnectivityBuilder<1>::GmshConnectivityBuilder(const GmshReader::GmshData&
descriptor.setCellTypeVector(cell_type_vector);
for (size_t j = 0; j < nb_cells; ++j) {
cell_number_vector[j] = gmsh_data.__edges_number[j];
cell_number_vector[j] = gmsh_data.edge_entity_number[j];
}
descriptor.setCellNumberVector(cell_number_vector);
std::map<unsigned int, std::vector<unsigned int>> ref_points_map;
for (unsigned int r = 0; r < gmsh_data.__points.size(); ++r) {
const unsigned int point_number = gmsh_data.__points[r];
for (auto ref : gmsh_data.__points_ref[r]) {
ref_points_map[ref].push_back(point_number);
}
for (unsigned int r = 0; r < gmsh_data.point_entity.size(); ++r) {
const unsigned int point_number = gmsh_data.point_entity[r];
ref_points_map[gmsh_data.point_entity_ref[r]].push_back(point_number);
}
Array<size_t> node_nb_cell(descriptor.nodeNumberVector().size());
......@@ -111,8 +109,8 @@ GmshConnectivityBuilder<1>::GmshConnectivityBuilder(const GmshReader::GmshData&
}();
auto ref_id = [&] {
if (auto i_physical_ref = gmsh_data.m_physical_ref_map.find(ref_point_list.first);
i_physical_ref != gmsh_data.m_physical_ref_map.end()) {
auto i_physical_ref = gmsh_data.m_physical_ref_map.find(ref_point_list.first);
if ((i_physical_ref != gmsh_data.m_physical_ref_map.end()) and (i_physical_ref->second.dimension() == 0)) {
return i_physical_ref->second.refId();
} else {
return RefId{ref_point_list.first};
......@@ -123,11 +121,9 @@ GmshConnectivityBuilder<1>::GmshConnectivityBuilder(const GmshReader::GmshData&
}
std::map<unsigned int, std::vector<unsigned int>> ref_cells_map;
for (unsigned int j = 0; j < gmsh_data.__edges_ref.size(); ++j) {
for (unsigned int j = 0; j < gmsh_data.edge_entity_ref.size(); ++j) {
const unsigned int elem_number = j;
for (auto ref : gmsh_data.__edges_ref[j]) {
ref_cells_map[ref].push_back(elem_number);
}
ref_cells_map[gmsh_data.edge_entity_ref[j]].push_back(elem_number);
}
for (const auto& ref_cell_list : ref_cells_map) {
......@@ -137,8 +133,8 @@ GmshConnectivityBuilder<1>::GmshConnectivityBuilder(const GmshReader::GmshData&
}
auto ref_id = [&] {
if (auto i_physical_ref = gmsh_data.m_physical_ref_map.find(ref_cell_list.first);
i_physical_ref != gmsh_data.m_physical_ref_map.end()) {
auto i_physical_ref = gmsh_data.m_physical_ref_map.find(ref_cell_list.first);
if ((i_physical_ref != gmsh_data.m_physical_ref_map.end()) and (i_physical_ref->second.dimension() == 1)) {
return i_physical_ref->second.refId();
} else {
return RefId{ref_cell_list.first};
......@@ -168,18 +164,18 @@ GmshConnectivityBuilder<2>::GmshConnectivityBuilder(const GmshReader::GmshData&
{
ConnectivityDescriptor descriptor;
descriptor.setNodeNumberVector(gmsh_data.__verticesNumbers);
descriptor.setNodeNumberVector(gmsh_data.node_number);
Array<CellType> cell_type_vector(nb_cells);
Array<int> cell_number_vector(nb_cells);
Array<unsigned int> cell_to_node_row(gmsh_data.__triangles.size() + gmsh_data.__quadrangles.size() + 1);
Array<unsigned int> cell_to_node_row(gmsh_data.triangle_entity.size() + gmsh_data.quadrangle_entity.size() + 1);
{
cell_to_node_row[0] = 0;
size_t i_cell = 0;
for (size_t i_triangle = 0; i_triangle < gmsh_data.__triangles.size(); ++i_triangle, ++i_cell) {
for (size_t i_triangle = 0; i_triangle < gmsh_data.triangle_entity.size(); ++i_triangle, ++i_cell) {
cell_to_node_row[i_cell + 1] = cell_to_node_row[i_cell] + 3;
}
for (size_t i_quadrangle = 0; i_quadrangle < gmsh_data.__quadrangles.size(); ++i_quadrangle, ++i_cell) {
for (size_t i_quadrangle = 0; i_quadrangle < gmsh_data.quadrangle_entity.size(); ++i_quadrangle, ++i_cell) {
cell_to_node_row[i_cell + 1] = cell_to_node_row[i_cell] + 4;
}
}
......@@ -187,50 +183,46 @@ GmshConnectivityBuilder<2>::GmshConnectivityBuilder(const GmshReader::GmshData&
Array<unsigned int> cell_to_node_list(cell_to_node_row[cell_to_node_row.size() - 1]);
{
size_t i_cell_node = 0;
for (size_t i_triangle = 0; i_triangle < gmsh_data.__triangles.size(); ++i_triangle) {
cell_to_node_list[i_cell_node++] = gmsh_data.__triangles[i_triangle][0];
cell_to_node_list[i_cell_node++] = gmsh_data.__triangles[i_triangle][1];
cell_to_node_list[i_cell_node++] = gmsh_data.__triangles[i_triangle][2];
for (size_t i_triangle = 0; i_triangle < gmsh_data.triangle_entity.size(); ++i_triangle) {
cell_to_node_list[i_cell_node++] = gmsh_data.triangle_entity[i_triangle][0];
cell_to_node_list[i_cell_node++] = gmsh_data.triangle_entity[i_triangle][1];
cell_to_node_list[i_cell_node++] = gmsh_data.triangle_entity[i_triangle][2];
}
for (size_t i_quadrangle = 0; i_quadrangle < gmsh_data.__quadrangles.size(); ++i_quadrangle) {
cell_to_node_list[i_cell_node++] = gmsh_data.__quadrangles[i_quadrangle][0];
cell_to_node_list[i_cell_node++] = gmsh_data.__quadrangles[i_quadrangle][1];
cell_to_node_list[i_cell_node++] = gmsh_data.__quadrangles[i_quadrangle][2];
cell_to_node_list[i_cell_node++] = gmsh_data.__quadrangles[i_quadrangle][3];
for (size_t i_quadrangle = 0; i_quadrangle < gmsh_data.quadrangle_entity.size(); ++i_quadrangle) {
cell_to_node_list[i_cell_node++] = gmsh_data.quadrangle_entity[i_quadrangle][0];
cell_to_node_list[i_cell_node++] = gmsh_data.quadrangle_entity[i_quadrangle][1];
cell_to_node_list[i_cell_node++] = gmsh_data.quadrangle_entity[i_quadrangle][2];
cell_to_node_list[i_cell_node++] = gmsh_data.quadrangle_entity[i_quadrangle][3];
}
}
descriptor.setCellToNodeMatrix(ConnectivityMatrix(cell_to_node_row, cell_to_node_list));
const size_t nb_triangles = gmsh_data.__triangles.size();
const size_t nb_triangles = gmsh_data.triangle_entity.size();
for (size_t i_triangle = 0; i_triangle < nb_triangles; ++i_triangle) {
cell_type_vector[i_triangle] = CellType::Triangle;
cell_number_vector[i_triangle] = gmsh_data.__triangles_number[i_triangle];
cell_number_vector[i_triangle] = gmsh_data.triangle_entity_number[i_triangle];
}
const size_t nb_quadrangles = gmsh_data.__quadrangles.size();
const size_t nb_quadrangles = gmsh_data.quadrangle_entity.size();
for (size_t i_quadrangle = 0; i_quadrangle < nb_quadrangles; ++i_quadrangle) {
const size_t i_cell = i_quadrangle + nb_triangles;
cell_type_vector[i_cell] = CellType::Quadrangle;
cell_number_vector[i_cell] = gmsh_data.__quadrangles_number[i_quadrangle];
cell_number_vector[i_cell] = gmsh_data.quadrangle_entity_number[i_quadrangle];
}
descriptor.setCellNumberVector(cell_number_vector);
descriptor.setCellTypeVector(cell_type_vector);
std::map<unsigned int, std::vector<unsigned int>> ref_cells_map;
for (unsigned int j = 0; j < gmsh_data.__triangles_ref.size(); ++j) {
for (unsigned int j = 0; j < gmsh_data.triangle_entity_ref.size(); ++j) {
const unsigned int elem_number = j;
for (auto ref : gmsh_data.__triangles_ref[j]) {
ref_cells_map[ref].push_back(elem_number);
}
ref_cells_map[gmsh_data.triangle_entity_ref[j]].push_back(elem_number);
}
for (unsigned int j = 0; j < gmsh_data.__quadrangles_ref.size(); ++j) {
for (unsigned int j = 0; j < gmsh_data.quadrangle_entity_ref.size(); ++j) {
const size_t elem_number = nb_triangles + j;
for (auto ref : gmsh_data.__quadrangles_ref[j]) {
ref_cells_map[ref].push_back(elem_number);
}
ref_cells_map[gmsh_data.quadrangle_entity_ref[j]].push_back(elem_number);
}
for (const auto& ref_cell_list : ref_cells_map) {
......@@ -240,8 +232,8 @@ GmshConnectivityBuilder<2>::GmshConnectivityBuilder(const GmshReader::GmshData&
}
auto ref_id = [&] {
if (auto i_physical_ref = gmsh_data.m_physical_ref_map.find(ref_cell_list.first);
i_physical_ref != gmsh_data.m_physical_ref_map.end()) {
auto i_physical_ref = gmsh_data.m_physical_ref_map.find(ref_cell_list.first);
if ((i_physical_ref != gmsh_data.m_physical_ref_map.end()) and (i_physical_ref->second.dimension() == 2)) {
return i_physical_ref->second.refId();
} else {
return RefId{ref_cell_list.first};
......@@ -263,12 +255,12 @@ GmshConnectivityBuilder<2>::GmshConnectivityBuilder(const GmshReader::GmshData&
return mutable_face_number_id_map;
}();
const auto& node_number_vector = descriptor.nodeNumberVector();
const auto& node_number = descriptor.nodeNumberVector();
const auto& node_to_face_matrix = descriptor.nodeToFaceMatrix();
const auto& face_to_node_matrix = descriptor.faceToNodeMatrix();
const auto find_face = [&](uint32_t node0, uint32_t node1) {
if (node_number_vector[node0] > node_number_vector[node1]) {
if (node_number[node0] > node_number[node1]) {
std::swap(node0, node1);
}
const auto& face_id_vector = node_to_face_matrix[node0];
......@@ -287,15 +279,13 @@ GmshConnectivityBuilder<2>::GmshConnectivityBuilder(const GmshReader::GmshData&
Array<int> face_number_vector = copy(descriptor.faceNumberVector());
std::map<unsigned int, std::vector<unsigned int>> ref_faces_map;
for (unsigned int e = 0; e < gmsh_data.__edges.size(); ++e) {
const unsigned int edge_id = find_face(gmsh_data.__edges[e][0], gmsh_data.__edges[e][1]);
for (unsigned int e = 0; e < gmsh_data.edge_entity.size(); ++e) {
const unsigned int edge_id = find_face(gmsh_data.edge_entity[e][0], gmsh_data.edge_entity[e][1]);
for (auto ref : gmsh_data.__edges_ref[e]) {
ref_faces_map[ref].push_back(edge_id);
}
ref_faces_map[gmsh_data.edge_entity_ref[e]].push_back(edge_id);
if (face_number_vector[edge_id] != gmsh_data.__edges_number[e]) {
if (auto i_face = face_number_id_map.find(gmsh_data.__edges_number[e]); i_face != face_number_id_map.end()) {
if (face_number_vector[edge_id] != gmsh_data.edge_entity_number[e]) {
if (auto i_face = face_number_id_map.find(gmsh_data.edge_entity_number[e]); i_face != face_number_id_map.end()) {
const int other_edge_id = i_face->second;
std::swap(face_number_vector[edge_id], face_number_vector[other_edge_id]);
......@@ -306,7 +296,7 @@ GmshConnectivityBuilder<2>::GmshConnectivityBuilder(const GmshReader::GmshData&
face_number_id_map[face_number_vector[other_edge_id]] = other_edge_id;
} else {
face_number_id_map.erase(face_number_vector[edge_id]);
face_number_vector[edge_id] = gmsh_data.__edges_number[e];
face_number_vector[edge_id] = gmsh_data.edge_entity_number[e];
face_number_id_map[face_number_vector[edge_id]] = edge_id;
}
}
......@@ -331,8 +321,8 @@ GmshConnectivityBuilder<2>::GmshConnectivityBuilder(const GmshReader::GmshData&
}
auto ref_id = [&] {
if (auto i_physical_ref = gmsh_data.m_physical_ref_map.find(ref_face_list.first);
i_physical_ref != gmsh_data.m_physical_ref_map.end()) {
auto i_physical_ref = gmsh_data.m_physical_ref_map.find(ref_face_list.first);
if ((i_physical_ref != gmsh_data.m_physical_ref_map.end()) and (i_physical_ref->second.dimension() == 1)) {
return i_physical_ref->second.refId();
} else {
return RefId{ref_face_list.first};
......@@ -368,7 +358,7 @@ GmshConnectivityBuilder<2>::GmshConnectivityBuilder(const GmshReader::GmshData&
descriptor.addRefItemList(RefFaceList{ref_id, face_list, list_type});
}
Array<bool> is_boundary_node(node_number_vector.size());
Array<bool> is_boundary_node(node_number.size());
is_boundary_node.fill(false);
for (size_t i_face = 0; i_face < face_nb_cell.size(); ++i_face) {
......@@ -381,11 +371,9 @@ GmshConnectivityBuilder<2>::GmshConnectivityBuilder(const GmshReader::GmshData&
}
std::map<unsigned int, std::vector<unsigned int>> ref_points_map;
for (unsigned int r = 0; r < gmsh_data.__points.size(); ++r) {
const unsigned int point_number = gmsh_data.__points[r];
for (auto ref : gmsh_data.__points_ref[r]) {
ref_points_map[ref].push_back(point_number);
}
for (unsigned int r = 0; r < gmsh_data.point_entity.size(); ++r) {
const unsigned int point_number = gmsh_data.point_entity[r];
ref_points_map[gmsh_data.point_entity_ref[r]].push_back(point_number);
}
for (const auto& ref_point_list : ref_points_map) {
......@@ -395,8 +383,8 @@ GmshConnectivityBuilder<2>::GmshConnectivityBuilder(const GmshReader::GmshData&
}
auto ref_id = [&] {
if (auto i_physical_ref = gmsh_data.m_physical_ref_map.find(ref_point_list.first);
i_physical_ref != gmsh_data.m_physical_ref_map.end()) {
auto i_physical_ref = gmsh_data.m_physical_ref_map.find(ref_point_list.first);
if ((i_physical_ref != gmsh_data.m_physical_ref_map.end()) and (i_physical_ref->second.dimension() == 0)) {
return i_physical_ref->second.refId();
} else {
return RefId{ref_point_list.first};
......@@ -434,7 +422,7 @@ GmshConnectivityBuilder<2>::GmshConnectivityBuilder(const GmshReader::GmshData&
}());
descriptor.setNodeOwnerVector([&] {
Array<int> node_owner_vector(node_number_vector.size());
Array<int> node_owner_vector(node_number.size());
node_owner_vector.fill(parallel::rank());
return node_owner_vector;
}());
......@@ -447,14 +435,14 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData&
{
ConnectivityDescriptor descriptor;
descriptor.setNodeNumberVector(gmsh_data.__verticesNumbers);
descriptor.setNodeNumberVector(gmsh_data.node_number);
Array<CellType> cell_type_vector(nb_cells);
Array<int> cell_number_vector(nb_cells);
const size_t nb_tetrahedra = gmsh_data.__tetrahedra.size();
const size_t nb_hexahedra = gmsh_data.__hexahedra.size();
const size_t nb_prisms = gmsh_data.__prisms.size();
const size_t nb_pyramids = gmsh_data.__pyramids.size();
const size_t nb_tetrahedra = gmsh_data.tetrahedron_entity.size();
const size_t nb_hexahedra = gmsh_data.hexahedron_entity.size();
const size_t nb_prisms = gmsh_data.prism_entity.size();
const size_t nb_pyramids = gmsh_data.pyramid_entity.size();
Array<unsigned int> cell_to_node_row(nb_cells + 1);
{
......@@ -478,62 +466,62 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData&
{
size_t i_cell_node = 0;
for (size_t i_tetrahedron = 0; i_tetrahedron < nb_tetrahedra; ++i_tetrahedron) {
cell_to_node_list[i_cell_node++] = gmsh_data.__tetrahedra[i_tetrahedron][0];
cell_to_node_list[i_cell_node++] = gmsh_data.__tetrahedra[i_tetrahedron][1];
cell_to_node_list[i_cell_node++] = gmsh_data.__tetrahedra[i_tetrahedron][2];
cell_to_node_list[i_cell_node++] = gmsh_data.__tetrahedra[i_tetrahedron][3];
cell_to_node_list[i_cell_node++] = gmsh_data.tetrahedron_entity[i_tetrahedron][0];
cell_to_node_list[i_cell_node++] = gmsh_data.tetrahedron_entity[i_tetrahedron][1];
cell_to_node_list[i_cell_node++] = gmsh_data.tetrahedron_entity[i_tetrahedron][2];
cell_to_node_list[i_cell_node++] = gmsh_data.tetrahedron_entity[i_tetrahedron][3];
}
for (size_t i_hexahedron = 0; i_hexahedron < nb_hexahedra; ++i_hexahedron) {
cell_to_node_list[i_cell_node++] = gmsh_data.__hexahedra[i_hexahedron][0];
cell_to_node_list[i_cell_node++] = gmsh_data.__hexahedra[i_hexahedron][1];
cell_to_node_list[i_cell_node++] = gmsh_data.__hexahedra[i_hexahedron][2];
cell_to_node_list[i_cell_node++] = gmsh_data.__hexahedra[i_hexahedron][3];
cell_to_node_list[i_cell_node++] = gmsh_data.__hexahedra[i_hexahedron][4];
cell_to_node_list[i_cell_node++] = gmsh_data.__hexahedra[i_hexahedron][5];
cell_to_node_list[i_cell_node++] = gmsh_data.__hexahedra[i_hexahedron][6];
cell_to_node_list[i_cell_node++] = gmsh_data.__hexahedra[i_hexahedron][7];
cell_to_node_list[i_cell_node++] = gmsh_data.hexahedron_entity[i_hexahedron][0];
cell_to_node_list[i_cell_node++] = gmsh_data.hexahedron_entity[i_hexahedron][1];
cell_to_node_list[i_cell_node++] = gmsh_data.hexahedron_entity[i_hexahedron][2];
cell_to_node_list[i_cell_node++] = gmsh_data.hexahedron_entity[i_hexahedron][3];
cell_to_node_list[i_cell_node++] = gmsh_data.hexahedron_entity[i_hexahedron][4];
cell_to_node_list[i_cell_node++] = gmsh_data.hexahedron_entity[i_hexahedron][5];
cell_to_node_list[i_cell_node++] = gmsh_data.hexahedron_entity[i_hexahedron][6];
cell_to_node_list[i_cell_node++] = gmsh_data.hexahedron_entity[i_hexahedron][7];
}
for (size_t i_prism = 0; i_prism < nb_prisms; ++i_prism) {
cell_to_node_list[i_cell_node++] = gmsh_data.__prisms[i_prism][0];
cell_to_node_list[i_cell_node++] = gmsh_data.__prisms[i_prism][1];
cell_to_node_list[i_cell_node++] = gmsh_data.__prisms[i_prism][2];
cell_to_node_list[i_cell_node++] = gmsh_data.__prisms[i_prism][3];
cell_to_node_list[i_cell_node++] = gmsh_data.__prisms[i_prism][4];
cell_to_node_list[i_cell_node++] = gmsh_data.__prisms[i_prism][5];
cell_to_node_list[i_cell_node++] = gmsh_data.prism_entity[i_prism][0];
cell_to_node_list[i_cell_node++] = gmsh_data.prism_entity[i_prism][1];
cell_to_node_list[i_cell_node++] = gmsh_data.prism_entity[i_prism][2];
cell_to_node_list[i_cell_node++] = gmsh_data.prism_entity[i_prism][3];
cell_to_node_list[i_cell_node++] = gmsh_data.prism_entity[i_prism][4];
cell_to_node_list[i_cell_node++] = gmsh_data.prism_entity[i_prism][5];
}
for (size_t i_pyramid = 0; i_pyramid < nb_pyramids; ++i_pyramid) {
cell_to_node_list[i_cell_node++] = gmsh_data.__pyramids[i_pyramid][0];
cell_to_node_list[i_cell_node++] = gmsh_data.__pyramids[i_pyramid][1];
cell_to_node_list[i_cell_node++] = gmsh_data.__pyramids[i_pyramid][2];
cell_to_node_list[i_cell_node++] = gmsh_data.__pyramids[i_pyramid][3];
cell_to_node_list[i_cell_node++] = gmsh_data.__pyramids[i_pyramid][4];
cell_to_node_list[i_cell_node++] = gmsh_data.pyramid_entity[i_pyramid][0];
cell_to_node_list[i_cell_node++] = gmsh_data.pyramid_entity[i_pyramid][1];
cell_to_node_list[i_cell_node++] = gmsh_data.pyramid_entity[i_pyramid][2];
cell_to_node_list[i_cell_node++] = gmsh_data.pyramid_entity[i_pyramid][3];
cell_to_node_list[i_cell_node++] = gmsh_data.pyramid_entity[i_pyramid][4];
}
}
for (size_t j = 0; j < nb_tetrahedra; ++j) {
cell_type_vector[j] = CellType::Tetrahedron;
cell_number_vector[j] = gmsh_data.__tetrahedra_number[j];
cell_number_vector[j] = gmsh_data.tetrahedron_entity_number[j];
}
for (size_t j = 0; j < nb_hexahedra; ++j) {
const size_t jh = nb_tetrahedra + j;
cell_type_vector[jh] = CellType::Hexahedron;
cell_number_vector[jh] = gmsh_data.__hexahedra_number[j];
cell_number_vector[jh] = gmsh_data.hexahedron_entity_number[j];
}
for (size_t j = 0; j < nb_prisms; ++j) {
const size_t jp = nb_tetrahedra + nb_hexahedra + j;
cell_type_vector[jp] = CellType::Prism;
cell_number_vector[jp] = gmsh_data.__prisms_number[j];
cell_number_vector[jp] = gmsh_data.prism_entity_number[j];
}
for (size_t j = 0; j < nb_pyramids; ++j) {
const size_t jh = nb_tetrahedra + nb_hexahedra + nb_prisms + j;
cell_type_vector[jh] = CellType::Pyramid;
cell_number_vector[jh] = gmsh_data.__pyramids_number[j];
cell_number_vector[jh] = gmsh_data.pyramid_entity_number[j];
}
descriptor.setCellNumberVector(cell_number_vector);
......@@ -541,32 +529,24 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData&
descriptor.setCellToNodeMatrix(ConnectivityMatrix(cell_to_node_row, cell_to_node_list));
std::map<unsigned int, std::vector<unsigned int>> ref_cells_map;
for (unsigned int j = 0; j < gmsh_data.__tetrahedra_ref.size(); ++j) {
for (unsigned int j = 0; j < gmsh_data.tetrahedron_entity_ref.size(); ++j) {
const unsigned int elem_number = j;
for (auto ref : gmsh_data.__tetrahedra_ref[j]) {
ref_cells_map[ref].push_back(elem_number);
}
ref_cells_map[gmsh_data.tetrahedron_entity_ref[j]].push_back(elem_number);
}
for (unsigned int j = 0; j < gmsh_data.__hexahedra_ref.size(); ++j) {
for (unsigned int j = 0; j < gmsh_data.hexahedron_entity_ref.size(); ++j) {
const size_t elem_number = nb_tetrahedra + j;
for (auto ref : gmsh_data.__hexahedra_ref[j]) {
ref_cells_map[ref].push_back(elem_number);
}
ref_cells_map[gmsh_data.hexahedron_entity_ref[j]].push_back(elem_number);
}
for (unsigned int j = 0; j < gmsh_data.__prisms_ref.size(); ++j) {
for (unsigned int j = 0; j < gmsh_data.prism_entity_ref.size(); ++j) {
const size_t elem_number = nb_tetrahedra + nb_hexahedra + j;
for (auto ref : gmsh_data.__prisms_ref[j]) {
ref_cells_map[ref].push_back(elem_number);
}
ref_cells_map[gmsh_data.prism_entity_ref[j]].push_back(elem_number);
}
for (unsigned int j = 0; j < gmsh_data.__pyramids_ref.size(); ++j) {
for (unsigned int j = 0; j < gmsh_data.pyramid_entity_ref.size(); ++j) {
const size_t elem_number = nb_tetrahedra + nb_hexahedra + nb_prisms + j;
for (auto ref : gmsh_data.__pyramids_ref[j]) {
ref_cells_map[ref].push_back(elem_number);
}
ref_cells_map[gmsh_data.pyramid_entity_ref[j]].push_back(elem_number);
}
for (const auto& ref_cell_list : ref_cells_map) {
......@@ -576,8 +556,8 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData&
}
auto ref_id = [&] {
if (auto i_physical_ref = gmsh_data.m_physical_ref_map.find(ref_cell_list.first);
i_physical_ref != gmsh_data.m_physical_ref_map.end()) {
auto i_physical_ref = gmsh_data.m_physical_ref_map.find(ref_cell_list.first);
if ((i_physical_ref != gmsh_data.m_physical_ref_map.end()) and (i_physical_ref->second.dimension() == 3)) {
return i_physical_ref->second.refId();
} else {
return RefId{ref_cell_list.first};
......@@ -589,7 +569,7 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData&
ConnectivityBuilderBase::_computeCellFaceAndFaceNodeConnectivities<3>(descriptor);
const auto& node_number_vector = descriptor.nodeNumberVector();
const auto& node_number = descriptor.nodeNumberVector();
Array<size_t> face_nb_cell(descriptor.faceNumberVector().size());
face_nb_cell.fill(0);
......@@ -617,10 +597,10 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData&
const auto& node_to_face_matrix = descriptor.nodeToFaceMatrix();
const auto& face_to_node_matrix = descriptor.faceToNodeMatrix();
const auto find_face = [&](std::vector<uint32_t> node_list) {
const auto find_face = [&](std::vector<unsigned int> node_list) {
size_t i_node_smallest_number = 0;
for (size_t i_node = 1; i_node < node_list.size(); ++i_node) {
if (node_number_vector[node_list[i_node]] < node_number_vector[node_list[i_node_smallest_number]]) {
if (node_number[node_list[i_node]] < node_number[node_list[i_node_smallest_number]]) {
i_node_smallest_number = i_node;
}
}
......@@ -639,7 +619,7 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData&
}
}
if (node_number_vector[node_list[1]] > node_number_vector[node_list[node_list.size() - 1]]) {
if (node_number[node_list[1]] > node_number[node_list[node_list.size() - 1]]) {
for (size_t i_node = 1; i_node <= (node_list.size() + 1) / 2 - 1; ++i_node) {
std::swap(node_list[i_node], node_list[node_list.size() - i_node]);
}
......@@ -672,16 +652,14 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData&
Array<int> face_number_vector = copy(descriptor.faceNumberVector());
std::map<unsigned int, std::vector<unsigned int>> ref_faces_map;
for (unsigned int f = 0; f < gmsh_data.__triangles.size(); ++f) {
for (unsigned int f = 0; f < gmsh_data.triangle_entity.size(); ++f) {
const unsigned int face_id =
find_face({gmsh_data.__triangles[f][0], gmsh_data.__triangles[f][1], gmsh_data.__triangles[f][2]});
find_face({gmsh_data.triangle_entity[f][0], gmsh_data.triangle_entity[f][1], gmsh_data.triangle_entity[f][2]});
for (auto ref : gmsh_data.__triangles_ref[f]) {
ref_faces_map[ref].push_back(face_id);
}
ref_faces_map[gmsh_data.triangle_entity_ref[f]].push_back(face_id);
if (face_number_vector[face_id] != gmsh_data.__triangles_number[f]) {
if (auto i_face = face_number_id_map.find(gmsh_data.__triangles_number[f]);
if (face_number_vector[face_id] != gmsh_data.triangle_entity_number[f]) {
if (auto i_face = face_number_id_map.find(gmsh_data.triangle_entity_number[f]);
i_face != face_number_id_map.end()) {
const int other_face_id = i_face->second;
std::swap(face_number_vector[face_id], face_number_vector[other_face_id]);
......@@ -693,21 +671,20 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData&
face_number_id_map[face_number_vector[other_face_id]] = other_face_id;
} else {
face_number_id_map.erase(face_number_vector[face_id]);
face_number_vector[face_id] = gmsh_data.__triangles_number[f];
face_number_vector[face_id] = gmsh_data.triangle_entity_number[f];
face_number_id_map[face_number_vector[face_id]] = face_id;
}
}
}
for (unsigned int f = 0; f < gmsh_data.__quadrangles.size(); ++f) {
const unsigned int face_id = find_face({gmsh_data.__quadrangles[f][0], gmsh_data.__quadrangles[f][1],
gmsh_data.__quadrangles[f][2], gmsh_data.__quadrangles[f][3]});
for (unsigned int f = 0; f < gmsh_data.quadrangle_entity.size(); ++f) {
const unsigned int face_id = find_face({gmsh_data.quadrangle_entity[f][0], gmsh_data.quadrangle_entity[f][1],
gmsh_data.quadrangle_entity[f][2], gmsh_data.quadrangle_entity[f][3]});
for (auto ref : gmsh_data.__quadrangles_ref[f]) {
ref_faces_map[ref].push_back(face_id);
}
if (face_number_vector[face_id] != gmsh_data.__quadrangles_number[f]) {
if (auto i_face = face_number_id_map.find(gmsh_data.__quadrangles_number[f]);
ref_faces_map[gmsh_data.quadrangle_entity_ref[f]].push_back(face_id);
if (face_number_vector[face_id] != gmsh_data.quadrangle_entity_number[f]) {
if (auto i_face = face_number_id_map.find(gmsh_data.quadrangle_entity_number[f]);
i_face != face_number_id_map.end()) {
const int other_face_id = i_face->second;
std::swap(face_number_vector[face_id], face_number_vector[other_face_id]);
......@@ -719,7 +696,7 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData&
face_number_id_map[face_number_vector[other_face_id]] = other_face_id;
} else {
face_number_id_map.erase(face_number_vector[face_id]);
face_number_vector[face_id] = gmsh_data.__quadrangles_number[f];
face_number_vector[face_id] = gmsh_data.quadrangle_entity_number[f];
face_number_id_map[face_number_vector[face_id]] = face_id;
}
}
......@@ -733,8 +710,8 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData&
}
auto ref_id = [&] {
if (auto i_physical_ref = gmsh_data.m_physical_ref_map.find(ref_face_list.first);
i_physical_ref != gmsh_data.m_physical_ref_map.end()) {
auto i_physical_ref = gmsh_data.m_physical_ref_map.find(ref_face_list.first);
if ((i_physical_ref != gmsh_data.m_physical_ref_map.end()) and (i_physical_ref->second.dimension() == 2)) {
return i_physical_ref->second.refId();
} else {
return RefId{ref_face_list.first};
......@@ -792,7 +769,7 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData&
const auto& edge_to_node_matrix = descriptor.edgeToNodeMatrix();
const auto find_edge = [&](uint32_t node0, uint32_t node1) {
if (node_number_vector[node0] > node_number_vector[node1]) {
if (node_number[node0] > node_number[node1]) {
std::swap(node0, node1);
}
const auto& edge_id_vector = node_to_edge_matrix[node0];
......@@ -821,15 +798,14 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData&
Array<int> edge_number_vector = copy(descriptor.edgeNumberVector());
std::map<unsigned int, std::vector<unsigned int>> ref_edges_map;
for (unsigned int e = 0; e < gmsh_data.__edges.size(); ++e) {
const unsigned int edge_id = find_edge(gmsh_data.__edges[e][0], gmsh_data.__edges[e][1]);
for (unsigned int e = 0; e < gmsh_data.edge_entity.size(); ++e) {
const unsigned int edge_id = find_edge(gmsh_data.edge_entity[e][0], gmsh_data.edge_entity[e][1]);
for (auto ref : gmsh_data.__edges_ref[e]) {
ref_edges_map[ref].push_back(edge_id);
}
ref_edges_map[gmsh_data.edge_entity_ref[e]].push_back(edge_id);
if (edge_number_vector[edge_id] != gmsh_data.__edges_number[e]) {
if (auto i_edge = edge_number_id_map.find(gmsh_data.__edges_number[e]); i_edge != edge_number_id_map.end()) {
if (edge_number_vector[edge_id] != gmsh_data.edge_entity_number[e]) {
if (auto i_edge = edge_number_id_map.find(gmsh_data.edge_entity_number[e]);
i_edge != edge_number_id_map.end()) {
const int other_edge_id = i_edge->second;
std::swap(edge_number_vector[edge_id], edge_number_vector[other_edge_id]);
......@@ -840,7 +816,7 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData&
edge_number_id_map[edge_number_vector[other_edge_id]] = other_edge_id;
} else {
edge_number_id_map.erase(edge_number_vector[edge_id]);
edge_number_vector[edge_id] = gmsh_data.__edges_number[e];
edge_number_vector[edge_id] = gmsh_data.edge_entity_number[e];
edge_number_id_map[edge_number_vector[edge_id]] = edge_id;
}
}
......@@ -854,8 +830,8 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData&
}
auto ref_id = [&] {
if (auto i_physical_ref = gmsh_data.m_physical_ref_map.find(ref_edge_list.first);
i_physical_ref != gmsh_data.m_physical_ref_map.end()) {
auto i_physical_ref = gmsh_data.m_physical_ref_map.find(ref_edge_list.first);
if ((i_physical_ref != gmsh_data.m_physical_ref_map.end()) and (i_physical_ref->second.dimension() == 1)) {
return i_physical_ref->second.refId();
} else {
return RefId{ref_edge_list.first};
......@@ -880,7 +856,7 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData&
}
}
Array<bool> is_boundary_node(node_number_vector.size());
Array<bool> is_boundary_node(node_number.size());
is_boundary_node.fill(false);
const auto& face_to_node_matrix = descriptor.faceToNodeMatrix();
......@@ -896,11 +872,9 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData&
}
std::map<unsigned int, std::vector<unsigned int>> ref_points_map;
for (unsigned int r = 0; r < gmsh_data.__points.size(); ++r) {
const unsigned int point_number = gmsh_data.__points[r];
for (auto ref : gmsh_data.__points_ref[r]) {
ref_points_map[ref].push_back(point_number);
}
for (unsigned int r = 0; r < gmsh_data.point_entity.size(); ++r) {
const unsigned int point_number = gmsh_data.point_entity[r];
ref_points_map[gmsh_data.point_entity_ref[r]].push_back(point_number);
}
for (const auto& ref_point_list : ref_points_map) {
......@@ -910,8 +884,8 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData&
}
auto ref_id = [&] {
if (auto i_physical_ref = gmsh_data.m_physical_ref_map.find(ref_point_list.first);
i_physical_ref != gmsh_data.m_physical_ref_map.end()) {
auto i_physical_ref = gmsh_data.m_physical_ref_map.find(ref_point_list.first);
if ((i_physical_ref != gmsh_data.m_physical_ref_map.end()) and (i_physical_ref->second.dimension() == 0)) {
return i_physical_ref->second.refId();
} else {
return RefId{ref_point_list.first};
......@@ -954,7 +928,7 @@ GmshConnectivityBuilder<3>::GmshConnectivityBuilder(const GmshReader::GmshData&
}());
descriptor.setNodeOwnerVector([&] {
Array<int> node_owner_vector(node_number_vector.size());
Array<int> node_owner_vector(node_number.size());
node_owner_vector.fill(parallel::rank());
return node_owner_vector;
}());
......@@ -1341,8 +1315,8 @@ struct actions<nodes_section>
throw NormalError(error_msg.str());
}
gmsh_data.__verticesNumbers = vertex_number_list;
gmsh_data.__vertices = vertex_list;
gmsh_data.node_number = vertex_number_list;
gmsh_data.node_coordinate = vertex_list;
}
};
......@@ -1385,9 +1359,9 @@ struct actions<elements_section>
const size_t number_of_elements = nb_elements;
std::cout << "- Number Of Elements: " << number_of_elements << '\n' << std::flush;
gmsh_data.__elementNumber.resize(number_of_elements);
gmsh_data.__elementType.resize(number_of_elements);
gmsh_data.__references.resize(number_of_elements);
gmsh_data.entity_number.resize(number_of_elements);
gmsh_data.entity_type.resize(number_of_elements);
gmsh_data.entity_reference.resize(number_of_elements);
gmsh_data.__elementVertices.resize(number_of_elements);
std::vector<std::string_view> tokens;
......@@ -1412,33 +1386,29 @@ struct actions<elements_section>
throw std::runtime_error("empty element line");
}
gmsh_data.__elementNumber[i_element] = parse_int(tokens[0]);
gmsh_data.entity_number[i_element] = parse_int(tokens[0]);
const int element_type = parse_int(tokens[1]) - 1;
if ((element_type < 0) or (element_type > 15)) {
const int entity_type = parse_int(tokens[1]) - 1;
if ((entity_type < 0) or (entity_type > 15)) {
throw std::runtime_error("unknown element type");
}
if (not supported_element[element_type]) {
throw std::runtime_error(std::string{element_name[element_type]} + " is not supported");
if (not supported_element[entity_type]) {
throw std::runtime_error(std::string{element_name[entity_type]} + " is not supported");
}
gmsh_data.__elementType[i_element] = element_type;
gmsh_data.entity_type[i_element] = entity_type;
int nb_tags = parse_int(tokens[2]);
if (nb_tags < 0) {
throw std::runtime_error(std::string{element_name[element_type]} + " number of tags is negative");
throw std::runtime_error(std::string{element_name[entity_type]} + " number of tags is negative");
}
const size_t number_of_tags = nb_tags;
gmsh_data.__references[i_element].reserve(number_of_tags);
for (size_t i_tag = 0; i_tag < number_of_tags; ++i_tag) {
gmsh_data.__references[i_element].push_back(parse_int(tokens[3 + i_tag]));
#warning store all referencess
break;
}
gmsh_data.entity_reference[i_element] = parse_int(tokens[3]);
// In Gmsh file format 2.2 the first tag IS the physical tag
// other are spurious construction element tags
const size_t first_node_index = 3 + number_of_tags;
const size_t number_of_nodes = primitive_number_of_nodes[element_type];
const size_t number_of_nodes = primitive_number_of_nodes[entity_type];
gmsh_data.__elementVertices[i_element].resize(number_of_nodes);
for (size_t i_node = 0; i_node < number_of_nodes; ++i_node) {
......@@ -1548,27 +1518,25 @@ GmshReader::GmshReader(const std::string& filename) : m_filename(filename)
void
GmshReader::__proceedData()
{
TinyVector<15, int> elementNumber = zero;
std::vector<std::set<size_t>> elementReferences(15);
TinyVector<15, int> element_number = zero;
std::vector<std::set<size_t>> element_reference(15);
for (size_t i = 0; i < m_mesh_data.__elementType.size(); ++i) {
const int elementType = m_mesh_data.__elementType[i];
elementNumber[elementType]++;
for (auto ref : m_mesh_data.__references[i]) {
elementReferences[elementType].insert(ref);
}
for (size_t i = 0; i < m_mesh_data.entity_type.size(); ++i) {
const int elementType = m_mesh_data.entity_type[i];
element_number[elementType]++;
element_reference[elementType].insert(m_mesh_data.entity_reference[i]);
}
for (size_t i = 0; i < elementNumber.dimension(); ++i) {
if (elementNumber[i] > 0) {
std::cout << " - Number of " << gmsh_parser::element_name[i] << ": " << elementNumber[i];
for (size_t i = 0; i < element_number.dimension(); ++i) {
if (element_number[i] > 0) {
std::cout << " - Number of " << gmsh_parser::element_name[i] << ": " << element_number[i];
if (not(gmsh_parser::supported_element[i])) {
std::cout << " [" << rang::fg::yellow << "IGNORED" << rang::style::reset << "]";
} else {
std::cout << " references: ";
for (std::set<size_t>::const_iterator iref = elementReferences[i].begin(); iref != elementReferences[i].end();
for (std::set<size_t>::const_iterator iref = element_reference[i].begin(); iref != element_reference[i].end();
++iref) {
if (iref != elementReferences[i].begin())
if (iref != element_reference[i].begin())
std::cout << ',';
std::cout << *iref;
}
......@@ -1576,51 +1544,51 @@ GmshReader::__proceedData()
switch (i) {
// Supported entities
case 0: { // edges
m_mesh_data.__edges = Array<GmshData::Edge>(elementNumber[i]);
m_mesh_data.__edges_ref.resize(elementNumber[i]);
m_mesh_data.__edges_number.resize(elementNumber[i]);
m_mesh_data.edge_entity = Array<GmshData::Edge>(element_number[i]);
m_mesh_data.edge_entity_ref = Array<int>(element_number[i]);
m_mesh_data.edge_entity_number = Array<int>(element_number[i]);
break;
}
case 1: { // triangles
m_mesh_data.__triangles = Array<GmshData::Triangle>(elementNumber[i]);
m_mesh_data.__triangles_ref.resize(elementNumber[i]);
m_mesh_data.__triangles_number.resize(elementNumber[i]);
m_mesh_data.triangle_entity = Array<GmshData::Triangle>(element_number[i]);
m_mesh_data.triangle_entity_ref = Array<int>(element_number[i]);
m_mesh_data.triangle_entity_number = Array<int>(element_number[i]);
break;
}
case 2: { // quadrangles
m_mesh_data.__quadrangles = Array<GmshData::Quadrangle>(elementNumber[i]);
m_mesh_data.__quadrangles_ref.resize(elementNumber[i]);
m_mesh_data.__quadrangles_number.resize(elementNumber[i]);
m_mesh_data.quadrangle_entity = Array<GmshData::Quadrangle>(element_number[i]);
m_mesh_data.quadrangle_entity_ref = Array<int>(element_number[i]);
m_mesh_data.quadrangle_entity_number = Array<int>(element_number[i]);
break;
}
case 3: { // tetrahedra
m_mesh_data.__tetrahedra = Array<GmshData::Tetrahedron>(elementNumber[i]);
m_mesh_data.__tetrahedra_ref.resize(elementNumber[i]);
m_mesh_data.__tetrahedra_number.resize(elementNumber[i]);
m_mesh_data.tetrahedron_entity = Array<GmshData::Tetrahedron>(element_number[i]);
m_mesh_data.tetrahedron_entity_ref = Array<int>(element_number[i]);
m_mesh_data.tetrahedron_entity_number = Array<int>(element_number[i]);
break;
}
case 4: { // hexahedra
m_mesh_data.__hexahedra = Array<GmshData::Hexahedron>(elementNumber[i]);
m_mesh_data.__hexahedra_ref.resize(elementNumber[i]);
m_mesh_data.__hexahedra_number.resize(elementNumber[i]);
m_mesh_data.hexahedron_entity = Array<GmshData::Hexahedron>(element_number[i]);
m_mesh_data.hexahedron_entity_ref = Array<int>(element_number[i]);
m_mesh_data.hexahedron_entity_number = Array<int>(element_number[i]);
break;
}
case 5: { // prism
m_mesh_data.__prisms = Array<GmshData::Prism>(elementNumber[i]);
m_mesh_data.__prisms_ref.resize(elementNumber[i]);
m_mesh_data.__prisms_number.resize(elementNumber[i]);
m_mesh_data.prism_entity = Array<GmshData::Prism>(element_number[i]);
m_mesh_data.prism_entity_ref = Array<int>(element_number[i]);
m_mesh_data.prism_entity_number = Array<int>(element_number[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]);
m_mesh_data.pyramid_entity = Array<GmshData::Pyramid>(element_number[i]);
m_mesh_data.pyramid_entity_ref = Array<int>(element_number[i]);
m_mesh_data.pyramid_entity_number = Array<int>(element_number[i]);
break;
}
case 14: { // point
m_mesh_data.__points = Array<GmshData::Point>(elementNumber[i]);
m_mesh_data.__points_ref.resize(elementNumber[i]);
m_mesh_data.__points_number.resize(elementNumber[i]);
m_mesh_data.point_entity = Array<GmshData::Point>(element_number[i]);
m_mesh_data.point_entity_ref = Array<int>(element_number[i]);
m_mesh_data.point_entity_number = Array<int>(element_number[i]);
break;
}
// Unsupported entities
......@@ -1640,156 +1608,135 @@ GmshReader::__proceedData()
}
}
const int minNumber = min(m_mesh_data.__verticesNumbers);
const int maxNumber = max(m_mesh_data.__verticesNumbers);
const int min_node_number = min(m_mesh_data.node_number);
const int max_node_number = max(m_mesh_data.node_number);
if (minNumber < 0) {
if (min_node_number < 0) {
throw NotImplementedError("reading file '" + m_filename + "': negative vertices number");
}
#warning consider the case of very sparse node lists (max-min>3*nb nodes). Uses unordered_map?
// A value of -1 means that the vertex is unknown
m_mesh_data.__verticesCorrepondance.resize(maxNumber + 1, -1);
m_mesh_data.gmsh_node_to_node_id_correspondance.resize(max_node_number + 1, -1);
for (size_t i = 0; i < m_mesh_data.__verticesNumbers.size(); ++i) {
m_mesh_data.__verticesCorrepondance[m_mesh_data.__verticesNumbers[i]] = i;
for (size_t i = 0; i < m_mesh_data.node_number.size(); ++i) {
m_mesh_data.gmsh_node_to_node_id_correspondance[m_mesh_data.node_number[i]] = i;
}
// reset element number to count them while filling structures
elementNumber = zero;
element_number = zero;
// Element structures filling
for (size_t i = 0; i < m_mesh_data.__elementType.size(); ++i) {
for (size_t i = 0; i < m_mesh_data.entity_type.size(); ++i) {
std::vector<int>& elementVertices = m_mesh_data.__elementVertices[i];
switch (m_mesh_data.__elementType[i]) {
switch (m_mesh_data.entity_type[i]) {
// Supported entities
case 0: { // edge
int& edgeNumber = elementNumber[0];
const int a = m_mesh_data.__verticesCorrepondance[elementVertices[0]];
const int b = m_mesh_data.__verticesCorrepondance[elementVertices[1]];
if ((a < 0) or (b < 0)) {
throw NormalError("reading file '" + m_filename + "': error reading element " + stringify(i) +
" [bad vertices definition]");
}
m_mesh_data.__edges[edgeNumber] = GmshData::Edge(a, b);
m_mesh_data.__edges_ref[edgeNumber] = std::move(m_mesh_data.__references[i]);
m_mesh_data.__edges_number[edgeNumber] = m_mesh_data.__elementNumber[i];
edgeNumber++; // one more edge
int& edge_number = element_number[0];
const auto a = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[0]];
const auto b = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[1]];
m_mesh_data.edge_entity[edge_number] = GmshData::Edge{a, b};
m_mesh_data.edge_entity_ref[edge_number] = m_mesh_data.entity_reference[i];
m_mesh_data.edge_entity_number[edge_number] = m_mesh_data.entity_number[i];
edge_number++; // one more edge
break;
}
case 1: { // triangles
int& triangleNumber = elementNumber[1];
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]];
if ((a < 0) or (b < 0) or (c < 0)) {
throw NormalError("reading file '" + m_filename + "': error reading element " + stringify(i) +
" [bad vertices definition]");
}
m_mesh_data.__triangles[triangleNumber] = GmshData::Triangle(a, b, c);
m_mesh_data.__triangles_ref[triangleNumber] = std::move(m_mesh_data.__references[i]);
m_mesh_data.__triangles_number[triangleNumber] = m_mesh_data.__elementNumber[i];
triangleNumber++; // one more triangle
int& triangle_number = element_number[1];
const auto a = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[0]];
const auto b = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[1]];
const auto c = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[2]];
m_mesh_data.triangle_entity[triangle_number] = GmshData::Triangle{a, b, c};
m_mesh_data.triangle_entity_ref[triangle_number] = m_mesh_data.entity_reference[i];
m_mesh_data.triangle_entity_number[triangle_number] = m_mesh_data.entity_number[i];
triangle_number++; // one more triangle
break;
}
case 2: { // quadrangle
int& quadrilateralNumber = elementNumber[2];
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]];
if ((a < 0) or (b < 0) or (c < 0) or (d < 0)) {
throw NormalError("reading file '" + m_filename + "': error reading element " + stringify(i) +
" [bad vertices definition]");
}
m_mesh_data.__quadrangles[quadrilateralNumber] = GmshData::Quadrangle(a, b, c, d);
m_mesh_data.__quadrangles_ref[quadrilateralNumber] = std::move(m_mesh_data.__references[i]);
m_mesh_data.__quadrangles_number[quadrilateralNumber] = m_mesh_data.__elementNumber[i];
quadrilateralNumber++; // one more quadrangle
int& quadrilateral_number = element_number[2];
const auto a = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[0]];
const auto b = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[1]];
const auto c = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[2]];
const auto d = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[3]];
m_mesh_data.quadrangle_entity[quadrilateral_number] = GmshData::Quadrangle{a, b, c, d};
m_mesh_data.quadrangle_entity_ref[quadrilateral_number] = m_mesh_data.entity_reference[i];
m_mesh_data.quadrangle_entity_number[quadrilateral_number] = m_mesh_data.entity_number[i];
quadrilateral_number++; // one more quadrangle
break;
}
case 3: { // tetrahedra
int& tetrahedronNumber = elementNumber[3];
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]];
if ((a < 0) or (b < 0) or (c < 0) or (d < 0)) {
throw NormalError("reading file '" + m_filename + "': error reading element " + stringify(i) +
" [bad vertices definition]");
}
m_mesh_data.__tetrahedra[tetrahedronNumber] = GmshData::Tetrahedron(a, b, c, d);
m_mesh_data.__tetrahedra_ref[tetrahedronNumber] = std::move(m_mesh_data.__references[i]);
m_mesh_data.__tetrahedra_number[tetrahedronNumber] = m_mesh_data.__elementNumber[i];
tetrahedronNumber++; // one more tetrahedron
int& tetrahedron_number = element_number[3];
const auto a = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[0]];
const auto b = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[1]];
const auto c = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[2]];
const auto d = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[3]];
m_mesh_data.tetrahedron_entity[tetrahedron_number] = GmshData::Tetrahedron{a, b, c, d};
m_mesh_data.tetrahedron_entity_ref[tetrahedron_number] = m_mesh_data.entity_reference[i];
m_mesh_data.tetrahedron_entity_number[tetrahedron_number] = m_mesh_data.entity_number[i];
tetrahedron_number++; // one more tetrahedron
break;
}
case 4: { // hexaredron
int& hexahedronNumber = elementNumber[4];
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]];
const int f = m_mesh_data.__verticesCorrepondance[elementVertices[5]];
const int g = m_mesh_data.__verticesCorrepondance[elementVertices[6]];
const int h = m_mesh_data.__verticesCorrepondance[elementVertices[7]];
if ((a < 0) or (b < 0) or (c < 0) or (d < 0) or (e < 0) or (f < 0) or (g < 0) or (h < 0)) {
throw NormalError("reading file '" + m_filename + "': error reading element " + stringify(i) +
" [bad vertices definition]");
}
m_mesh_data.__hexahedra[hexahedronNumber] = GmshData::Hexahedron(a, b, c, d, e, f, g, h);
m_mesh_data.__hexahedra_ref[hexahedronNumber] = std::move(m_mesh_data.__references[i]);
m_mesh_data.__hexahedra_number[hexahedronNumber] = m_mesh_data.__elementNumber[i];
hexahedronNumber++; // one more hexahedron
int& hexahedron_number = element_number[4];
const auto a = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[0]];
const auto b = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[1]];
const auto c = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[2]];
const auto d = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[3]];
const auto e = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[4]];
const auto f = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[5]];
const auto g = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[6]];
const auto h = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[7]];
m_mesh_data.hexahedron_entity[hexahedron_number] = GmshData::Hexahedron{a, b, c, d, e, f, g, h};
m_mesh_data.hexahedron_entity_ref[hexahedron_number] = m_mesh_data.entity_reference[i];
m_mesh_data.hexahedron_entity_number[hexahedron_number] = m_mesh_data.entity_number[i];
hexahedron_number++; // one more hexahedron
break;
}
case 5: { // prism
int& prism_number = elementNumber[5];
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]];
const int f = m_mesh_data.__verticesCorrepondance[elementVertices[5]];
if ((a < 0) or (b < 0) or (c < 0) or (d < 0) or (e < 0) or (f < 0)) {
throw NormalError("reading file '" + m_filename + "': error reading element " + stringify(i) +
" [bad vertices definition]");
}
m_mesh_data.__prisms[prism_number] = GmshData::Prism(a, b, c, d, e, f);
m_mesh_data.__prisms_ref[prism_number] = std::move(m_mesh_data.__references[i]);
m_mesh_data.__prisms_number[prism_number] = m_mesh_data.__elementNumber[i];
int& prism_number = element_number[5];
const auto a = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[0]];
const auto b = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[1]];
const auto c = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[2]];
const auto d = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[3]];
const auto e = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[4]];
const auto f = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[5]];
m_mesh_data.prism_entity[prism_number] = GmshData::Prism{a, b, c, d, e, f};
m_mesh_data.prism_entity_ref[prism_number] = m_mesh_data.entity_reference[i];
m_mesh_data.prism_entity_number[prism_number] = m_mesh_data.entity_number[i];
prism_number++; // one more prism
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 " + stringify(i) +
" [bad vertices definition]");
}
m_mesh_data.__pyramids[pyramid_number] = GmshData::Pyramid(a, b, c, d, e);
m_mesh_data.__pyramids_ref[pyramid_number] = std::move(m_mesh_data.__references[i]);
m_mesh_data.__pyramids_number[pyramid_number] = m_mesh_data.__elementNumber[i];
int& pyramid_number = element_number[6];
const auto a = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[0]];
const auto b = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[1]];
const auto c = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[2]];
const auto d = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[3]];
const auto e = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[4]];
m_mesh_data.pyramid_entity[pyramid_number] = GmshData::Pyramid{a, b, c, d, e};
m_mesh_data.pyramid_entity_ref[pyramid_number] = m_mesh_data.entity_reference[i];
m_mesh_data.pyramid_entity_number[pyramid_number] = m_mesh_data.entity_number[i];
pyramid_number++; // one more hexahedron
break;
}
// Unsupported entities
case 14: { // point
int& point_number = elementNumber[14];
const int a = m_mesh_data.__verticesCorrepondance[elementVertices[0]];
m_mesh_data.__points[point_number] = a;
m_mesh_data.__points_ref[point_number] = std::move(m_mesh_data.__references[i]);
m_mesh_data.__points_number[point_number] = m_mesh_data.__elementNumber[i];
int& point_number = element_number[14];
m_mesh_data.point_entity[point_number] = m_mesh_data.gmsh_node_to_node_id_correspondance[elementVertices[0]];
m_mesh_data.point_entity_ref[point_number] = m_mesh_data.entity_reference[i];
m_mesh_data.point_entity_number[point_number] = m_mesh_data.entity_number[i];
point_number++;
break;
}
......@@ -1802,7 +1749,7 @@ GmshReader::__proceedData()
case 13: { // second order pyramid
}
default: {
throw NormalError("reading file '" + m_filename + "': ff3d cannot read this kind of element");
throw NormalError("reading file '" + m_filename + "': cannot read this kind of element");
}
}
}
......@@ -1827,8 +1774,8 @@ GmshReader::__proceedData()
dimension3_mask[12] = 1;
dimension3_mask[13] = 1;
if (dot(dimension3_mask, elementNumber) > 0) {
const size_t nb_cells = dot(dimension3_mask, elementNumber);
if (dot(dimension3_mask, element_number) > 0) {
const size_t nb_cells = dot(dimension3_mask, element_number);
GmshConnectivityBuilder<3> connectivity_builder(m_mesh_data, nb_cells);
......@@ -1839,17 +1786,12 @@ GmshReader::__proceedData()
using MeshType = Mesh<3>;
using Rd = TinyVector<3, double>;
NodeValue<Rd> xr(connectivity);
for (NodeId i = 0; i < m_mesh_data.__vertices.size(); ++i) {
xr[i][0] = m_mesh_data.__vertices[i][0];
xr[i][1] = m_mesh_data.__vertices[i][1];
xr[i][2] = m_mesh_data.__vertices[i][2];
}
NodeValue<Rd> xr(connectivity, m_mesh_data.node_coordinate);
m_mesh = std::make_shared<MeshVariant>(std::make_shared<const MeshType>(p_connectivity, xr));
this->_checkMesh<3>();
} else if (dot(dimension2_mask, elementNumber) > 0) {
const size_t nb_cells = dot(dimension2_mask, elementNumber);
} else if (dot(dimension2_mask, element_number) > 0) {
const size_t nb_cells = dot(dimension2_mask, element_number);
GmshConnectivityBuilder<2> connectivity_builder(m_mesh_data, nb_cells);
......@@ -1861,15 +1803,15 @@ GmshReader::__proceedData()
using Rd = TinyVector<2, double>;
NodeValue<Rd> xr(connectivity);
for (NodeId i = 0; i < m_mesh_data.__vertices.size(); ++i) {
xr[i][0] = m_mesh_data.__vertices[i][0];
xr[i][1] = m_mesh_data.__vertices[i][1];
for (NodeId i = 0; i < m_mesh_data.node_coordinate.size(); ++i) {
xr[i][0] = m_mesh_data.node_coordinate[i][0];
xr[i][1] = m_mesh_data.node_coordinate[i][1];
}
m_mesh = std::make_shared<MeshVariant>(std::make_shared<const MeshType>(p_connectivity, xr));
this->_checkMesh<2>();
} else if (dot(dimension1_mask, elementNumber) > 0) {
const size_t nb_cells = dot(dimension1_mask, elementNumber);
} else if (dot(dimension1_mask, element_number) > 0) {
const size_t nb_cells = dot(dimension1_mask, element_number);
GmshConnectivityBuilder<1> connectivity_builder(m_mesh_data, nb_cells);
......@@ -1881,8 +1823,8 @@ GmshReader::__proceedData()
using Rd = TinyVector<1, double>;
NodeValue<Rd> xr(connectivity);
for (NodeId i = 0; i < m_mesh_data.__vertices.size(); ++i) {
xr[i][0] = m_mesh_data.__vertices[i][0];
for (NodeId i = 0; i < m_mesh_data.node_coordinate.size(); ++i) {
xr[i][0] = m_mesh_data.node_coordinate[i][0];
}
m_mesh = std::make_shared<MeshVariant>(std::make_shared<const MeshType>(p_connectivity, xr));
......
src/mesh/GmshReader.hpp
View file @ dbc6302d
......@@ -55,10 +55,7 @@ class GmshReader : public MeshBuilderBase
PhysicalRefId& operator=(const PhysicalRefId&) = default;
PhysicalRefId& operator=(PhysicalRefId&&) = default;
PhysicalRefId(int dimension, const RefId& ref_id) : m_dimension(dimension), m_ref_id(ref_id)
{
;
}
PhysicalRefId(int dimension, const RefId& ref_id) : m_dimension(dimension), m_ref_id(ref_id) {}
PhysicalRefId() = default;
PhysicalRefId(const PhysicalRefId&) = default;
PhysicalRefId(PhysicalRefId&&) = default;
......@@ -81,70 +78,70 @@ class GmshReader : public MeshBuilderBase
* Gmsh format provides a numbered, none ordrered and none dense
* vertices list, this stores the number of read vertices.
*/
Array<int> __verticesNumbers;
Array<int> node_number;
Array<R3> __vertices;
Array<R3> node_coordinate;
using Point = unsigned int;
Array<Point> __points;
std::vector<std::vector<int>> __points_ref;
std::vector<int> __points_number;
using Edge = TinyVector<2, unsigned int>;
Array<Edge> __edges;
std::vector<std::vector<int>> __edges_ref;
std::vector<int> __edges_number;
using Triangle = TinyVector<3, unsigned int>;
Array<Triangle> __triangles;
std::vector<std::vector<int>> __triangles_ref;
std::vector<int> __triangles_number;
using Quadrangle = TinyVector<4, unsigned int>;
Array<Quadrangle> __quadrangles;
std::vector<std::vector<int>> __quadrangles_ref;
std::vector<int> __quadrangles_number;
using Tetrahedron = TinyVector<4, unsigned int>;
Array<Tetrahedron> __tetrahedra;
std::vector<std::vector<int>> __tetrahedra_ref;
std::vector<int> __tetrahedra_number;
using Pyramid = TinyVector<5, unsigned int>;
Array<Pyramid> __pyramids;
std::vector<std::vector<int>> __pyramids_ref;
std::vector<int> __pyramids_number;
using Prism = TinyVector<6, unsigned int>;
Array<Prism> __prisms;
std::vector<std::vector<int>> __prisms_ref;
std::vector<int> __prisms_number;
using Hexahedron = TinyVector<8, unsigned int>;
Array<Hexahedron> __hexahedra;
std::vector<std::vector<int>> __hexahedra_ref;
std::vector<int> __hexahedra_number;
Array<Point> point_entity;
Array<int> point_entity_ref;
Array<int> point_entity_number;
using Edge = std::array<unsigned int, 2>;
Array<Edge> edge_entity;
Array<int> edge_entity_ref;
Array<int> edge_entity_number;
using Triangle = std::array<unsigned int, 3>;
Array<Triangle> triangle_entity;
Array<int> triangle_entity_ref;
Array<int> triangle_entity_number;
using Quadrangle = std::array<unsigned int, 4>;
Array<Quadrangle> quadrangle_entity;
Array<int> quadrangle_entity_ref;
Array<int> quadrangle_entity_number;
using Tetrahedron = std::array<unsigned int, 4>;
Array<Tetrahedron> tetrahedron_entity;
Array<int> tetrahedron_entity_ref;
Array<int> tetrahedron_entity_number;
using Pyramid = std::array<unsigned int, 5>;
Array<Pyramid> pyramid_entity;
Array<int> pyramid_entity_ref;
Array<int> pyramid_entity_number;
using Prism = std::array<unsigned int, 6>;
Array<Prism> prism_entity;
Array<int> prism_entity_ref;
Array<int> prism_entity_number;
using Hexahedron = std::array<unsigned int, 8>;
Array<Hexahedron> hexahedron_entity;
Array<int> hexahedron_entity_ref;
Array<int> hexahedron_entity_number;
/**
* Gmsh format provides a numbered, none ordrered and none dense
* vertices list, this provides vertices renumbering correspondence
*/
std::vector<int> __verticesCorrepondance;
std::vector<unsigned int> gmsh_node_to_node_id_correspondance;
/**
* elements types
*/
std::vector<int> __elementNumber;
std::vector<int> entity_number;
/**
* elements types
*/
std::vector<short> __elementType;
std::vector<short> entity_type;
/**
* References
*/
std::vector<std::vector<int>> __references;
std::vector<int> entity_reference;
/**
* References
......