From 362adc0288e3a889e6a5af18cc4de9b8b47b4073 Mon Sep 17 00:00:00 2001
From: Stephane Del Pino <stephane.delpino44@gmail.com>
Date: Fri, 18 May 2018 18:10:34 +0200
Subject: [PATCH] Set of improvements
- GmshReader:
- now knows about vertices
- '"' are filtered from physical names
- uses RefId to store entity lists
- RefId: simple definition of references. only ref<->name association
- EntityListManager: begining of entity sets handling
---
src/mesh/Connectivity2D.hpp | 206 +++++---
src/mesh/GmshReader.cpp | 966 +++++++++++++++++++-----------------
src/mesh/GmshReader.hpp | 158 ++++--
src/mesh/RefId.hpp | 69 +++
4 files changed, 811 insertions(+), 588 deletions(-)
create mode 100644 src/mesh/RefId.hpp
diff --git a/src/mesh/Connectivity2D.hpp b/src/mesh/Connectivity2D.hpp
index 47d19c024..1603b386c 100644
--- a/src/mesh/Connectivity2D.hpp
+++ b/src/mesh/Connectivity2D.hpp
@@ -9,15 +9,86 @@
#include <map>
#include <algorithm>
+#include <RefId.hpp>
+
+class EntityListManager
+{
+ public:
+ typedef Kokkos::View<const unsigned int*> CellIdArray;
+ typedef Kokkos::View<const unsigned int*> FaceIdArray;
+ typedef Kokkos::View<const unsigned int*> EdgeIdArray;
+ typedef Kokkos::View<const unsigned int*> NodeIdArray;
+
+ typedef std::map<RefId, CellIdArray> CellIdArraySet;
+ typedef std::map<RefId, FaceIdArray> FaceIdArraySet;
+ typedef std::map<RefId, EdgeIdArray> EdgeIdArraySet;
+ typedef std::map<RefId, NodeIdArray> NodeIdArraySet;
+
+ private:
+ CellIdArraySet m_cell_id_array_set;
+ FaceIdArraySet m_face_id_array_set;
+ EdgeIdArraySet m_edge_id_array_set;
+ NodeIdArraySet m_node_id_array_set;
+
+ public:
+
+ void addCellIdArray(const RefId& set_id, const CellIdArray& cell_id_array)
+ {
+ m_cell_id_array_set.insert(std::make_pair(set_id, cell_id_array));
+ }
+
+ void addFaceIdArray(const RefId& set_id, const FaceIdArray& face_id_array)
+ {
+ m_face_id_array_set.insert(std::make_pair(set_id, face_id_array));
+ }
+
+ void addEdgeIdArray(const RefId& set_id, const EdgeIdArray& edge_id_array)
+ {
+ m_edge_id_array_set.insert(std::make_pair(set_id, edge_id_array));
+ }
+
+ void addNodeIdArray(const RefId& set_id, const NodeIdArray& node_id_array)
+ {
+ m_node_id_array_set.insert(std::make_pair(set_id, node_id_array));
+ }
+
+ const CellIdArraySet& getCellIdArraySet() const
+ {
+ return m_cell_id_array_set;
+ }
+
+ const FaceIdArraySet& getFaceIdArraySet() const
+ {
+ return m_face_id_array_set;
+ }
+
+ const EdgeIdArraySet& getEdgeIdArraySet() const
+ {
+ return m_edge_id_array_set;
+ }
+
+ const NodeIdArraySet& getNodeIdArraySet() const
+ {
+ return m_node_id_array_set;
+ }
+
+ EntityListManager& operator=(const EntityListManager&) = delete;
+ EntityListManager() = default;
+ EntityListManager(const EntityListManager&) = delete;
+ ~EntityListManager() = default;
+};
+
class Connectivity2D
{
public:
static constexpr size_t dimension = 2;
typedef std::pair<unsigned int, Kokkos::View<const unsigned int*>> FacesBoundary;
- typedef std::pair<unsigned int, Kokkos::View<const unsigned int*>> NodesBoundary;
+ typedef std::pair<RefId, Kokkos::View<const unsigned int*>> NodesBoundary;
typedef std::vector<FacesBoundary> FacesBoundaryList;
typedef std::vector<NodesBoundary> NodesBoundaryList;
+
+ EntityListManager m_entity_list_manager;
private:
const size_t m_number_of_cells;
size_t m_number_of_faces;
@@ -33,7 +104,7 @@ private:
Kokkos::View<unsigned short*> m_node_nb_cells;
Kokkos::View<unsigned int**> m_node_cells;
Kokkos::View<unsigned short**> m_node_cell_local_node;
-
+
Kokkos::View<unsigned short*> m_face_nb_cells;
Kokkos::View<unsigned int**> m_face_cells;
Kokkos::View<unsigned short**> m_face_cell_local_face;
@@ -56,8 +127,8 @@ private:
bool operator<(const Face& f) const
{
return ((m_node0_id<f.m_node0_id) or
- ((m_node0_id == f.m_node0_id) and
- (m_node1_id<f.m_node1_id)));
+ ((m_node0_id == f.m_node0_id) and
+ (m_node1_id<f.m_node1_id)));
}
KOKKOS_INLINE_FUNCTION
@@ -74,9 +145,9 @@ private:
KOKKOS_INLINE_FUNCTION
Face(unsigned int node0_id,
- unsigned int node1_id)
- : m_node0_id(node0_id),
- m_node1_id(node1_id)
+ unsigned int node1_id)
+ : m_node0_id(node0_id),
+ m_node1_id(node1_id)
{
;
}
@@ -93,19 +164,19 @@ private:
for (unsigned int j=0; j<m_number_of_cells; ++j) {
const unsigned short cell_nb_nodes = m_cell_nb_nodes[j];
for (unsigned short r=0; r<cell_nb_nodes; ++r) {
- unsigned int node0_id = m_cell_nodes(j,r);
- unsigned int node1_id = m_cell_nodes(j,(r+1)%cell_nb_nodes);
- if (node1_id<node0_id) {
- std::swap(node0_id, node1_id);
- }
- face_cells_map[Face(node0_id, node1_id)].push_back(std::make_pair(j, r));
+ unsigned int node0_id = m_cell_nodes(j,r);
+ unsigned int node1_id = m_cell_nodes(j,(r+1)%cell_nb_nodes);
+ if (node1_id<node0_id) {
+ std::swap(node0_id, node1_id);
+ }
+ face_cells_map[Face(node0_id, node1_id)].push_back(std::make_pair(j, r));
}
}
m_number_of_faces = face_cells_map.size();
Kokkos::View<unsigned short*> face_nb_nodes("face_nb_nodes", m_number_of_faces);
Kokkos::parallel_for(m_number_of_faces, KOKKOS_LAMBDA(const unsigned int& l) {
- face_nb_nodes[l] = 2;
+ face_nb_nodes[l] = 2;
});
m_face_nb_nodes = face_nb_nodes;
@@ -113,10 +184,10 @@ private:
{
int l=0;
for (const auto& face_cells_vector : face_cells_map) {
- const Face& face = face_cells_vector.first;
- face_nodes(l,0) = face.m_node0_id;
- face_nodes(l,1) = face.m_node1_id;
- ++l;
+ const Face& face = face_cells_vector.first;
+ face_nodes(l,0) = face.m_node0_id;
+ face_nodes(l,1) = face.m_node1_id;
+ ++l;
}
}
m_face_nodes = face_nodes;
@@ -125,9 +196,9 @@ private:
{
int l=0;
for (const auto& face_cells_vector : face_cells_map) {
- const auto& cells_vector = face_cells_vector.second;
- face_nb_cells[l] = cells_vector.size();
- ++l;
+ const auto& cells_vector = face_cells_vector.second;
+ face_nb_cells[l] = cells_vector.size();
+ ++l;
}
}
m_face_nb_cells = face_nb_cells;
@@ -136,12 +207,12 @@ private:
{
int l=0;
for (const auto& face_cells_vector : face_cells_map) {
- const auto& cells_vector = face_cells_vector.second;
- for (unsigned short lj=0; lj<cells_vector.size(); ++lj) {
- unsigned int cell_number = cells_vector[lj].first;
- face_cells(l,lj) = cell_number;
- }
- ++l;
+ const auto& cells_vector = face_cells_vector.second;
+ for (unsigned short lj=0; lj<cells_vector.size(); ++lj) {
+ unsigned int cell_number = cells_vector[lj].first;
+ face_cells(l,lj) = cell_number;
+ }
+ ++l;
}
}
m_face_cells = face_cells;
@@ -152,35 +223,34 @@ private:
{
int l=0;
for (const auto& face_cells_vector : face_cells_map) {
- const auto& cells_vector = face_cells_vector.second;
- for (unsigned short lj=0; lj<cells_vector.size(); ++lj) {
- unsigned int cell_number = cells_vector[lj].first;
- unsigned short cell_local_face = cells_vector[lj].second;
- cell_faces(cell_number,cell_local_face) = l;
- }
- ++l;
+ const auto& cells_vector = face_cells_vector.second;
+ for (unsigned short lj=0; lj<cells_vector.size(); ++lj) {
+ unsigned int cell_number = cells_vector[lj].first;
+ unsigned short cell_local_face = cells_vector[lj].second;
+ cell_faces(cell_number,cell_local_face) = l;
+ }
+ ++l;
}
}
m_cell_faces = cell_faces;
Kokkos::View<unsigned short**> face_cell_local_face("face_cell_local_face",
- m_number_of_faces, m_max_nb_node_per_cell);
+ m_number_of_faces, m_max_nb_node_per_cell);
{
int l=0;
for (const auto& face_cells_vector : face_cells_map) {
- const auto& cells_vector = face_cells_vector.second;
- for (unsigned short lj=0; lj<cells_vector.size(); ++lj) {
- unsigned short cell_local_face = cells_vector[lj].second;
- face_cell_local_face(l,lj) = cell_local_face;
- }
- ++l;
+ const auto& cells_vector = face_cells_vector.second;
+ for (unsigned short lj=0; lj<cells_vector.size(); ++lj) {
+ unsigned short cell_local_face = cells_vector[lj].second;
+ face_cell_local_face(l,lj) = cell_local_face;
+ }
+ ++l;
}
}
m_face_cell_local_face = face_cell_local_face;
}
-
-public:
+ public:
const size_t& numberOfNodes() const
{
return m_number_of_nodes;
@@ -200,7 +270,7 @@ public:
{
return m_max_nb_node_per_cell;
}
-
+
const Kokkos::View<const unsigned int**> cellNodes() const
{
return m_cell_nodes;
@@ -262,7 +332,7 @@ public:
}
unsigned int getFaceNumber(const unsigned int node0_id,
- const unsigned int node1_id) const
+ const unsigned int node1_id) const
{
#warning rewrite
const unsigned int n0_id = std::min(node0_id, node1_id);
@@ -270,7 +340,7 @@ public:
// worst code ever
for (unsigned int l=0; l<m_number_of_faces; ++l) {
if ((m_face_nodes(l,0) == n0_id) and (m_face_nodes(l,1) == n1_id)) {
- return l;
+ return l;
}
}
std::cerr << "Face not found!\n";
@@ -278,27 +348,17 @@ public:
return 0;
}
- size_t numberOfFaceBoundaries() const
+ EntityListManager& entityListManager()
{
- return m_faces_boundary_list.size();
+ return m_entity_list_manager;
}
- const FacesBoundary& facesBoundary(const size_t& i) const
+ const EntityListManager& entityListManager() const
{
- return m_faces_boundary_list[i];
- }
-
- void addFaceBoundary(const unsigned int& ref,
- const std::vector<unsigned int>& face_id_vector)
- {
- Kokkos::View<unsigned int*> faces_boundary("faces_boundary", face_id_vector.size());
- Kokkos::parallel_for(face_id_vector.size(), KOKKOS_LAMBDA(const int& l) {
- faces_boundary[l] = face_id_vector[l];
- });
- m_faces_boundary_list.push_back(std::make_pair(ref, faces_boundary));
+ return m_entity_list_manager;
}
- size_t numberOfNodeBoundaries() const
+size_t numberOfNodeBoundaries() const
{
return m_nodes_boundary_list.size();
}
@@ -308,8 +368,8 @@ public:
return m_nodes_boundary_list[i];
}
- void addNodeBoundary(const unsigned int& ref,
- const Kokkos::View<const unsigned int*> nodes_boundary)
+ void addNodeBoundary(const RefId& ref,
+ const Kokkos::View<const unsigned int*> nodes_boundary)
{
m_nodes_boundary_list.push_back(std::make_pair(ref, nodes_boundary));
}
@@ -317,7 +377,7 @@ public:
Connectivity2D(const Connectivity2D&) = delete;
Connectivity2D(const Kokkos::View<const unsigned short*> cell_nb_nodes,
- const Kokkos::View<const unsigned int**> cell_nodes)
+ const Kokkos::View<const unsigned int**> cell_nodes)
: m_number_of_cells (cell_nodes.extent(0)),
m_cell_nb_nodes(cell_nb_nodes),
m_cell_nodes (cell_nodes)
@@ -325,21 +385,21 @@ public:
{
Kokkos::View<double*> inv_cell_nb_nodes("inv_cell_nb_nodes", m_number_of_cells);
Kokkos::parallel_for(m_number_of_cells, KOKKOS_LAMBDA(const int&j){
- inv_cell_nb_nodes[j] = 1./m_cell_nb_nodes[j];
- });
+ inv_cell_nb_nodes[j] = 1./m_cell_nb_nodes[j];
+ });
m_inv_cell_nb_nodes = inv_cell_nb_nodes;
}
assert(m_number_of_cells>0);
ConnectivityUtils utils;
utils.computeNodeCellConnectivity(m_cell_nodes,
- m_cell_nb_nodes,
- m_number_of_cells,
- m_max_nb_node_per_cell,
- m_number_of_nodes,
- m_node_nb_cells,
- m_node_cells,
- m_node_cell_local_node);
+ m_cell_nb_nodes,
+ m_number_of_cells,
+ m_max_nb_node_per_cell,
+ m_number_of_nodes,
+ m_node_nb_cells,
+ m_node_cells,
+ m_node_cell_local_node);
this->_computeFaceCellConnectivities();
//this->_computeNodeFaceConnectivities();
diff --git a/src/mesh/GmshReader.cpp b/src/mesh/GmshReader.cpp
index ba398607b..53c704e13 100644
--- a/src/mesh/GmshReader.cpp
+++ b/src/mesh/GmshReader.cpp
@@ -9,8 +9,9 @@
#include <Connectivity2D.hpp>
#include <Mesh.hpp>
#include <map>
+#include <regex>
-template<typename T>
+template<typename T>
inline std::string stringify(const T & t)
{
std::ostringstream oss;
@@ -18,7 +19,7 @@ inline std::string stringify(const T & t)
return oss.str();
}
-template<>
+template<>
inline std::string stringify<std::string>(const std::string& t)
{
return t;
@@ -26,7 +27,7 @@ inline std::string stringify<std::string>(const std::string& t)
class ErrorHandler
{
-public:
+ public:
enum Type {
asked, /**< execution request by the user*/
compilation, /**< syntax error in a language */
@@ -34,49 +35,49 @@ public:
unexpected /**< Unexpected execution error */
};
-private:
+ private:
const std::string __filename; /**< The source file name where the error occured */
const size_t __lineNumber; /**< The line number where exception was raised */
const std::string __errorMessage; /**< The reporting message */
const Type __type; /**< the type of the error */
-public:
- /**
+ public:
+ /**
* Prints the error message
- *
+ *
*/
virtual void writeErrorMessage();
- /**
+ /**
* The copy constructor
- *
+ *
* @param e an handled error
*/
ErrorHandler(const ErrorHandler& e)
- : __filename(e.__filename),
- __lineNumber(e.__lineNumber),
- __errorMessage(e.__errorMessage),
- __type(e.__type)
+ : __filename(e.__filename),
+ __lineNumber(e.__lineNumber),
+ __errorMessage(e.__errorMessage),
+ __type(e.__type)
{
;
}
- /**
+ /**
* Constructor
- *
+ *
* @param filename the source file name
* @param lineNumber the source line
* @param errorMessage the reported message
* @param type the type of the error
*/
ErrorHandler(const std::string& filename,
- const size_t& lineNumber,
- const std::string& errorMessage,
- const Type& type);
+ const size_t& lineNumber,
+ const std::string& errorMessage,
+ const Type& type);
- /**
+ /**
* The destructor
- *
+ *
*/
virtual ~ErrorHandler()
{
@@ -86,205 +87,205 @@ public:
void ErrorHandler::writeErrorMessage()
{
switch(__type) {
- case asked: {
- std::cerr << "\nremark: exit command explicitly called\n";
- }
- case normal: {
- std::cerr << '\n' << __filename << ':' << __lineNumber
- << ":remark: emitted the following message\n";
- std::cerr << "error: " << __errorMessage << '\n';
- break;
- }
- case compilation: {
- std::cerr << "\nline " << __lineNumber << ':' << __errorMessage << '\n';
- break;
- }
- case unexpected: {
- std::cerr << '\n' << __filename << ':' << __lineNumber << ":\n" << __errorMessage << '\n';
- std::cerr << "\nUNEXPECTED ERROR: this should not occure, please report it\n";
- std::cerr << "\nBUG REPORT: Please send bug reports to:\n"
- << " ff3d-dev@nongnu.org or freefem@ann.jussieu.fr\n"
- << "or better, use the Bug Tracking System:\n"
- << " http://savannah.nongnu.org/bugs/?group=ff3d\n";
- break;
- }
- default: {
- std::cerr << __filename << ':' << __lineNumber << ": " << __errorMessage << '\n';
- std::cerr << __FILE__ << ':' << __LINE__ << ":remark: error type not implemented!\n";
- }
+ case asked: {
+ std::cerr << "\nremark: exit command explicitly called\n";
+ }
+ case normal: {
+ std::cerr << '\n' << __filename << ':' << __lineNumber
+ << ":remark: emitted the following message\n";
+ std::cerr << "error: " << __errorMessage << '\n';
+ break;
+ }
+ case compilation: {
+ std::cerr << "\nline " << __lineNumber << ':' << __errorMessage << '\n';
+ break;
+ }
+ case unexpected: {
+ std::cerr << '\n' << __filename << ':' << __lineNumber << ":\n" << __errorMessage << '\n';
+ std::cerr << "\nUNEXPECTED ERROR: this should not occure, please report it\n";
+ std::cerr << "\nBUG REPORT: Please send bug reports to:\n"
+ << " ff3d-dev@nongnu.org or freefem@ann.jussieu.fr\n"
+ << "or better, use the Bug Tracking System:\n"
+ << " http://savannah.nongnu.org/bugs/?group=ff3d\n";
+ break;
+ }
+ default: {
+ std::cerr << __filename << ':' << __lineNumber << ": " << __errorMessage << '\n';
+ std::cerr << __FILE__ << ':' << __LINE__ << ":remark: error type not implemented!\n";
+ }
}
}
ErrorHandler::
ErrorHandler(const std::string& filename,
- const size_t& lineNumber,
- const std::string& errorMessage,
- const Type& type)
- : __filename(filename),
- __lineNumber(lineNumber),
- __errorMessage(errorMessage),
- __type(type)
+ const size_t& lineNumber,
+ const std::string& errorMessage,
+ const Type& type)
+ : __filename(filename),
+ __lineNumber(lineNumber),
+ __errorMessage(errorMessage),
+ __type(type)
{
;
}
GmshReader::GmshReader(const std::string& filename)
- : m_filename(filename)
+ : m_filename(filename)
{
try {
- m_fin.open(m_filename);
- if (not m_fin) {
- std::cerr << "cannot read file '" << m_filename << "'\n";
- std::exit(0);
- }
-
- // gmsh 2.2 format keywords
- __keywordList["$MeshFormat"] = MESHFORMAT;
- __keywordList["$EndMeshFormat"] = ENDMESHFORMAT;
- __keywordList["$Nodes"] = NODES;
- __keywordList["$EndNodes"] = ENDNODES;
- __keywordList["$Elements"] = ELEMENTS;
- __keywordList["$EndElements"] = ENDELEMENTS;
- __keywordList["$PhysicalNames"] = PHYSICALNAMES;
- __keywordList["$EndPhysicalNames"] = ENDPHYSICALNAMES;
-
- __numberOfPrimitiveNodes.resize(16);
- __numberOfPrimitiveNodes[ 0] = 2; // edge
- __numberOfPrimitiveNodes[ 1] = 3; // triangle
- __numberOfPrimitiveNodes[ 2] = 4; // quadrangle
- __numberOfPrimitiveNodes[ 3] = 4; // Tetrahedron
- __numberOfPrimitiveNodes[ 4] = 8; // Hexaredron
- __numberOfPrimitiveNodes[ 5] = 6; // Prism
- __numberOfPrimitiveNodes[ 6] = 5; // Pyramid
- __numberOfPrimitiveNodes[ 7] = 3; // second order edge
- __numberOfPrimitiveNodes[ 8] = 6; // second order triangle
- __numberOfPrimitiveNodes[ 9] = 9; // second order quadrangle
- __numberOfPrimitiveNodes[10] = 10; // second order tetrahedron
- __numberOfPrimitiveNodes[11] = 27; // second order hexahedron
- __numberOfPrimitiveNodes[12] = 18; // second order prism
- __numberOfPrimitiveNodes[13] = 14; // second order pyramid
- __numberOfPrimitiveNodes[14] = 1; // point
-
- __primitivesNames[0] = "edges";
- __supportedPrimitives[0] = true;
- __primitivesNames[1] = "triangles";
- __supportedPrimitives[1] = true;
- __primitivesNames[2] = "quadrangles";
- __supportedPrimitives[2] = true;
- __primitivesNames[3] = "tetrahedra";
- __supportedPrimitives[3] = true;
- __primitivesNames[4] = "hexahedra";
- __supportedPrimitives[4] = true;
- __primitivesNames[5] = "prisms";
- __supportedPrimitives[5] = false;
- __primitivesNames[6] = "pyramids";
- __supportedPrimitives[6] = false;
- __primitivesNames[7] = "second order edges";
- __supportedPrimitives[7] = false;
- __primitivesNames[8] = "second order triangles";
- __supportedPrimitives[8] = false;
- __primitivesNames[9] = "second order quadrangles";
- __supportedPrimitives[9] = false;
- __primitivesNames[10] = "second order tetrahedra";
- __supportedPrimitives[10]= false;
- __primitivesNames[11] = "second order hexahedra";
- __supportedPrimitives[11]= false;
- __primitivesNames[12] = "second order prisms";
- __supportedPrimitives[12]= false;
- __primitivesNames[13] = "second order pyramids";
- __supportedPrimitives[13]= false;
- __primitivesNames[14] = "point";
- __supportedPrimitives[14]= false;
-
- std::cout << "Reading file '" << m_filename << "'\n";
-
- // Getting vertices list
- GmshReader::Keyword kw = this->__nextKeyword();
- switch(kw.second) {
- // case NOD: {
- // this->__readGmsh1();
- // break;
- // }
- case MESHFORMAT: {
- double fileVersion = this->_getReal();
- if (fileVersion != 2.2) {
- throw ErrorHandler(__FILE__,__LINE__,
- "Cannot read Gmsh format '"+stringify(fileVersion)+"'",
- ErrorHandler::normal);
- }
- int fileType = this->_getInteger();
- __binary = (fileType == 1);
- if ((fileType < 0) or (fileType > 1)) {
- throw ErrorHandler(__FILE__,__LINE__,
- "Cannot read Gmsh file type '"+stringify(fileType)+"'",
- ErrorHandler::normal);
- }
-
- int dataSize = this->_getInteger();
- if (dataSize != sizeof(double)) {
- throw ErrorHandler(__FILE__,__LINE__,
- "Data size not supported '"+stringify(dataSize)+"'",
- ErrorHandler::normal);
- }
-
- if (__binary) {
-// int one=0;
-
-// fseek(__ifh,1L,SEEK_CUR);
-// fread(reinterpret_cast<char*>(&one),sizeof(int),1,__ifh);
-
-// if (one == 1) {
-// __convertEndian = false;
-// } else {
-// invertEndianess(one);
-
-// if (one == 1) {
-// __convertEndian = true;
-// } else {
-// throw ErrorHandler(__FILE__,__LINE__,
-// "Cannot determine endianess",
-// ErrorHandler::normal);
-// }
-// }
-
-// std::cout << "- Binary format: ";
-// #ifdef WORDS_BIGENDIAN
-// if (not __convertEndian) {
-// std::cout << "big endian\n";
-// } else {
-// std::cout << "little endian\n";
-// }
-// #else // WORDS_BIGENDIAN
-// if (not __convertEndian) {
-// std::cout << "little endian\n";
-// } else {
-// std::cout << "big endian\n";
-// }
-// #endif // WORDS_BIGENDIAN
- }
-
- kw = this->__nextKeyword();
- if (kw.second != ENDMESHFORMAT) {
- throw ErrorHandler(__FILE__,__LINE__,
- "reading file '"+m_filename
- +"': expecting $EndMeshFormat, '"+kw.first+"' was found",
- ErrorHandler::normal);
+ m_fin.open(m_filename);
+ if (not m_fin) {
+ std::cerr << "cannot read file '" << m_filename << "'\n";
+ std::exit(0);
+ }
+
+ // gmsh 2.2 format keywords
+ __keywordList["$MeshFormat"] = MESHFORMAT;
+ __keywordList["$EndMeshFormat"] = ENDMESHFORMAT;
+ __keywordList["$Nodes"] = NODES;
+ __keywordList["$EndNodes"] = ENDNODES;
+ __keywordList["$Elements"] = ELEMENTS;
+ __keywordList["$EndElements"] = ENDELEMENTS;
+ __keywordList["$PhysicalNames"] = PHYSICALNAMES;
+ __keywordList["$EndPhysicalNames"] = ENDPHYSICALNAMES;
+
+ __numberOfPrimitiveNodes.resize(16);
+ __numberOfPrimitiveNodes[ 0] = 2; // edge
+ __numberOfPrimitiveNodes[ 1] = 3; // triangle
+ __numberOfPrimitiveNodes[ 2] = 4; // quadrangle
+ __numberOfPrimitiveNodes[ 3] = 4; // Tetrahedron
+ __numberOfPrimitiveNodes[ 4] = 8; // Hexaredron
+ __numberOfPrimitiveNodes[ 5] = 6; // Prism
+ __numberOfPrimitiveNodes[ 6] = 5; // Pyramid
+ __numberOfPrimitiveNodes[ 7] = 3; // second order edge
+ __numberOfPrimitiveNodes[ 8] = 6; // second order triangle
+ __numberOfPrimitiveNodes[ 9] = 9; // second order quadrangle
+ __numberOfPrimitiveNodes[10] = 10; // second order tetrahedron
+ __numberOfPrimitiveNodes[11] = 27; // second order hexahedron
+ __numberOfPrimitiveNodes[12] = 18; // second order prism
+ __numberOfPrimitiveNodes[13] = 14; // second order pyramid
+ __numberOfPrimitiveNodes[14] = 1; // point
+
+ __primitivesNames[0] = "edges";
+ __supportedPrimitives[0] = true;
+ __primitivesNames[1] = "triangles";
+ __supportedPrimitives[1] = true;
+ __primitivesNames[2] = "quadrangles";
+ __supportedPrimitives[2] = true;
+ __primitivesNames[3] = "tetrahedra";
+ __supportedPrimitives[3] = true;
+ __primitivesNames[4] = "hexahedra";
+ __supportedPrimitives[4] = true;
+ __primitivesNames[5] = "prisms";
+ __supportedPrimitives[5] = false;
+ __primitivesNames[6] = "pyramids";
+ __supportedPrimitives[6] = false;
+ __primitivesNames[7] = "second order edges";
+ __supportedPrimitives[7] = false;
+ __primitivesNames[8] = "second order triangles";
+ __supportedPrimitives[8] = false;
+ __primitivesNames[9] = "second order quadrangles";
+ __supportedPrimitives[9] = false;
+ __primitivesNames[10] = "second order tetrahedra";
+ __supportedPrimitives[10]= false;
+ __primitivesNames[11] = "second order hexahedra";
+ __supportedPrimitives[11]= false;
+ __primitivesNames[12] = "second order prisms";
+ __supportedPrimitives[12]= false;
+ __primitivesNames[13] = "second order pyramids";
+ __supportedPrimitives[13]= false;
+ __primitivesNames[14] = "point";
+ __supportedPrimitives[14]= true;
+
+ std::cout << "Reading file '" << m_filename << "'\n";
+
+ // Getting vertices list
+ GmshReader::Keyword kw = this->__nextKeyword();
+ switch(kw.second) {
+ // case NOD: {
+ // this->__readGmsh1();
+ // break;
+ // }
+ case MESHFORMAT: {
+ double fileVersion = this->_getReal();
+ if (fileVersion != 2.2) {
+ throw ErrorHandler(__FILE__,__LINE__,
+ "Cannot read Gmsh format '"+stringify(fileVersion)+"'",
+ ErrorHandler::normal);
+ }
+ int fileType = this->_getInteger();
+ __binary = (fileType == 1);
+ if ((fileType < 0) or (fileType > 1)) {
+ throw ErrorHandler(__FILE__,__LINE__,
+ "Cannot read Gmsh file type '"+stringify(fileType)+"'",
+ ErrorHandler::normal);
+ }
+
+ int dataSize = this->_getInteger();
+ if (dataSize != sizeof(double)) {
+ throw ErrorHandler(__FILE__,__LINE__,
+ "Data size not supported '"+stringify(dataSize)+"'",
+ ErrorHandler::normal);
+ }
+
+ if (__binary) {
+ // int one=0;
+
+ // fseek(__ifh,1L,SEEK_CUR);
+ // fread(reinterpret_cast<char*>(&one),sizeof(int),1,__ifh);
+
+ // if (one == 1) {
+ // __convertEndian = false;
+ // } else {
+ // invertEndianess(one);
+
+ // if (one == 1) {
+ // __convertEndian = true;
+ // } else {
+ // throw ErrorHandler(__FILE__,__LINE__,
+ // "Cannot determine endianess",
+ // ErrorHandler::normal);
+ // }
+ // }
+
+ // std::cout << "- Binary format: ";
+ // #ifdef WORDS_BIGENDIAN
+ // if (not __convertEndian) {
+ // std::cout << "big endian\n";
+ // } else {
+ // std::cout << "little endian\n";
+ // }
+ // #else // WORDS_BIGENDIAN
+ // if (not __convertEndian) {
+ // std::cout << "little endian\n";
+ // } else {
+ // std::cout << "big endian\n";
+ // }
+ // #endif // WORDS_BIGENDIAN
+ }
+
+ kw = this->__nextKeyword();
+ if (kw.second != ENDMESHFORMAT) {
+ throw ErrorHandler(__FILE__,__LINE__,
+ "reading file '"+m_filename
+ +"': expecting $EndMeshFormat, '"+kw.first+"' was found",
+ ErrorHandler::normal);
+ }
+
+ this->__readGmshFormat2_2();
+
+ break;
+ }
+ default: {
+ throw ErrorHandler(__FILE__,__LINE__,
+ "cannot determine format version of '"+m_filename+"'",
+ ErrorHandler::normal);
+ }
}
- this->__readGmshFormat2_2();
-
- break;
- }
- default: {
- throw ErrorHandler(__FILE__,__LINE__,
- "cannot determine format version of '"+m_filename+"'",
- ErrorHandler::normal);
- }
- }
-
- this->__proceedData();
- // this->__createMesh();
+ this->__proceedData();
+ // this->__createMesh();
}
catch(ErrorHandler e) {
e.writeErrorMessage();
@@ -298,14 +299,14 @@ void GmshReader::__readVertices()
std::cout << "- Number Of Vertices: " << numberOfVerices << '\n';
if (numberOfVerices<0) {
throw ErrorHandler(__FILE__,__LINE__,
- "reading file '"+this->m_filename
- +"': number of vertices is negative",
- ErrorHandler::normal);
+ "reading file '"+this->m_filename
+ +"': number of vertices is negative",
+ ErrorHandler::normal);
}
__verticesNumbers.resize(numberOfVerices);
__vertices = Kokkos::View<R3*>("vertices",numberOfVerices);
-
+
if (not __binary) {
for (int i=0; i<numberOfVerices; ++i) {
__verticesNumbers[i] = this->_getInteger();
@@ -313,7 +314,7 @@ void GmshReader::__readVertices()
const double y = this->_getReal();
const double z = this->_getReal();
__vertices[i] = TinyVector<3, double>(x,y,z);
- }
+ }
} else {
// fseek(m_fin.file()__ifh,1L,SEEK_CUR);
// for (int i=0; i<numberOfVerices; ++i) {
@@ -394,9 +395,9 @@ GmshReader::__readElements2_2()
std::cout << "- Number Of Elements: " << numberOfElements << '\n';
if (numberOfElements<0) {
throw ErrorHandler(__FILE__,__LINE__,
- "reading file '"+m_filename
- +"': number of elements is negative",
- ErrorHandler::normal);
+ "reading file '"+m_filename
+ +"': number of elements is negative",
+ ErrorHandler::normal);
}
__elementType.resize(numberOfElements);
@@ -409,22 +410,22 @@ GmshReader::__readElements2_2()
const int elementType = this->_getInteger()-1;
if ((elementType < 0) or (elementType > 14)) {
- throw ErrorHandler(__FILE__,__LINE__,
- "reading file '"+m_filename
- +"': unknown element type '"+stringify(elementType)+"'",
- ErrorHandler::normal);
+ throw ErrorHandler(__FILE__,__LINE__,
+ "reading file '"+m_filename
+ +"': unknown element type '"+stringify(elementType)+"'",
+ ErrorHandler::normal);
}
__elementType[i] = elementType;
const int numberOfTags = this->_getInteger();
__references[i] = this->_getInteger(); // physical reference
for (int tag=1; tag<numberOfTags; ++tag) {
- this->_getInteger(); // drops remaining tags
+ this->_getInteger(); // drops remaining tags
}
const int numberOfVertices = __numberOfPrimitiveNodes[elementType];
__elementVertices[i].resize(numberOfVertices);
for (int j=0; j<numberOfVertices; ++j) {
- __elementVertices[i][j] = this->_getInteger();
+ __elementVertices[i][j] = this->_getInteger();
}
}
} else {
@@ -502,7 +503,19 @@ __readPhysicalNames2_2()
const int physical_number = this->_getInteger();
std::string physical_name;
m_fin >> physical_name;
- std::cout << "- " << rang::style::bold << physical_name << " <-> " << physical_number << " for dimension " << physical_dimension << rang::style::reset << " [IGNORED]\n";
+ physical_name = std::regex_replace(physical_name, std::regex("(\")"), "");
+ std::cout << "- " << rang::style::bold << physical_name
+ << " <-> " << physical_number << " for dimension "
+ << physical_dimension << rang::style::reset
+ << " [" << rang::fg::yellow << "IGNORED" << rang::style::reset << "]\n";
+ PhysicalRefId physical_ref_id(physical_dimension, RefId(physical_number, physical_name));
+ auto searched_physical_ref_id = m_physical_ref_map.find(physical_number);
+ if (searched_physical_ref_id != m_physical_ref_map.end()) {
+ std::cerr << "Physical reference id '" << physical_ref_id
+ << "' already defined as '" << searched_physical_ref_id->second << "'!";
+ std::exit(0);
+ }
+ m_physical_ref_map[physical_number] = physical_ref_id;
}
}
@@ -521,76 +534,74 @@ GmshReader::__proceedData()
for (size_t i=0; i<elementNumber.dimension(); ++i) {
if (elementNumber[i]>0) {
std::cout << " - Number of "
- << __primitivesNames[i]
- << ": " << elementNumber[i];
+ << __primitivesNames[i]
+ << ": " << elementNumber[i];
if (not(__supportedPrimitives[i])) {
- std::cout << " [IGNORED]";
+ 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(); ++iref) {
- if (iref != elementReferences[i].begin()) std::cout << ',';
- std::cout << *iref;
- }
-
- switch (i) {
- // Supported entities
- case 0: { // edges
- __edges
- = Kokkos::View<Edge*>("edges",
- elementNumber[i]);
- __edges_ref.resize(elementNumber[i]);
- break;
- }
- case 1: { // triangles
- __triangles
- = Kokkos::View<Triangle*>("triangles",
- elementNumber[i]);
- __triangles_ref.resize(elementNumber[i]);
- break;
- }
- case 2: { // quadrangles
- __quadrangles
- = Kokkos::View<Quadrangle*>("quadrangles",
- elementNumber[i]);
- __quadrangles_ref.resize(elementNumber[i]);
- break;
- }
- case 3: { // tetrahedra
- __tetrahedra
- = Kokkos::View<Tetrahedron*>("tetrahedra",
- elementNumber[i]);
- __tetrahedra_ref.resize(elementNumber[i]);
- break;
- }
- case 4: {// hexahedra
- __hexahedra
- = Kokkos::View<Hexahedron*>("hexahedra",
- elementNumber[i]);
- __hexahedra_ref.resize(elementNumber[i]);
- }
- // Ignored entities
- case 14: {// point
- break;
- }
- // Unsupported entities
- case 5: // prism
- case 6: // pyramid
- case 7: // second order edge
- case 8: // second order triangle
- case 9: // second order quadrangle
- case 10: // second order tetrahedron
- case 11: // second order hexahedron
- case 12: // second order prism
- case 13: // second order pyramid
- default: {
- throw ErrorHandler(__FILE__,__LINE__,
- "reading file '"+m_filename
- +"': ff3d cannot read this kind of element",
- ErrorHandler::normal);
- }
- }
+ std::cout << " references: ";
+ for(std::set<size_t>::const_iterator iref
+ = elementReferences[i].begin();
+ iref != elementReferences[i].end(); ++iref) {
+ if (iref != elementReferences[i].begin()) std::cout << ',';
+ std::cout << *iref;
+ }
+
+ switch (i) {
+ // Supported entities
+ case 0: { // edges
+ __edges = Kokkos::View<Edge*>("edges",
+ elementNumber[i]);
+ __edges_ref.resize(elementNumber[i]);
+ break;
+ }
+ case 1: { // triangles
+ __triangles = Kokkos::View<Triangle*>("triangles",
+ elementNumber[i]);
+ __triangles_ref.resize(elementNumber[i]);
+ break;
+ }
+ case 2: { // quadrangles
+ __quadrangles = Kokkos::View<Quadrangle*>("quadrangles",
+ elementNumber[i]);
+ __quadrangles_ref.resize(elementNumber[i]);
+ break;
+ }
+ case 3: { // tetrahedra
+ __tetrahedra = Kokkos::View<Tetrahedron*>("tetrahedra",
+ elementNumber[i]);
+ __tetrahedra_ref.resize(elementNumber[i]);
+ break;
+ }
+ case 4: {// hexahedra
+ __hexahedra = Kokkos::View<Hexahedron*>("hexahedra",
+ elementNumber[i]);
+ __hexahedra_ref.resize(elementNumber[i]);
+ }
+ // Ignored entities
+ case 14: {// point
+ __points = Kokkos::View<Point*>("points",
+ elementNumber[i]);
+ __points_ref.resize(elementNumber[i]);
+ break;
+ }
+ // Unsupported entities
+ case 5: // prism
+ case 6: // pyramid
+ case 7: // second order edge
+ case 8: // second order triangle
+ case 9: // second order quadrangle
+ case 10: // second order tetrahedron
+ case 11: // second order hexahedron
+ case 12: // second order prism
+ case 13: // second order pyramid
+ default: {
+ throw ErrorHandler(__FILE__,__LINE__,
+ "reading file '"+m_filename
+ +"': ff3d cannot read this kind of element",
+ ErrorHandler::normal);
+ }
+ }
}
std::cout << '\n';
}
@@ -607,16 +618,16 @@ GmshReader::__proceedData()
if (minNumber<0) {
throw ErrorHandler(__FILE__,__LINE__,
- "reading file '"+m_filename
- +"': ff3d does not implement negative vertices number",
- ErrorHandler::normal);
+ "reading file '"+m_filename
+ +"': ff3d does not implement negative vertices number",
+ ErrorHandler::normal);
}
- // A value of -1 means that the vertex is unknown
+ // A value of -1 means that the vertex is unknown
__verticesCorrepondance.resize(maxNumber+1,-1);
for (size_t i=0; i<__verticesNumbers.size(); ++i) {
-
+
__verticesCorrepondance[__verticesNumbers[i]] = i;
}
@@ -628,125 +639,130 @@ GmshReader::__proceedData()
std::vector<int>& elementVertices = __elementVertices[i];
switch (__elementType[i]) {
// Supported entities
- case 0: { // edge
- int& edgeNumber = elementNumber[0];
- const int a = __verticesCorrepondance[elementVertices[0]];
- const int b = __verticesCorrepondance[elementVertices[1]];
- if ((a<0) or (b<0)) {
- throw ErrorHandler(__FILE__,__LINE__,
- "reading file '"+m_filename
- +"': error reading element "+stringify(i)
- +" [bad vertices definition]",
- ErrorHandler::normal);
+ case 0: { // edge
+ int& edgeNumber = elementNumber[0];
+ const int a = __verticesCorrepondance[elementVertices[0]];
+ const int b = __verticesCorrepondance[elementVertices[1]];
+ if ((a<0) or (b<0)) {
+ throw ErrorHandler(__FILE__,__LINE__,
+ "reading file '"+m_filename
+ +"': error reading element "+stringify(i)
+ +" [bad vertices definition]",
+ ErrorHandler::normal);
+ }
+ __edges[edgeNumber]
+ = Edge(a, b);
+ __edges_ref[edgeNumber] = __references[i];
+ edgeNumber++; // one more edge
+ break;
}
- __edges[edgeNumber]
- = Edge(a, b);
- __edges_ref[edgeNumber] = __references[i];
- edgeNumber++; // one more edge
- break;
- }
- case 1: { // triangles
- int& triangleNumber = elementNumber[1];
-
- const int a = __verticesCorrepondance[elementVertices[0]];
- const int b = __verticesCorrepondance[elementVertices[1]];
- const int c = __verticesCorrepondance[elementVertices[2]];
- if ((a<0) or (b<0) or (c<0)) {
- throw ErrorHandler(__FILE__,__LINE__,
- "reading file '"+m_filename
- +"': error reading element "+stringify(i)
- +" [bad vertices definition]",
- ErrorHandler::normal);
+ case 1: { // triangles
+ int& triangleNumber = elementNumber[1];
+
+ const int a = __verticesCorrepondance[elementVertices[0]];
+ const int b = __verticesCorrepondance[elementVertices[1]];
+ const int c = __verticesCorrepondance[elementVertices[2]];
+ if ((a<0) or (b<0) or (c<0)) {
+ throw ErrorHandler(__FILE__,__LINE__,
+ "reading file '"+m_filename
+ +"': error reading element "+stringify(i)
+ +" [bad vertices definition]",
+ ErrorHandler::normal);
+ }
+ __triangles[triangleNumber]
+ = Triangle(a, b, c);
+ __triangles_ref[triangleNumber] = __references[i];
+ triangleNumber++; // one more triangle
+ break;
}
- __triangles[triangleNumber]
- = Triangle(a, b, c);
- __triangles_ref[triangleNumber] = __references[i];
- triangleNumber++; // one more triangle
- break;
- }
- case 2: { // quadrangle
- int& quadrilateralNumber = elementNumber[2];
-
- const int a = __verticesCorrepondance[elementVertices[0]];
- const int b = __verticesCorrepondance[elementVertices[1]];
- const int c = __verticesCorrepondance[elementVertices[2]];
- const int d = __verticesCorrepondance[elementVertices[3]];
- if ((a<0) or (b<0) or (c<0) or (d<0)) {
- throw ErrorHandler(__FILE__,__LINE__,
- "reading file '"+m_filename
- +"': error reading element "+stringify(i)
- +" [bad vertices definition]",
- ErrorHandler::normal);
+ case 2: { // quadrangle
+ int& quadrilateralNumber = elementNumber[2];
+
+ const int a = __verticesCorrepondance[elementVertices[0]];
+ const int b = __verticesCorrepondance[elementVertices[1]];
+ const int c = __verticesCorrepondance[elementVertices[2]];
+ const int d = __verticesCorrepondance[elementVertices[3]];
+ if ((a<0) or (b<0) or (c<0) or (d<0)) {
+ throw ErrorHandler(__FILE__,__LINE__,
+ "reading file '"+m_filename
+ +"': error reading element "+stringify(i)
+ +" [bad vertices definition]",
+ ErrorHandler::normal);
+ }
+ __quadrangles[quadrilateralNumber] = Quadrangle(a,b,c,d);
+ __quadrangles_ref[quadrilateralNumber] = __references[i];
+ quadrilateralNumber++; // one more quadrangle
+ break;
}
- __quadrangles[quadrilateralNumber] = Quadrangle(a,b,c,d);
- __quadrangles_ref[quadrilateralNumber] = __references[i];
- quadrilateralNumber++; // one more quadrangle
- break;
- }
- case 3: { // tetrahedra
- int& tetrahedronNumber = elementNumber[3];
-
- const int a = __verticesCorrepondance[elementVertices[0]];
- const int b = __verticesCorrepondance[elementVertices[1]];
- const int c = __verticesCorrepondance[elementVertices[2]];
- const int d = __verticesCorrepondance[elementVertices[3]];
- if ((a<0) or (b<0) or (c<0) or (d<0)) {
- throw ErrorHandler(__FILE__,__LINE__,
- "reading file '"+m_filename
- +"': error reading element "+stringify(i)
- +" [bad vertices definition]",
- ErrorHandler::normal);
+ case 3: { // tetrahedra
+ int& tetrahedronNumber = elementNumber[3];
+
+ const int a = __verticesCorrepondance[elementVertices[0]];
+ const int b = __verticesCorrepondance[elementVertices[1]];
+ const int c = __verticesCorrepondance[elementVertices[2]];
+ const int d = __verticesCorrepondance[elementVertices[3]];
+ if ((a<0) or (b<0) or (c<0) or (d<0)) {
+ throw ErrorHandler(__FILE__,__LINE__,
+ "reading file '"+m_filename
+ +"': error reading element "+stringify(i)
+ +" [bad vertices definition]",
+ ErrorHandler::normal);
+ }
+ __tetrahedra[tetrahedronNumber] = Tetrahedron(a, b, c, d);
+ __tetrahedra_ref[tetrahedronNumber] = __references[i];
+ tetrahedronNumber++; // one more tetrahedron
+ break;
}
- __tetrahedra[tetrahedronNumber] = Tetrahedron(a, b, c, d);
- __tetrahedra_ref[tetrahedronNumber] = __references[i];
- tetrahedronNumber++; // one more tetrahedron
- break;
- }
- case 4: { // hexaredron
- int& hexahedronNumber = elementNumber[4];
-
- const int a = __verticesCorrepondance[elementVertices[0]];
- const int b = __verticesCorrepondance[elementVertices[1]];
- const int c = __verticesCorrepondance[elementVertices[2]];
- const int d = __verticesCorrepondance[elementVertices[3]];
- const int e = __verticesCorrepondance[elementVertices[4]];
- const int f = __verticesCorrepondance[elementVertices[5]];
- const int g = __verticesCorrepondance[elementVertices[6]];
- const int h = __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 ErrorHandler(__FILE__,__LINE__,
- "reading file '"+m_filename
- +"': error reading element "+stringify(i)
- +" [bad vertices definition]",
- ErrorHandler::normal);
+ case 4: { // hexaredron
+ int& hexahedronNumber = elementNumber[4];
+
+ const int a = __verticesCorrepondance[elementVertices[0]];
+ const int b = __verticesCorrepondance[elementVertices[1]];
+ const int c = __verticesCorrepondance[elementVertices[2]];
+ const int d = __verticesCorrepondance[elementVertices[3]];
+ const int e = __verticesCorrepondance[elementVertices[4]];
+ const int f = __verticesCorrepondance[elementVertices[5]];
+ const int g = __verticesCorrepondance[elementVertices[6]];
+ const int h = __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 ErrorHandler(__FILE__,__LINE__,
+ "reading file '"+m_filename
+ +"': error reading element "+stringify(i)
+ +" [bad vertices definition]",
+ ErrorHandler::normal);
+ }
+ __hexahedra[hexahedronNumber]
+ = Hexahedron(a, b, c, d,
+ e, f, g, h);
+ __hexahedra_ref[hexahedronNumber] = __references[i];
+ hexahedronNumber++; // one more hexahedron
+ break;
+ }
+ // Unsupported entities
+ case 14:{// point
+ int& point_number = elementNumber[14];
+ const int a = __verticesCorrepondance[elementVertices[0]];
+ __points[point_number] = a;
+ __points_ref[point_number] = __references[i];
+ point_number++;
+ break;
+ }
+ case 5: // prism
+ case 6: // pyramid
+ case 7: // second order edge
+ case 8: // second order triangle
+ case 9: // second order quadrangle
+ case 10: // second order tetrahedron
+ case 11: // second order hexahedron
+ case 12: // second order prism
+ case 13:{// second order pyramid
+ }
+ default: {
+ throw ErrorHandler(__FILE__,__LINE__,
+ "reading file '"+m_filename
+ +"': ff3d cannot read this kind of element",
+ ErrorHandler::normal);
}
- __hexahedra[hexahedronNumber]
- = Hexahedron(a, b, c, d,
- e, f, g, h);
- __hexahedra_ref[hexahedronNumber] = __references[i];
- hexahedronNumber++; // one more hexahedron
- break;
- }
- // Unsupported entities
- case 14:{// point
- break;
- }
- case 5: // prism
- case 6: // pyramid
- case 7: // second order edge
- case 8: // second order triangle
- case 9: // second order quadrangle
- case 10: // second order tetrahedron
- case 11: // second order hexahedron
- case 12: // second order prism
- case 13:{// second order pyramid
- }
- default: {
- throw ErrorHandler(__FILE__,__LINE__,
- "reading file '"+m_filename
- +"': ff3d cannot read this kind of element",
- ErrorHandler::normal);
- }
}
}
@@ -790,7 +806,7 @@ GmshReader::__proceedData()
cell_nodes(j,1) = __triangles[j][1];
cell_nodes(j,2) = __triangles[j][2];
}
-
+
const size_t nb_quadrangles = __quadrangles.extent(0);
for (size_t j=0; j<nb_quadrangles; ++j) {
const size_t jq = j+nb_triangles;
@@ -812,43 +828,63 @@ GmshReader::__proceedData()
std::map<unsigned int, std::vector<unsigned int>> ref_faces_map;
for (unsigned int e=0; e<__edges.extent(0); ++e) {
const unsigned int edge_number = connectivity.getFaceNumber(__edges[e][0], __edges[e][1]);
- const unsigned int ref = __edges_ref[e];
+ const unsigned int& ref = __edges_ref[e];
ref_faces_map[ref].push_back(edge_number);
}
for (const auto& ref_face_list : ref_faces_map) {
- connectivity.addFaceBoundary(ref_face_list.first,
- ref_face_list.second);
+ Kokkos::View<unsigned int*> face_list("face_list", ref_face_list.second.size());
+ for (size_t j=0; j<ref_face_list.second.size(); ++j) {
+ face_list[j]=ref_face_list.second[j];
+ }
+ const PhysicalRefId& physical_ref_id = m_physical_ref_map.at(ref_face_list.first);
+ connectivity.entityListManager().addFaceIdArray(physical_ref_id.refId(), face_list);
}
const Kokkos::View<const unsigned int**> face_nodes = connectivity.faceNodes();
- for (size_t i=0; i<connectivity.numberOfFaceBoundaries(); ++i) {
- const Connectivity2D::FacesBoundary& faces_boundary
- = connectivity.facesBoundary(i);
+ std::cout << rang::fgB::red << "Creating nodes from faces for BC should not change connectivity definition" << rang::style::reset << '\n';
+
+ const EntityListManager& entity_list_manager = connectivity.entityListManager();
+ for (const auto& id_face_array : entity_list_manager.getFaceIdArraySet()) {
+ const RefId& ref = id_face_array.first;
+ const auto& face_ids = id_face_array.second;
- const unsigned int ref = faces_boundary.first;
std::set<unsigned int> node_id_set;
- const Kokkos::View<const unsigned int*> face_ids
- = faces_boundary.second;
-
+
for (unsigned int l=0; l<face_ids.extent(0); ++l) {
- for (unsigned short r=0; r<2; ++r) {
- node_id_set.insert(face_nodes(face_ids[l],r));
- }
+ for (unsigned short r=0; r<2; ++r) {
+ node_id_set.insert(face_nodes(face_ids[l],r));
+ }
}
Kokkos::View<unsigned int*> node_id_list("node_ids", node_id_set.size());
{
- int r=0;
- for (auto node_id : node_id_set) {
- node_id_list[r] = node_id;
- ++r;
- }
+ int r=0;
+ for (auto node_id : node_id_set) {
+ node_id_list[r] = node_id;
+ ++r;
+ }
}
connectivity.addNodeBoundary(ref, node_id_list);
}
+ std::map<unsigned int, std::vector<unsigned int>> ref_points_map;
+ for (unsigned int r=0; r<__points.extent(0); ++r) {
+ const unsigned int point_number = __points[r];
+ const unsigned int& ref = __points_ref[r];
+ ref_points_map[ref].push_back(point_number);
+ }
+
+ for (const auto& ref_point_list : ref_points_map) {
+ Kokkos::View<unsigned int*> point_list("point_list", ref_point_list.second.size());
+ for (size_t j=0; j<ref_point_list.second.size(); ++j) {
+ point_list[j]=ref_point_list.second[j];
+ }
+ const PhysicalRefId& physical_ref_id = m_physical_ref_map.at(ref_point_list.first);
+ connectivity.entityListManager().addNodeIdArray(physical_ref_id.refId(), point_list);
+ }
+
typedef Mesh<Connectivity2D> MeshType;
typedef TinyVector<2, double> Rd;
@@ -913,9 +949,9 @@ GmshReader::__nextKeyword()
}
throw ErrorHandler(__FILE__,__LINE__,
- "reading file '"+m_filename
- +"': unknown keyword '"+aKeyword+"'",
- ErrorHandler::normal);
+ "reading file '"+m_filename
+ +"': unknown keyword '"+aKeyword+"'",
+ ErrorHandler::normal);
kw.first = aKeyword;
kw.second = Unknown;
@@ -930,47 +966,47 @@ __readGmshFormat2_2()
while (kw.second != EndOfFile) {
kw = this->__nextKeyword();
switch (kw.second) {
- case NODES: {
- this->__readVertices();
- if (this->__nextKeyword().second != ENDNODES) {
- throw ErrorHandler(__FILE__,__LINE__,
- "reading file '"+m_filename
- +"': expecting $EndNodes, '"+kw.first+"' was found",
- ErrorHandler::normal);
+ case NODES: {
+ this->__readVertices();
+ if (this->__nextKeyword().second != ENDNODES) {
+ throw ErrorHandler(__FILE__,__LINE__,
+ "reading file '"+m_filename
+ +"': expecting $EndNodes, '"+kw.first+"' was found",
+ ErrorHandler::normal);
+ }
+ break;
}
- break;
- }
- case ELEMENTS: {
- this->__readElements2_2();
- kw = this->__nextKeyword();
- if (kw.second != ENDELEMENTS) {
- throw ErrorHandler(__FILE__,__LINE__,
- "reading file '"+m_filename
- +"': expecting $EndElements, '"+kw.first+"' was found",
- ErrorHandler::normal);
+ case ELEMENTS: {
+ this->__readElements2_2();
+ kw = this->__nextKeyword();
+ if (kw.second != ENDELEMENTS) {
+ throw ErrorHandler(__FILE__,__LINE__,
+ "reading file '"+m_filename
+ +"': expecting $EndElements, '"+kw.first+"' was found",
+ ErrorHandler::normal);
+ }
+ break;
}
- break;
- }
- case PHYSICALNAMES: {
- this->__readPhysicalNames2_2();
- if (this->__nextKeyword().second != ENDPHYSICALNAMES) {
- throw ErrorHandler(__FILE__,__LINE__,
- "reading file '"+m_filename
- +"': expecting $EndNodes, '"+kw.first+"' was found",
- ErrorHandler::normal);
+ case PHYSICALNAMES: {
+ this->__readPhysicalNames2_2();
+ if (this->__nextKeyword().second != ENDPHYSICALNAMES) {
+ throw ErrorHandler(__FILE__,__LINE__,
+ "reading file '"+m_filename
+ +"': expecting $EndNodes, '"+kw.first+"' was found",
+ ErrorHandler::normal);
+ }
+ break;
}
- break;
- }
- case EndOfFile: {
- break;
- }
- default: {
- throw ErrorHandler(__FILE__,__LINE__,
- "reading file '"+m_filename
- +"': unexpected '"+kw.first+"'",
- ErrorHandler::normal);
- }
+ case EndOfFile: {
+ break;
+ }
+ default: {
+ throw ErrorHandler(__FILE__,__LINE__,
+ "reading file '"+m_filename
+ +"': unexpected '"+kw.first+"'",
+ ErrorHandler::normal);
+ }
}
}
}
diff --git a/src/mesh/GmshReader.hpp b/src/mesh/GmshReader.hpp
index fb41a5c35..f02cfa16f 100644
--- a/src/mesh/GmshReader.hpp
+++ b/src/mesh/GmshReader.hpp
@@ -1,42 +1,79 @@
#ifndef GMSH_READER_HPP
#define GMSH_READER_HPP
-#include <string>
-#include <array>
-#include <vector>
-#include <map>
-#include <fstream>
#include <Kokkos_Core.hpp>
#include <TinyVector.hpp>
#include <Connectivity2D.hpp>
#include <Mesh.hpp>
+#include <RefId.hpp>
+
+#include <string>
+#include <array>
+#include <vector>
+#include <map>
+#include <fstream>
+
class GmshReader
{
public:
typedef TinyVector<3, double> R3;
-private:
- std::ifstream m_fin;
- const std::string m_filename;
-
- int _getInteger()
- {
- int i;
- m_fin >> i;
- return std::move(i);
- }
-
- double _getReal()
+ class PhysicalRefId
{
- double d;
- m_fin >> d;
+ private:
+ int m_dimension;
+ RefId m_ref_id;
+
+ public:
+ friend std::ostream& operator<<(std::ostream& os,
+ const PhysicalRefId& physical_ref_id)
+ {
+ os << physical_ref_id.m_ref_id << "[dimension " << physical_ref_id.m_dimension << "]";
+ return os;
+ }
+
+ bool operator<(const PhysicalRefId& physical_ref_id) const
+ {
+ return m_ref_id.tagNumber()<physical_ref_id.m_ref_id.tagNumber();
+ }
+
+ bool operator==(const PhysicalRefId& physical_ref_id) const
+ {
+ return m_ref_id.tagNumber()==physical_ref_id.m_ref_id.tagNumber();
+ }
+
+ const int& dimension() const
+ {
+ return m_dimension;
+ }
+
+ const RefId& refId() const
+ {
+ return m_ref_id;
+ }
+
+ PhysicalRefId& operator=(const PhysicalRefId&) = default;
+ PhysicalRefId& operator=(PhysicalRefId&&) =default;
+ PhysicalRefId(const int& dimension,
+ const RefId& ref_id)
+ : m_dimension(dimension),
+ m_ref_id(ref_id)
+ {
+ ;
+ }
+ PhysicalRefId() = default;
+ PhysicalRefId(const PhysicalRefId&) = default;
+ PhysicalRefId(PhysicalRefId&&) = default;
+ ~PhysicalRefId() = default;
+ };
- return std::move(d);
- }
+private:
+ std::ifstream m_fin;
+ const std::string m_filename;
- /**
+ /**
* Gmsh format provides a numbered, none ordrered and none dense
* vertices list, this stores the number of read vertices.
*/
@@ -44,6 +81,9 @@ private:
Kokkos::View<R3*> __vertices;
+ typedef unsigned int Point;
+ Kokkos::View<Point*> __points;
+ std::vector<int> __points_ref;
typedef TinyVector<2,unsigned int> Edge;
Kokkos::View<Edge*> __edges;
std::vector<int> __edges_ref;
@@ -60,63 +100,81 @@ private:
Kokkos::View<Hexahedron*> __hexahedra;
std::vector<int> __hexahedra_ref;
- /**
+ /**
* Gmsh format provides a numbered, none ordrered and none dense
* vertices list, this provides vertices renumbering correspondance
*/
std::vector<int> __verticesCorrepondance;
- /**
+ /**
* elements types
*/
std::vector<short> __elementType;
- /**
+ /**
* References
*/
std::vector<int> __references;
- /**
+ /**
* References
*/
std::vector<std::vector<int> > __elementVertices;
- /**
+ /**
* Stores the number of nodes associated to each primitive
- *
+ *
*/
std::vector<int> __numberOfPrimitiveNodes;
- /**
+ /**
* Array of boolean describing gmsh supported primitives by ff3d
- *
+ *
*/
std::array<bool, 15> __supportedPrimitives;
- /**
+ /**
* Primitives names according to there number
- *
+ *
*/
std::map<int, std::string> __primitivesNames;
bool __binary; /**< true if binary format */
bool __convertEndian; /**< true if needs to adapt endianess */
- /**
- * Adapts gmsh data to ff3d's structures
- *
+ std::map<unsigned int, PhysicalRefId> m_physical_ref_map;
+
+ int _getInteger()
+ {
+ int i;
+ m_fin >> i;
+
+ return std::move(i);
+ }
+
+ double _getReal()
+ {
+ double d;
+ m_fin >> d;
+
+ return std::move(d);
+ }
+
+ /**
+ * Adapts gmsh data to ff3d's structures
+ *
*/
void __proceedData();
- /**
+ /**
* Reads data in format 2.2
- *
+ *
*/
void __readGmshFormat2_2();
/**
* List of allowed keyword in mesh file
- *
+ *
*/
enum KeywordType {
// gmsh 2.2
@@ -135,7 +193,7 @@ private:
/**
* List of known keywords
- *
+ *
*/
typedef std::map<std::string, int> KeywordList;
@@ -143,44 +201,44 @@ private:
/**
* Type for keyword
- *
+ *
*/
typedef std::pair<std::string, int> Keyword;
- /**
+ /**
* Skips next comments if exists
- *
+ *
*/
void __skipComments();
- /**
+ /**
* Reads the next keyword and returns its token
- *
+ *
* @return KeywordToken
*/
Keyword __nextKeyword();
- /**
+ /**
* get list of vertices
- *
+ *
*/
void __getVertices();
/**
* Read list of vertices
- *
+ *
*/
void __readVertices();
/**
* Read all elements in format 2.0
- *
+ *
*/
void __readElements2_2();
/**
* Reads physical names
- *
+ *
*/
void __readPhysicalNames2_2();
@@ -195,7 +253,7 @@ public:
GmshReader(const std::string& filename);
~GmshReader() = default;
-
+
};
#endif // GMSH_READER_HPP
diff --git a/src/mesh/RefId.hpp b/src/mesh/RefId.hpp
new file mode 100644
index 000000000..775b1a713
--- /dev/null
+++ b/src/mesh/RefId.hpp
@@ -0,0 +1,69 @@
+#ifndef REF_ID_HPP
+#define REF_ID_HPP
+
+#include <string>
+#include <iostream>
+
+class RefId
+{
+ private:
+ unsigned int m_tag_number;
+ std::string m_tag_name;
+
+ public:
+ friend std::ostream& operator<<(std::ostream& os, const RefId& ref_id)
+ {
+ if (ref_id.m_tag_name.size()>0) {
+ os << ref_id.m_tag_name << '(' << ref_id.m_tag_number << ')';
+ } else {
+ os << ref_id.m_tag_number;
+ }
+ return os;
+ }
+
+ bool operator==(const RefId& ref_id) const
+ {
+ return ((m_tag_number==ref_id.m_tag_number) and
+ (m_tag_name==ref_id.m_tag_name));
+ }
+
+ bool operator<(const RefId& ref_id) const
+ {
+ return ((m_tag_number<ref_id.m_tag_number) or
+ ((m_tag_number==ref_id.m_tag_number) and
+ (m_tag_name<ref_id.m_tag_name)));
+ }
+
+ const unsigned int& tagNumber() const
+ {
+ return m_tag_number;
+ }
+
+ const std::string& tagName() const
+ {
+ return m_tag_name;
+ }
+
+ RefId& operator=(const RefId&) = default;
+ RefId& operator=(RefId&&) = default;
+ RefId() = default;
+ RefId(const RefId&) = default;
+ RefId(RefId&&) = default;
+ explicit RefId(const unsigned int& tag_number,
+ const std::string& tag_name)
+ : m_tag_number(tag_number),
+ m_tag_name(tag_name)
+ {
+ ;
+ }
+
+ explicit RefId(const unsigned int& tag_number)
+ : m_tag_number(tag_number)
+ {
+ ;
+ }
+
+ ~RefId() = default;
+};
+
+#endif // REF_ID_HPP
--
GitLab