Branch data Line data Source code
1 : : // *****************************************************************************
2 : : /*!
3 : : \file src/IO/GmshMeshReader.cpp
4 : : \copyright 2012-2015 J. Bakosi,
5 : : 2016-2018 Los Alamos National Security, LLC.,
6 : : 2019-2021 Triad National Security, LLC.
7 : : All rights reserved. See the LICENSE file for details.
8 : : \brief Gmsh mesh reader class definition
9 : : \details Gmsh mesh reader class definition. Currently, this class supports
10 : : line, triangle, tetrahedron, and point Gmsh element types.
11 : : */
12 : : // *****************************************************************************
13 : :
14 : : #include <limits>
15 : : #include <array>
16 : : #include <cmath>
17 : : #include <cstddef>
18 : : #include <istream>
19 : : #include <string>
20 : : #include <utility>
21 : : #include <vector>
22 : :
23 : : #include "UnsMesh.hpp"
24 : : #include "GmshMeshReader.hpp"
25 : : #include "GmshMeshIO.hpp"
26 : : #include "Reorder.hpp"
27 : : #include "PrintUtil.hpp"
28 : :
29 : : using tk::GmshMeshReader;
30 : :
31 : : void
32 : 2 : GmshMeshReader::readMesh( UnsMesh& mesh )
33 : : // *****************************************************************************
34 : : // Public interface for read a Gmsh mesh from file
35 : : //! \param[in] mesh Unstructured mesh object
36 : : // *****************************************************************************
37 : : {
38 : : // Read in mandatory "$MeshFormat" section
39 : 2 : readMeshFormat();
40 : :
41 : : // Keep reading in sections until end of file. These sections can be in
42 : : // arbitrary order, hence a while loop.
43 [ + + ]: 8 : while ( !m_inFile.eof() ) {
44 : : std::string s;
45 [ + - ]: 6 : getline( m_inFile, s );
46 [ + + ]: 6 : if ( s == "$Nodes" )
47 [ + - ]: 2 : readNodes( mesh );
48 [ + + ]: 4 : else if ( s == "$Elements" )
49 [ + - ]: 2 : readElements( mesh );
50 [ - + ]: 2 : else if ( s == "$PhysicalNames" )
51 : 0 : readPhysicalNames();
52 : : }
53 : 2 : }
54 : :
55 : : void
56 [ + - ]: 2 : GmshMeshReader::readMeshFormat()
57 : : // *****************************************************************************
58 : : // Read mandatory "$MeshFormat--$EndMeshFormat" section
59 : : // *****************************************************************************
60 : : {
61 : : using tk::operator<<;
62 : : std::string s;
63 : :
64 : : // Read in beginning of header: $MeshFormat
65 [ + - ]: 2 : getline( m_inFile, s );
66 [ - + ][ - - ]: 2 : ErrChk( s == "$MeshFormat",
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ][ - - ]
67 : : std::string("Unsupported mesh format '") + s + "' in file " +
68 : : m_filename );
69 : :
70 : : // Read in "version-number file-type data-size"
71 : : int type;
72 [ + - ][ + - ]: 2 : m_inFile >> m_version >> type >> m_datasize;
[ + - ]
73 [ + + ]: 2 : if (type == 0 )
74 : 1 : m_type = GmshFileType::ASCII;
75 [ + - ]: 1 : else if (type == 1 )
76 : 1 : m_type = GmshFileType::BINARY;
77 : :
78 [ - + ][ - - ]: 2 : ErrChk( ( fabs(m_version-2.2) < std::numeric_limits<tk::real>::epsilon() ||
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ]
79 : : fabs(m_version-2.0) < std::numeric_limits<tk::real>::epsilon() ),
80 : : std::string("Unsupported mesh version '") << m_version <<
81 : : "' in file " << m_filename );
82 : :
83 [ - + ][ - - ]: 2 : ErrChk( ( m_type == GmshFileType::ASCII || m_type == GmshFileType::BINARY ),
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ]
84 : : std::string("Unsupported mesh type '") << m_type << "' in file " <<
85 : : m_filename );
86 : :
87 [ - + ][ - - ]: 2 : ErrChk( m_datasize == sizeof(tk::real),
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ][ - - ]
88 : : std::string("Unsupported mesh datasize '") << m_datasize <<
89 : : "' in file " << m_filename );
90 : :
91 [ + - ]: 2 : getline( m_inFile, s ); // finish reading the line
92 : :
93 : : // if file is binary, binary-read in binary "one"
94 [ + + ]: 2 : if ( isBinary() ) {
95 : : int one;
96 [ + - ]: 1 : m_inFile.read( reinterpret_cast<char*>(&one), sizeof(int) );
97 : : #ifdef __bg__
98 : : one = tk::swap_endian< int >( one );
99 : : #endif
100 [ - + ][ - - ]: 1 : ErrChk( one == 1, "Endianness does not match in file " + m_filename );
[ - - ][ - - ]
[ - - ][ - - ]
101 [ + - ]: 1 : getline( m_inFile, s ); // finish reading the line
102 : : }
103 : :
104 : : // Read in end of header: $EndMeshFormat
105 [ + - ]: 2 : getline( m_inFile, s );
106 [ - + ][ - - ]: 2 : ErrChk( s == "$EndMeshFormat",
[ - - ][ - - ]
[ - - ][ - - ]
107 : : "'$EndMeshFormat' keyword is missing in file " + m_filename );
108 : 2 : }
109 : :
110 : : void
111 : 2 : GmshMeshReader::readNodes( UnsMesh& mesh )
112 : : // *****************************************************************************
113 : : // Read "$Nodes--$EndNodes" section
114 : : //! \param[in] mesh Unstructured mesh object
115 : : // *****************************************************************************
116 : : {
117 : : // Read in number of nodes in this node set
118 : : std::size_t nnode;
119 : 2 : m_inFile >> nnode;
120 [ - + ][ - - ]: 2 : ErrChk( nnode > 0,
[ - - ][ - - ]
[ - - ][ - - ]
121 : : "Number of nodes must be greater than zero in file " + m_filename );
122 : : std::string s;
123 [ + + ][ + - ]: 2 : if (isBinary()) getline( m_inFile, s ); // finish reading the line
124 : :
125 : : // Read in node ids and coordinates: node-number x-coord y-coord z-coord
126 [ + + ]: 30 : for ( std::size_t i=0; i<nnode; ++i ) {
127 : : int id;
128 : : std::array< tk::real, 3 > coord;
129 : :
130 [ + + ]: 28 : if (isASCII()) {
131 [ + - ]: 14 : m_inFile >> id >> coord[0] >> coord[1] >> coord[2];
132 : : } else {
133 [ + - ]: 14 : m_inFile.read( reinterpret_cast<char*>(&id), sizeof(int) );
134 : : #ifdef __bg__
135 : : id = tk::swap_endian< int >( id );
136 : : #endif
137 [ + - ]: 14 : m_inFile.read( reinterpret_cast<char*>(coord.data()), 3*sizeof(double) );
138 : : #ifdef __bg__
139 : : coord[0] = tk::swap_endian< double >( coord[0] );
140 : : coord[1] = tk::swap_endian< double >( coord[1] );
141 : : coord[2] = tk::swap_endian< double >( coord[2] );
142 : : #endif
143 : : }
144 : :
145 [ + - ]: 28 : mesh.x().push_back( coord[0] );
146 [ + - ]: 28 : mesh.y().push_back( coord[1] );
147 [ + - ]: 28 : mesh.z().push_back( coord[2] );
148 : : }
149 [ + - ]: 2 : getline( m_inFile, s ); // finish reading the last line
150 : :
151 : : // Read in end of header: $EndNodes
152 [ + - ]: 2 : getline( m_inFile, s );
153 [ - + ][ - - ]: 2 : ErrChk( s == "$EndNodes",
[ - - ][ - - ]
[ - - ][ - - ]
154 : : "'$EndNodes' keyword is missing in file" + m_filename );
155 : 2 : }
156 : :
157 : : void
158 [ + - ]: 2 : GmshMeshReader::readElements( UnsMesh& mesh )
159 : : // *****************************************************************************
160 : : // Read "$Elements--$EndElements" section
161 : : //! \param[in] mesh Unstructured mesh object
162 : : // *****************************************************************************
163 : : {
164 : : using tk::operator<<;
165 : :
166 : : std::string s;
167 : :
168 : : // Read in number of elements in this element set
169 : : int nel;
170 [ + - ]: 2 : m_inFile >> nel;
171 [ - + ][ - - ]: 2 : ErrChk( nel > 0, "Number of elements must be greater than zero in file " +
[ - - ][ - - ]
[ - - ][ - - ]
172 : : m_filename );
173 [ + - ]: 2 : getline( m_inFile, s ); // finish reading the last line
174 : :
175 : : // Read in element ids, tags, and element connectivity (node list)
176 : 2 : int n=1;
177 [ + + ]: 27 : for (int i=0; i<nel; i+=n) {
178 : : int id, elmtype, ntags;
179 : :
180 [ + + ]: 25 : if (isASCII()) {
181 : : // elm-number elm-type number-of-tags < tag > ... node-number-list
182 [ + - ][ + - ]: 24 : m_inFile >> id >> elmtype >> ntags;
[ + - ]
183 : : } else {
184 : : // elm-type num-of-elm-follow number-of-tags
185 [ + - ]: 1 : m_inFile.read( reinterpret_cast<char*>(&elmtype), sizeof(int) );
186 [ + - ]: 1 : m_inFile.read( reinterpret_cast<char*>(&n), sizeof(int) );
187 [ + - ]: 1 : m_inFile.read( reinterpret_cast<char*>(&ntags), sizeof(int) );
188 : : #ifdef __bg__
189 : : elmtype = tk::swap_endian< int >( elmtype );
190 : : n = tk::swap_endian< int >( n );
191 : : ntags = tk::swap_endian< int >( ntags );
192 : : #endif
193 : : }
194 : :
195 : : // Find element type, throw exception if not supported
196 : : const auto it = m_elemNodes.find( elmtype );
197 [ + - ][ - - ]: 25 : ErrChk( it != m_elemNodes.end(),
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ][ - - ]
198 : : std::string("Unsupported element type ") << elmtype <<
199 : : " in mesh file: " << m_filename );
200 : :
201 [ + + ]: 73 : for (int e=0; e<n; ++e) {
202 : : // Read element id if binary
203 [ + + ]: 48 : if (isBinary()) {
204 [ + - ]: 24 : m_inFile.read( reinterpret_cast<char*>(&id), sizeof(int) );
205 : : #ifdef __bg__
206 : : id = tk::swap_endian< int >( id );
207 : : #endif
208 : : }
209 : :
210 : : // Read and ignore element tags
211 [ + - ][ - - ]: 48 : std::vector< int > tags( static_cast<std::size_t>(ntags), 0 );
212 [ + + ]: 48 : if (isASCII()) {
213 [ + + ]: 48 : for (std::size_t j=0; j<static_cast<std::size_t>(ntags); j++)
214 [ + - ]: 24 : m_inFile >> tags[j];
215 : : } else {
216 : : m_inFile.read(
217 : 24 : reinterpret_cast<char*>(tags.data()),
218 : : static_cast<std::streamsize>(
219 [ + - ]: 24 : static_cast<std::size_t>(ntags) * sizeof(int) ) );
220 : : #ifdef __bg__
221 : : for (auto& t : tags) t = tk::swap_endian< int >( t );
222 : : #endif
223 : : }
224 : :
225 : : // Read and add element node list (i.e. connectivity)
226 : 48 : std::size_t nnode = static_cast< std::size_t >( it->second );
227 [ + - ][ - - ]: 48 : std::vector< std::size_t > nodes( nnode, 0 );
228 [ + + ]: 48 : if (isASCII()) {
229 [ + + ]: 120 : for (std::size_t j=0; j<nnode; j++)
230 [ + - ]: 96 : m_inFile >> nodes[j];
231 : : } else {
232 [ + - ][ - - ]: 24 : std::vector< int > nds( nnode, 0 );
233 : : m_inFile.read(
234 : 24 : reinterpret_cast< char* >( nds.data() ),
235 [ + - ]: 24 : static_cast< std::streamsize >( nnode * sizeof(int) ) );
236 : : #ifdef __bg__
237 : : for (auto& j : nds) j = tk::swap_endian< int >( j );
238 : : #endif
239 [ + + ]: 120 : for (std::size_t j=0; j<nnode; j++)
240 : 96 : nodes[j] = static_cast< std::size_t >( nds[j] );
241 : : }
242 : : // Put in element connectivity for different types of elements
243 [ - - ][ + - ]: 48 : switch ( elmtype ) {
[ - ]
244 : 0 : case GmshElemType::LIN:
245 [ - - ][ - - ]: 0 : for (const auto& j : nodes) mesh.lininpoel().push_back( j );
246 : : break;
247 : 0 : case GmshElemType::TRI:
248 [ - - ][ - - ]: 0 : for (const auto& j : nodes) mesh.triinpoel().push_back( j );
249 : : break;
250 : 48 : case GmshElemType::TET:
251 [ + + ][ + - ]: 240 : for (const auto& j : nodes) mesh.tetinpoel().push_back( j );
252 : : break;
253 : : case GmshElemType::PNT:
254 : : break; // ignore 1-node 'point element' type
255 [ - - ][ - - ]: 0 : default: Throw( std::string("Unsupported element type ") << elmtype <<
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ]
256 : : " in mesh file: " << m_filename );
257 : : }
258 : : }
259 : : }
260 [ + - ]: 2 : getline( m_inFile, s ); // finish reading the last line
261 : :
262 : : // Shift node IDs to start from zero (gmsh likes one-based node ids)
263 [ + - ]: 2 : shiftToZero( mesh.lininpoel() );
264 [ + - ]: 2 : shiftToZero( mesh.triinpoel() );
265 [ + - ]: 2 : shiftToZero( mesh.tetinpoel() );
266 : :
267 : : // Read in end of header: $EndNodes
268 [ + - ]: 2 : getline( m_inFile, s );
269 [ - + ][ - - ]: 2 : ErrChk( s == "$EndElements",
[ - - ][ - - ]
[ - - ][ - - ]
270 : : "'$EndElements' keyword is missing in file" + m_filename );
271 : 2 : }
272 : :
273 : : void
274 : 0 : GmshMeshReader::readPhysicalNames()
275 : : // *****************************************************************************
276 : : // Read "$PhysicalNames--$EndPhysicalNames" section
277 : : // *****************************************************************************
278 : : {
279 [ - - ][ - - ]: 0 : Throw( "Mesh section '$PhysicalNames -- $EndPhysicalNames' not implemented" );
[ - - ][ - - ]
[ - - ][ - - ]
280 : : }
|
