Branch data Line data Source code
1 : : // *****************************************************************************
2 : : /*!
3 : : \file src/IO/MeshFactory.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 Unstructured mesh reader and writer factory
9 : : \details Unstructured mesh reader and writer factory.
10 : : */
11 : : // *****************************************************************************
12 : :
13 : : #include <string>
14 : :
15 : : #include "MeshFactory.hpp"
16 : : #include "MeshDetect.hpp"
17 : : #include "Timer.hpp"
18 : : #include "Reader.hpp"
19 : : #include "GmshMeshReader.hpp"
20 : : #include "NetgenMeshReader.hpp"
21 : : #include "ExodusIIMeshReader.hpp"
22 : : #include "HyperMeshReader.hpp"
23 : : #include "UGRIDMeshReader.hpp"
24 : : #include "RDGFLOMeshReader.hpp"
25 : : #include "ASCMeshReader.hpp"
26 : : #include "NetgenMeshWriter.hpp"
27 : : #include "GmshMeshWriter.hpp"
28 : : #include "ExodusIIMeshWriter.hpp"
29 : : #include "DerivedData.hpp"
30 : : #include "Reorder.hpp"
31 : : #include "QuinoaBuildConfig.hpp"
32 : :
33 : : #ifdef HAS_OMEGA_H
34 : : #include "Omega_h_MeshReader.hpp"
35 : : #endif
36 : :
37 : : namespace tk {
38 : :
39 : : UnsMesh
40 : 47 : readUnsMesh( const tk::Print& print,
41 : : const std::string& filename,
42 : : std::pair< std::string, tk::real >& timestamp )
43 : : // *****************************************************************************
44 : : // Read unstructured mesh from file
45 : : //! \param[in] print Pretty printer
46 : : //! \param[in] filename Filename to read mesh from
47 : : //! \param[out] timestamp A time stamp consisting of a timer label (a string),
48 : : //! and a time state (a tk::real in seconds) measuring the mesh read time
49 : : //! \return Unstructured mesh object
50 : : // *****************************************************************************
51 : : {
52 [ + - ][ + - ]: 94 : print.diagstart( "Reading mesh from file '" + filename + "' ..." );
[ - + ][ - - ]
53 : :
54 : : // Read in mesh
55 : : tk::Timer t;
56 : :
57 : : //! Create unstructured mesh to store mesh
58 : 47 : UnsMesh mesh;
59 : :
60 [ + - ]: 47 : const auto meshtype = detectInput( filename );
61 : :
62 [ + + ]: 47 : if (meshtype == MeshReaderType::GMSH)
63 [ + - ][ + - ]: 10 : GmshMeshReader( filename ).readMesh( mesh );
64 [ + + ]: 42 : else if (meshtype == MeshReaderType::NETGEN)
65 [ + - ]: 4 : NetgenMeshReader( filename ).readMesh( mesh );
66 [ + + ]: 40 : else if (meshtype == MeshReaderType::EXODUSII)
67 [ + - ][ + - ]: 33 : ExodusIIMeshReader( filename ).readMesh( mesh );
68 [ + + ]: 7 : else if (meshtype == MeshReaderType::ASC)
69 [ + - ]: 8 : ASCMeshReader( filename ).readMesh( mesh );
70 [ + + ]: 3 : else if (meshtype == MeshReaderType::UGRID)
71 [ + - ]: 2 : UGRIDMeshReader( filename ).readMesh( mesh );
72 [ + + ]: 2 : else if (meshtype == MeshReaderType::RDGFLO)
73 [ + - ]: 2 : RDGFLOMeshReader( filename ).readMesh( mesh );
74 [ + - ]: 1 : else if (meshtype == MeshReaderType::HYPER)
75 [ + - ]: 2 : HyperMeshReader( filename ).readMesh( mesh );
76 : :
77 : : timestamp =
78 [ + - ][ + - ]: 94 : std::make_pair( "Read mesh from file '" + filename + '\'', t.dsec() );
[ - + ][ - + ]
79 : :
80 [ + - ][ + - ]: 47 : print.diagend( "done" );
81 : :
82 : : // Return (move out) mesh object
83 : 47 : return mesh;
84 : : }
85 : :
86 : : std::vector< std::pair< std::string, tk::real > >
87 : 18 : writeUnsMesh( const tk::Print& print,
88 : : const std::string& filename,
89 : : UnsMesh& mesh,
90 : : bool reorder )
91 : : // *****************************************************************************
92 : : // Write unstructured mesh to file
93 : : //! \param[in] print Pretty printer
94 : : //! \param[in] filename Filename to write mesh to
95 : : //! \param[in] mesh Unstructured mesh object to write from
96 : : //! \param[in] reorder Whether to also reorder mesh nodes
97 : : //! \return Vector of time stamps consisting of a timer label (a string), and a
98 : : //! time state (a tk::real in seconds) measuring the renumber and the mesh
99 : : //! write time
100 : : // *****************************************************************************
101 : : {
102 : : std::vector< std::pair< std::string, tk::real > > times;
103 : :
104 : : tk::Timer t;
105 : :
106 : : // If mesh has tetrahedra but no triangles, generate triangle connectivity
107 [ + + ][ + + ]: 18 : if (!mesh.tetinpoel().empty() && mesh.triinpoel().empty()) {
108 [ + - ][ + - ]: 12 : print.diagstart( "Generating missing surface mesh ..." );
[ + - ]
109 : :
110 : : const auto& inpoel = mesh.tetinpoel(); // get tet connectivity
111 [ + - ]: 12 : auto esup = tk::genEsup( inpoel, 4 ); // elements surrounding points
112 [ + - ]: 6 : auto esuel = tk::genEsuelTet( inpoel, esup ); // elems surrounding elements
113 : : auto& triinpoel = mesh.triinpoel();
114 : : // collect boundary faces
115 [ + + ]: 426426 : for (std::size_t e=0; e<esuel.size()/4; ++e) {
116 : 426420 : auto mark = e*4;
117 [ + + ]: 2132100 : for (std::size_t f=0; f<4; ++f) {
118 [ + + ]: 1705680 : if (esuel[mark+f] == -1) {
119 : : // extract triangle element connectivity from tetrahedron
120 [ + - ]: 46712 : triinpoel.push_back( inpoel[ mark+tk::lpofa[f][0] ] );
121 [ + - ]: 46712 : triinpoel.push_back( inpoel[ mark+tk::lpofa[f][1] ] );
122 [ + - ]: 46712 : triinpoel.push_back( inpoel[ mark+tk::lpofa[f][2] ] );
123 : : }
124 : : }
125 : : }
126 : :
127 [ + - ][ + - ]: 12 : print.diagend( "done" );
128 [ + - ][ - - ]: 6 : times.emplace_back( "Generate surface mesh", t.dsec() );
129 : : t.zero();
130 : : }
131 : :
132 [ + + ]: 18 : if (reorder) {
133 [ + - ][ + - ]: 2 : print.diagstart( "Reordering mesh nodes ..." );
[ + - ]
134 : :
135 : : // If mesh has tetrahedra elements, reorder based on those
136 [ + - ]: 1 : if (!mesh.tetinpoel().empty()) {
137 : :
138 : : auto& inpoel = mesh.tetinpoel();
139 [ + - ][ + - ]: 2 : const auto psup = tk::genPsup( inpoel, 4, tk::genEsup( inpoel, 4 ) );
140 [ + - ]: 1 : auto map = tk::renumber( psup );
141 [ + - ]: 1 : tk::remap( inpoel, map );
142 [ + - ]: 1 : tk::remap( mesh.triinpoel(), map );
143 [ + - ]: 1 : tk::remap( mesh.x(), map );
144 [ + - ]: 1 : tk::remap( mesh.y(), map );
145 [ + - ]: 1 : tk::remap( mesh.z(), map );
146 : :
147 : : // If mesh has no tetrahedra elements, reorder based on triangle mesh if any
148 [ - - ]: 0 : } else if (!mesh.triinpoel().empty()) {
149 : :
150 : : auto& inpoel = mesh.triinpoel();
151 [ - - ][ - - ]: 0 : const auto psup = tk::genPsup( inpoel, 3, tk::genEsup( inpoel, 3 ) );
152 [ - - ]: 0 : auto map = tk::renumber( psup );
153 [ - - ]: 0 : tk::remap( inpoel, map );
154 [ - - ]: 0 : tk::remap( mesh.x(), map );
155 [ - - ]: 0 : tk::remap( mesh.y(), map );
156 [ - - ]: 0 : tk::remap( mesh.z(), map );
157 : : }
158 : :
159 [ + - ][ + - ]: 2 : print.diagend( "done" );
160 [ + - ]: 1 : times.emplace_back( "Reorder mesh", t.dsec() );
161 : : t.zero();
162 : : }
163 : :
164 [ + - ][ + - ]: 36 : print.diagstart( "Writing mesh to file '" + filename + "' ..." );
[ + - ][ - + ]
[ - - ]
165 : :
166 [ + - ]: 18 : const auto meshtype = pickOutput( filename );
167 : :
168 [ + + ]: 18 : if (meshtype == MeshWriterType::GMSH)
169 [ + - ][ + - ]: 8 : GmshMeshWriter( filename ).writeMesh( mesh );
170 [ + + ]: 14 : else if (meshtype == MeshWriterType::NETGEN)
171 [ + - ]: 8 : NetgenMeshWriter( filename ).writeMesh( mesh );
172 [ + - ]: 10 : else if (meshtype== MeshWriterType::EXODUSII)
173 [ + - ][ + - ]: 10 : ExodusIIMeshWriter( filename, ExoWriter::CREATE ).writeMesh( mesh );
174 : :
175 [ + - ][ + - ]: 36 : print.diagend( "done" );
176 [ + - ]: 18 : times.emplace_back( "Write mesh to file", t.dsec() );
177 : :
178 : 18 : return times;
179 : : }
180 : :
181 : : } // tk::
|