From de194e5db07bd726a02776ddad7e1d1257d34d75 Mon Sep 17 00:00:00 2001
From: Stephane Del Pino <stephane.delpino44@gmail.com>
Date: Wed, 10 Oct 2018 14:49:58 +0200
Subject: [PATCH] Read elements number in Gmsh file
This number is just an identifier and should not be used for anything else
Also prepared global index item structure. This has to be a contiguous numbering
starting from 0. This will have to take into account ghost cells.
---
src/mesh/Connectivity.cpp | 36 +++++++++++++++++++++++++-------
src/mesh/Connectivity.hpp | 19 ++++++++++++-----
src/mesh/GmshReader.cpp | 44 +++++++++++++++++++++++++++++----------
src/mesh/GmshReader.hpp | 11 ++++++++++
4 files changed, 87 insertions(+), 23 deletions(-)
diff --git a/src/mesh/Connectivity.cpp b/src/mesh/Connectivity.cpp
index 58649ddc9..cec70576f 100644
--- a/src/mesh/Connectivity.cpp
+++ b/src/mesh/Connectivity.cpp
@@ -216,16 +216,16 @@ void Connectivity<2>::_computeCellFaceAndFaceNodeConnectivities()
template<size_t Dimension>
Connectivity<Dimension>::
Connectivity(const std::vector<std::vector<unsigned int>>& cell_by_node_vector,
- const std::vector<CellType>& cell_type_vector)
+ const std::vector<CellType>& cell_type_vector,
+ const std::vector<int>& cell_number_vector)
{
Assert(cell_by_node_vector.size() == cell_type_vector.size());
+ Assert(cell_number_vector.size() == cell_type_vector.size());
auto& cell_to_node_matrix
= m_item_to_item_matrix[itemTId(ItemType::cell)][itemTId(ItemType::node)];
cell_to_node_matrix = cell_by_node_vector;
- Assert(this->numberOfCells()>0);
-
{
CellValue<CellType> cell_type(*this);
parallel_for(this->numberOfCells(), PASTIS_LAMBDA(const CellId& j){
@@ -233,9 +233,28 @@ Connectivity(const std::vector<std::vector<unsigned int>>& cell_by_node_vector,
});
m_cell_type = cell_type;
}
+
+ {
+ CellValue<int> cell_number(*this);
+ cell_number = convert_to_array(cell_number_vector);
+ m_cell_number = cell_number;
+ }
+
+ {
+ CellValue<int> cell_global_index(*this);
+#warning index must start accounting number of global indices of other procs
+#warning must take care of ghost cells
+ int first_index = 0;
+ parallel_for(this->numberOfCells(), PASTIS_LAMBDA(const CellId& j) {
+ cell_global_index[j] = first_index+j;
+ });
+ m_cell_global_index = cell_global_index;
+ }
+
+
{
CellValue<double> inv_cell_nb_nodes(*this);
- parallel_for(this->numberOfCells(), PASTIS_LAMBDA(const CellId& j){
+ parallel_for(this->numberOfCells(), PASTIS_LAMBDA(const CellId& j) {
const auto& cell_nodes = cell_to_node_matrix.rowConst(j);
inv_cell_nb_nodes[j] = 1./cell_nodes.length;
});
@@ -250,12 +269,15 @@ Connectivity(const std::vector<std::vector<unsigned int>>& cell_by_node_vector,
template Connectivity1D::
Connectivity(const std::vector<std::vector<unsigned int>>& cell_by_node_vector,
- const std::vector<CellType>& cell_type_vector);
+ const std::vector<CellType>& cell_type_vector,
+ const std::vector<int>& cell_number_vector);
template Connectivity2D::
Connectivity(const std::vector<std::vector<unsigned int>>& cell_by_node_vector,
- const std::vector<CellType>& cell_type_vector);
+ const std::vector<CellType>& cell_type_vector,
+ const std::vector<int>& cell_number_vector);
template Connectivity3D::
Connectivity(const std::vector<std::vector<unsigned int>>& cell_by_node_vector,
- const std::vector<CellType>& cell_type_vector);
+ const std::vector<CellType>& cell_type_vector,
+ const std::vector<int>& cell_number_vector);
diff --git a/src/mesh/Connectivity.hpp b/src/mesh/Connectivity.hpp
index 1a445be0d..7c990dc5c 100644
--- a/src/mesh/Connectivity.hpp
+++ b/src/mesh/Connectivity.hpp
@@ -246,6 +246,8 @@ class Connectivity final
private:
ConnectivityMatrix m_item_to_item_matrix[Dimension+1][Dimension+1];
CellValue<const CellType> m_cell_type;
+ CellValue<const int> m_cell_global_index;
+ CellValue<const int> m_cell_number;
FaceValuePerCell<const bool> m_cell_face_is_reversed;
@@ -316,6 +318,12 @@ class Connectivity final
return m_cell_type;
};
+ PASTIS_INLINE
+ const CellValue<const int>& cellNumber() const
+ {
+ return m_cell_number;
+ };
+
PASTIS_INLINE
const bool& isConnectivityMatrixBuilt(const ItemType& item_type_0,
const ItemType& item_type_1) const
@@ -548,8 +556,8 @@ class Connectivity final
for (FaceId l=0; l<this->numberOfFaces(); ++l) {
const auto& face_to_cell = face_to_cell_matrix[l];
if (face_to_cell.size() == 2) {
- const int cell_0 = face_to_cell[0];
- const int cell_1 = face_to_cell[1];
+ const CellId cell_0 = face_to_cell[0];
+ const CellId cell_1 = face_to_cell[1];
cell_cells[cell_0].insert(cell_1);
cell_cells[cell_1].insert(cell_0);
@@ -565,8 +573,8 @@ class Connectivity final
{
size_t k=0;
for (size_t j=0; j<this->numberOfCells(); ++j) {
- for (auto cell_id : cell_cells[j]) {
- neighbors[k] = cell_id;
+ for (CellId cell_id : cell_cells[j]) {
+ neighbors[k] = m_cell_global_index[cell_id];
++k;
}
}
@@ -627,7 +635,8 @@ class Connectivity final
Connectivity(const Connectivity&) = delete;
Connectivity(const std::vector<std::vector<unsigned int>>& cell_by_node_vector,
- const std::vector<CellType>& cell_type_vector);
+ const std::vector<CellType>& cell_type_vector,
+ const std::vector<int>& cell_number_vector);
~Connectivity()
{
diff --git a/src/mesh/GmshReader.cpp b/src/mesh/GmshReader.cpp
index 6b071de30..84aae75e8 100644
--- a/src/mesh/GmshReader.cpp
+++ b/src/mesh/GmshReader.cpp
@@ -150,7 +150,7 @@ void GmshReader::_dispatch()
using MeshType = Mesh<ConnectivityType>;
if (not m_mesh) {
- std::shared_ptr<ConnectivityType> connectivity(new ConnectivityType({},{}));
+ std::shared_ptr<ConnectivityType> connectivity(new ConnectivityType({},{},{}));
NodeValue<Rd> xr;
m_mesh = std::shared_ptr<IMesh>(new MeshType(connectivity, xr));
}
@@ -192,6 +192,7 @@ void GmshReader::_dispatch()
return std::move(cell_to_send_by_proc);
} ();
+ const auto& cell_number = mesh.connectivity().cellNumber();
const std::vector<Array<const int>> cell_number_to_send_by_proc
= [&] () {
std::vector<Array<const int>> cell_number_to_send_by_proc(commSize());
@@ -199,8 +200,7 @@ void GmshReader::_dispatch()
const Array<const CellId>& cell_list = cell_to_send_by_proc[i];
Array<int> cell_number_list(cell_list.size());
parallel_for (cell_list.size(), PASTIS_LAMBDA(const CellId& cell_id) {
-#warning must use table of cell numbers
- cell_number_list[cell_id] = cell_list[cell_id];
+ cell_number_list[cell_id] = cell_number[cell_list[cell_id]];
});
cell_number_to_send_by_proc[i] = cell_number_list;
}
@@ -214,9 +214,8 @@ void GmshReader::_dispatch()
Array<int> nb_cell_to_recv_by_proc = allToAll(nb_cell_to_send_by_proc);
- const size_t cell_number = mesh.numberOfCells();
const size_t new_cell_number
- = cell_number
+ = mesh.numberOfCells()
+ Sum(nb_cell_to_recv_by_proc)
- Sum(nb_cell_to_send_by_proc);
@@ -248,7 +247,7 @@ void GmshReader::_dispatch()
std::cout << '\n' << std::flush;
}
}
- barrier();
+ // barrier();
}
Messenger::destroy();
@@ -573,13 +572,14 @@ GmshReader::__readElements2_2()
ErrorHandler::normal);
}
+ __elementNumber.resize(numberOfElements);
__elementType.resize(numberOfElements);
__references.resize(numberOfElements);
__elementVertices.resize(numberOfElements);
if (not __binary) {
for (int i=0; i<numberOfElements; ++i) {
- this->_getInteger(); // drops element number
+ __elementNumber[i] = this->_getInteger();
const int elementType = this->_getInteger()-1;
if ((elementType < 0) or (elementType > 14)) {
@@ -722,31 +722,36 @@ GmshReader::__proceedData()
case 0: { // edges
__edges = Array<Edge>(elementNumber[i]);
__edges_ref.resize(elementNumber[i]);
+ __edges_number.resize(elementNumber[i]);
break;
}
case 1: { // triangles
__triangles = Array<Triangle>(elementNumber[i]);
__triangles_ref.resize(elementNumber[i]);
+ __triangles_number.resize(elementNumber[i]);
break;
}
case 2: { // quadrangles
__quadrangles = Array<Quadrangle>(elementNumber[i]);
__quadrangles_ref.resize(elementNumber[i]);
+ __quadrangles_number.resize(elementNumber[i]);
break;
}
case 3: { // tetrahedra
__tetrahedra = Array<Tetrahedron>(elementNumber[i]);
__tetrahedra_ref.resize(elementNumber[i]);
+ __tetrahedra_number.resize(elementNumber[i]);
break;
}
case 4: {// hexahedra
__hexahedra = Array<Hexahedron>(elementNumber[i]);
__hexahedra_ref.resize(elementNumber[i]);
+ __hexahedra_number.resize(elementNumber[i]);
}
- // Ignored entities
case 14: {// point
__points = Array<Point>(elementNumber[i]);
__points_ref.resize(elementNumber[i]);
+ __points_number.resize(elementNumber[i]);
break;
}
// Unsupported entities
@@ -787,6 +792,7 @@ GmshReader::__proceedData()
ErrorHandler::normal);
}
+#warning should use an unordered_map
// A value of -1 means that the vertex is unknown
__verticesCorrepondance.resize(maxNumber+1,-1);
@@ -836,6 +842,7 @@ GmshReader::__proceedData()
__triangles[triangleNumber]
= Triangle(a, b, c);
__triangles_ref[triangleNumber] = __references[i];
+ __triangles_number[triangleNumber] = __elementNumber[i];
triangleNumber++; // one more triangle
break;
}
@@ -855,6 +862,7 @@ GmshReader::__proceedData()
}
__quadrangles[quadrilateralNumber] = Quadrangle(a,b,c,d);
__quadrangles_ref[quadrilateralNumber] = __references[i];
+ __quadrangles_number[quadrilateralNumber] = __elementNumber[i];
quadrilateralNumber++; // one more quadrangle
break;
}
@@ -874,6 +882,7 @@ GmshReader::__proceedData()
}
__tetrahedra[tetrahedronNumber] = Tetrahedron(a, b, c, d);
__tetrahedra_ref[tetrahedronNumber] = __references[i];
+ __tetrahedra_number[tetrahedronNumber] = __elementNumber[i];
tetrahedronNumber++; // one more tetrahedron
break;
}
@@ -899,6 +908,7 @@ GmshReader::__proceedData()
= Hexahedron(a, b, c, d,
e, f, g, h);
__hexahedra_ref[hexahedronNumber] = __references[i];
+ __hexahedra_number[hexahedronNumber] = __elementNumber[i];
hexahedronNumber++; // one more hexahedron
break;
}
@@ -908,6 +918,7 @@ GmshReader::__proceedData()
const int a = __verticesCorrepondance[elementVertices[0]];
__points[point_number] = a;
__points_ref[point_number] = __references[i];
+ __points_number[point_number] = __elementNumber[i];
point_number++;
break;
}
@@ -958,6 +969,7 @@ GmshReader::__proceedData()
const size_t nb_cells = (dimension3_mask, elementNumber);
std::vector<CellType> cell_type_vector(nb_cells);
+ std::vector<int> cell_number_vector(nb_cells);
std::vector<std::vector<unsigned int>> cell_by_node_vector(nb_cells);
const size_t nb_tetrahedra = __tetrahedra.size();
@@ -967,6 +979,7 @@ GmshReader::__proceedData()
cell_by_node_vector[j][r] = __tetrahedra[j][r];
}
cell_type_vector[j] = CellType::Tetrahedron;
+ cell_number_vector[j] = __tetrahedra_number[j];
}
const size_t nb_hexahedra = __hexahedra.size();
for (size_t j=0; j<nb_hexahedra; ++j) {
@@ -976,10 +989,12 @@ GmshReader::__proceedData()
cell_by_node_vector[jh][r] = __hexahedra[j][r];
}
cell_type_vector[jh] = CellType::Hexahedron;
+ cell_number_vector[jh] = __hexahedra_number[j];
}
std::shared_ptr<Connectivity3D> p_connectivity(new Connectivity3D(cell_by_node_vector,
- cell_type_vector));
+ cell_type_vector,
+ cell_number_vector));
Connectivity3D& connectivity = *p_connectivity;
std::map<unsigned int, std::vector<unsigned int>> ref_faces_map;
@@ -1020,6 +1035,7 @@ GmshReader::__proceedData()
const size_t nb_cells = (dimension2_mask, elementNumber);
std::vector<CellType> cell_type_vector(nb_cells);
+ std::vector<int> cell_number_vector(nb_cells);
std::vector<std::vector<unsigned int>> cell_by_node_vector(nb_cells);
const size_t nb_triangles = __triangles.size();
@@ -1029,6 +1045,7 @@ GmshReader::__proceedData()
cell_by_node_vector[j][r] = __triangles[j][r];
}
cell_type_vector[j] = CellType::Triangle;
+ cell_number_vector[j] = __triangles_number[j];
}
const size_t nb_quadrangles = __quadrangles.size();
@@ -1039,10 +1056,12 @@ GmshReader::__proceedData()
cell_by_node_vector[jq][r] = __quadrangles[j][r];
}
cell_type_vector[jq] = CellType::Quadrangle;
+ cell_number_vector[jq] = __quadrangles_number[j];
}
std::shared_ptr<Connectivity2D> p_connectivity(new Connectivity2D(cell_by_node_vector,
- cell_type_vector));
+ cell_type_vector,
+ cell_number_vector));
Connectivity2D& connectivity = *p_connectivity;
std::map<unsigned int, std::vector<unsigned int>> ref_faces_map;
@@ -1091,6 +1110,7 @@ GmshReader::__proceedData()
const size_t nb_cells = (dimension1_mask, elementNumber);
std::vector<CellType> cell_type_vector(nb_cells);
+ std::vector<int> cell_number_vector(nb_cells);
std::vector<std::vector<unsigned int>> cell_by_node_vector(nb_cells);
for (size_t j=0; j<nb_cells; ++j) {
@@ -1099,10 +1119,12 @@ GmshReader::__proceedData()
cell_by_node_vector[j][r] = __edges[j][r];
}
cell_type_vector[j] = CellType::Line;
+ cell_number_vector[j] = __edges_number[j];
}
std::shared_ptr<Connectivity1D> p_connectivity(new Connectivity1D(cell_by_node_vector,
- cell_type_vector));
+ cell_type_vector,
+ cell_number_vector));
Connectivity1D& connectivity = *p_connectivity;
std::map<unsigned int, std::vector<unsigned int>> ref_points_map;
diff --git a/src/mesh/GmshReader.hpp b/src/mesh/GmshReader.hpp
index 7364878a1..20db0d321 100644
--- a/src/mesh/GmshReader.hpp
+++ b/src/mesh/GmshReader.hpp
@@ -83,26 +83,32 @@ private:
using Point = unsigned int;
Array<Point> __points;
std::vector<int> __points_ref;
+ std::vector<int> __points_number;
using Edge = TinyVector<2,unsigned int>;
Array<Edge> __edges;
std::vector<int> __edges_ref;
+ std::vector<int> __edges_number;
using Triangle = TinyVector<3,unsigned int>;
Array<Triangle> __triangles;
std::vector<int> __triangles_ref;
+ std::vector<int> __triangles_number;
using Quadrangle = TinyVector<4,unsigned int>;
Array<Quadrangle> __quadrangles;
std::vector<int> __quadrangles_ref;
+ std::vector<int> __quadrangles_number;
using Tetrahedron = TinyVector<4,unsigned int>;
Array<Tetrahedron> __tetrahedra;
std::vector<int> __tetrahedra_ref;
+ std::vector<int> __tetrahedra_number;
using Hexahedron = TinyVector<8,unsigned int>;
Array<Hexahedron> __hexahedra;
std::vector<int> __hexahedra_ref;
+ std::vector<int> __hexahedra_number;
/**
* Gmsh format provides a numbered, none ordrered and none dense
@@ -110,6 +116,11 @@ private:
*/
std::vector<int> __verticesCorrepondance;
+ /**
+ * elements types
+ */
+ std::vector<int> __elementNumber;
+
/**
* elements types
*/
--
GitLab