Branch data Line data Source code
1 : : // *****************************************************************************
2 : : /*!
3 : : \file src/Inciter/LaxCG.hpp
4 : : \copyright 2012-2015 J. Bakosi,
5 : : 2016-2018 Los Alamos National Security, LLC.,
6 : : 2019-2021 Triad National Security, LLC.,
7 : : 2022-2024 J. Bakosi
8 : : All rights reserved. See the LICENSE file for details.
9 : : \brief LaxCG: Time-derivative preconditinoing for all Ma
10 : : \see Luo, Baum, Lohner, "Extension of Harten-Lax-van Leer Scheme for
11 : : Flows at All Speeds", AIAA Journal, Vol. 43, No. 6, 2005
12 : : \see Weiss & Smith, "Preconditioning Applied to Variable and Constant
13 : : Density Time-Accurate Flows on Unstructured Meshes", AIAA Journal,
14 : : Vol. 33, No. 11, 1995, pp. 2050-2057.
15 : : */
16 : : // *****************************************************************************
17 : :
18 : : #pragma once
19 : :
20 : : #include <vector>
21 : : #include <map>
22 : :
23 : : #include "Types.hpp"
24 : : #include "Fields.hpp"
25 : : #include "Table.hpp"
26 : : #include "DerivedData.hpp"
27 : : #include "NodeDiagnostics.hpp"
28 : :
29 : : #include "NoWarning/laxcg.decl.h"
30 : :
31 : : namespace inciter {
32 : :
33 : : //! LaxCG Charm++ chare array used to advance PDEs in time with LaxCG
34 : : class LaxCG : public CBase_LaxCG {
35 : :
36 : : public:
37 : : #if defined(__clang__)
38 : : #pragma clang diagnostic push
39 : : #pragma clang diagnostic ignored "-Wunused-parameter"
40 : : #pragma clang diagnostic ignored "-Wdeprecated-declarations"
41 : : #elif defined(STRICT_GNUC)
42 : : #pragma GCC diagnostic push
43 : : #pragma GCC diagnostic ignored "-Wunused-parameter"
44 : : #pragma GCC diagnostic ignored "-Wdeprecated-declarations"
45 : : #elif defined(__INTEL_COMPILER)
46 : : #pragma warning( push )
47 : : #pragma warning( disable: 1478 )
48 : : #endif
49 : : // Include Charm++ SDAG code. See http://charm.cs.illinois.edu/manuals/html/
50 : : // charm++/manual.html, Sec. "Structured Control Flow: Structured Dagger".
51 : : LaxCG_SDAG_CODE
52 : : #if defined(__clang__)
53 : : #pragma clang diagnostic pop
54 : : #elif defined(STRICT_GNUC)
55 : : #pragma GCC diagnostic pop
56 : : #elif defined(__INTEL_COMPILER)
57 : : #pragma warning( pop )
58 : : #endif
59 : :
60 : : //! Constructor
61 : : explicit LaxCG( const CProxy_Discretization& disc,
62 : : const std::map< int, std::vector< std::size_t > >& bface,
63 : : const std::map< int, std::vector< std::size_t > >& bnode,
64 : : const std::vector< std::size_t >& triinpoel );
65 : :
66 : : #if defined(__clang__)
67 : : #pragma clang diagnostic push
68 : : #pragma clang diagnostic ignored "-Wundefined-func-template"
69 : : #endif
70 : : //! Migrate constructor
71 : : // cppcheck-suppress uninitMemberVar
72 : 687 : explicit LaxCG( CkMigrateMessage* m ) : CBase_LaxCG( m ) {}
73 : : #if defined(__clang__)
74 : : #pragma clang diagnostic pop
75 : : #endif
76 : :
77 : : //! Configure Charm++ custom reduction types initiated from this chare array
78 : : static void registerReducers();
79 : :
80 : : //! Return from migration
81 : : void ResumeFromSync() override;
82 : :
83 : : //! Start setup for solution
84 : : void setup( tk::real v );
85 : :
86 : : // Start time stepping
87 : : void start();
88 : :
89 : : //! Advance equations to next time step
90 : : void advance( tk::real newdt );
91 : :
92 : : //! Start (re-)computing domain and boundary integrals
93 : : void feop();
94 : :
95 : : //! Receive contributions to boundary point normals on chare-boundaries
96 : : void comnorm( const std::unordered_map< int,
97 : : std::unordered_map< std::size_t, std::array<tk::real,4> > >& inbnd );
98 : :
99 : : //! Receive contributions to node gradients on chare-boundaries
100 : : void comgrad( const std::unordered_map< std::size_t,
101 : : std::vector< tk::real > >& ingrad );
102 : :
103 : : //! Receive contributions to right-hand side vector on chare-boundaries
104 : : void comrhs( const std::unordered_map< std::size_t,
105 : : std::vector< tk::real > >& inrhs );
106 : :
107 : : //! Evaluate residuals
108 : : void evalres( const std::vector< tk::real >& l2res );
109 : :
110 : : //! Receive new mesh from Refiner
111 : : void resizePostAMR(
112 : : const std::vector< std::size_t >& ginpoel,
113 : : const tk::UnsMesh::Chunk& chunk,
114 : : const tk::UnsMesh::Coords& coord,
115 : : const std::unordered_map< std::size_t, tk::UnsMesh::Edge >& addedNodes,
116 : : const std::unordered_map< std::size_t, std::size_t >& addedTets,
117 : : const std::set< std::size_t >& removedNodes,
118 : : const std::unordered_map< int, std::unordered_set< std::size_t > >&
119 : : nodeCommMap,
120 : : const std::map< int, std::vector< std::size_t > >& bface,
121 : : const std::map< int, std::vector< std::size_t > >& bnode,
122 : : const std::vector< std::size_t >& triinpoel );
123 : :
124 : : //! Const-ref access to current solution
125 : : //! \return Const-ref to current solution
126 : 0 : const tk::Fields& solution() const { return m_u; }
127 : :
128 : : //! Compute integral quantities for output
129 : : void integrals();
130 : :
131 : : //! Evaluate whether to continue with next time step
132 : : void step();
133 : :
134 : : // Evaluate whether to do load balancing
135 : : void evalLB( int nrestart );
136 : :
137 : : //! Evaluate whether to continue with next time step stage
138 : : void stage();
139 : :
140 : : /** @name Charm++ pack/unpack serializer member functions */
141 : : ///@{
142 : : //! \brief Pack/Unpack serialize member function
143 : : //! \param[in,out] p Charm++'s PUP::er serializer object reference
144 : 2085 : void pup( PUP::er &p ) override {
145 : 2085 : p | m_disc;
146 : 2085 : p | m_nrhs;
147 : 2085 : p | m_nnorm;
148 : 2085 : p | m_nbpint;
149 : 2085 : p | m_nbeint;
150 : 2085 : p | m_ndeint;
151 : 2085 : p | m_ngrad;
152 : 2085 : p | m_bnode;
153 : 2085 : p | m_bface;
154 : 2085 : p | m_triinpoel;
155 : 2085 : p | m_u;
156 : : // do not pup these, will recompute after migration anyway
157 [ + + ]: 2085 : if (p.isUnpacking()) {
158 : 687 : m_un.resize( m_u.nunk(), m_u.nprop() );
159 : 687 : m_rhs.resize( m_u.nunk(), m_u.nprop() );
160 : 687 : m_grad.resize( m_u.nunk(), m_u.nprop()*3 );
161 : : }
162 : 2085 : p | m_rhsc;
163 : 2085 : p | m_gradc;
164 : 2085 : p | m_diag;
165 : 2085 : p | m_bnorm;
166 : 2085 : p | m_bnormc;
167 : 2085 : p | m_bndpoinint;
168 : 2085 : p | m_domedgeint;
169 : 2085 : p | m_bpint;
170 : 2085 : p | m_dsupedge;
171 : 2085 : p | m_dsupint;
172 : 2085 : p | m_besym;
173 : 2085 : p | m_dirbcmasks;
174 : 2085 : p | m_prebcnodes;
175 : 2085 : p | m_prebcvals;
176 : 2085 : p | m_symbcnodeset;
177 : 2085 : p | m_symbcnodes;
178 : 2085 : p | m_symbcnorms;
179 : 2085 : p | m_farbcnodeset;
180 : 2085 : p | m_farbcnodes;
181 : 2085 : p | m_farbcnorms;
182 : 2085 : p | m_surfint;
183 : 2085 : p | m_stage;
184 : 2085 : p | m_dtp;
185 : 2085 : p | m_tp;
186 : 2085 : p | m_finished;
187 : 2085 : }
188 : : //! \brief Pack/Unpack serialize operator|
189 : : //! \param[in,out] p Charm++'s PUP::er serializer object reference
190 : : //! \param[in,out] i LaxCG object reference
191 : : friend void operator|( PUP::er& p, LaxCG& i ) { i.pup(p); }
192 : : //@}
193 : :
194 : : private:
195 : : //! Discretization proxy
196 : : CProxy_Discretization m_disc;
197 : : //! Counter for right-hand side vector nodes updated
198 : : std::size_t m_nrhs;
199 : : //! Counter for receiving boundary point normals
200 : : std::size_t m_nnorm;
201 : : //! Counter for receiving boundary point integrals
202 : : std::size_t m_nbpint;
203 : : //! Counter for receiving boundary edge integrals
204 : : std::size_t m_nbeint;
205 : : //! Counter for receiving domain edge integrals
206 : : std::size_t m_ndeint;
207 : : //! Counter for receiving gradients
208 : : std::size_t m_ngrad;
209 : : //! Boundary node lists mapped to side set ids used in the input file
210 : : std::map< int, std::vector< std::size_t > > m_bnode;
211 : : //! Boundary face lists mapped to side set ids used in the input file
212 : : std::map< int, std::vector< std::size_t > > m_bface;
213 : : //! Boundary triangle face connecitivity where BCs are set by user
214 : : std::vector< std::size_t > m_triinpoel;
215 : : //! Unknown/solution vector at mesh nodes
216 : : tk::Fields m_u;
217 : : //! Unknown/solution vector at mesh nodes at previous time
218 : : tk::Fields m_un;
219 : : //! Right-hand side vector (for the high order system)
220 : : tk::Fields m_rhs;
221 : : //! Receive buffer for communication of the right hand side
222 : : //! \details Key: global node id, value: rhs for all scalar components per
223 : : //! node.
224 : : std::unordered_map< std::size_t, std::vector< tk::real > > m_rhsc;
225 : : //! Diagnostics object
226 : : NodeDiagnostics m_diag;
227 : : //! Boundary point normals
228 : : //! \details Outer key: side set id. Inner key: global node id of boundary
229 : : //! point, value: weighted normals, inverse distance square, nodal area.
230 : : std::unordered_map< int,
231 : : std::unordered_map< std::size_t, std::array<tk::real,4> > > m_bnorm;
232 : : //! Boundary point normals receive buffer
233 : : //! \details Outer key: side set id. Inner key: global node id of boundary
234 : : //! point, value: weighted normals and inverse distance square.
235 : : decltype(m_bnorm) m_bnormc;
236 : : //! Boundary point integrals
237 : : //! \details Key: global node id of boundary point, value: boundary point
238 : : //! integral contributions.
239 : : std::unordered_map< std::size_t, std::array<tk::real,3> > m_bndpoinint;
240 : : //! Domain edge integrals
241 : : std::unordered_map< tk::UnsMesh::Edge, std::array< tk::real, 3 >,
242 : : tk::UnsMesh::Hash<2>, tk::UnsMesh::Eq<2> > m_domedgeint;
243 : : //! Streamable boundary point integrals
244 : : std::vector< tk::real > m_bpint;
245 : : //! Superedge (tet, face, edge) end points with local ids for domain edges
246 : : std::array< std::vector< std::size_t >, 3 > m_dsupedge;
247 : : //! Superedge (tet, face, edge) domain edge integrals
248 : : std::array< std::vector< tk::real >, 3 > m_dsupint;
249 : : //! Streamable boundary element symmetry BC flags
250 : : std::vector< std::uint8_t > m_besym;
251 : : //! Gradients in mesh nodes
252 : : tk::Fields m_grad;
253 : : //! Gradients receive buffer
254 : : std::unordered_map< std::size_t, std::vector< tk::real > > m_gradc;
255 : : //! Nodes and their Dirichlet BC masks
256 : : std::vector< std::size_t > m_dirbcmasks;
257 : : //! Nodes at pressure BCs
258 : : std::vector< std::size_t > m_prebcnodes;
259 : : //! Density and pressure values at pressure BCs
260 : : std::vector< tk::real > m_prebcvals;
261 : : //! Unique set of ordered nodes at which symmetry BCs are set
262 : : std::set< std::size_t > m_symbcnodeset;
263 : : //! Streamable nodes at which symmetry BCs are set
264 : : std::vector< std::size_t > m_symbcnodes;
265 : : //! Streamable normals at nodes at which symmetry BCs are set
266 : : std::vector< tk::real > m_symbcnorms;
267 : : //! Unique set of ordered nodes at which farfield BCs are set
268 : : std::set< std::size_t > m_farbcnodeset;
269 : : //! Streamable nodes at which farfield BCs are set
270 : : std::vector< std::size_t > m_farbcnodes;
271 : : //! Streamable normals at nodes at which farfield BCs are set
272 : : std::vector< tk::real > m_farbcnorms;
273 : : //! Streamable surface integral nodes and normals * dA on surfaces
274 : : std::map< int, std::pair< std::vector< std::size_t >,
275 : : std::vector< tk::real > > > m_surfint;
276 : : //! Runge-Kutta stage counter
277 : : std::size_t m_stage;
278 : : //! Time step size for each mesh node
279 : : std::vector< tk::real > m_dtp;
280 : : //! Physical time for each mesh node
281 : : std::vector< tk::real > m_tp;
282 : : //! True in the last time step
283 : : int m_finished;
284 : :
285 : : //! Access bound Discretization class pointer
286 : 301107 : Discretization* Disc() const {
287 [ + - ][ + - ]: 301107 : Assert( m_disc[ thisIndex ].ckLocal() != nullptr, "ckLocal() null" );
[ - + ][ - - ]
[ - - ][ - - ]
288 [ + - ][ + - ]: 301107 : return m_disc[ thisIndex ].ckLocal();
289 : : }
290 : :
291 : : //! Convert from conservative to primitive variables
292 : : void primitive( tk::Fields& U );
293 : :
294 : : //! Convert from primitive to conservative variables
295 : : void conservative( tk::Fields& U );
296 : :
297 : : //! Compute the inverse time-derivative preconditioning matrix
298 : : std::array< tk::real, 5*5 > precond( const tk::Fields& U, std::size_t i );
299 : :
300 : : // Compute characteristic velocity of the preconditioned system at a point
301 : : tk::real charvel( std::size_t i );
302 : :
303 : : //! Prepare boundary condition data structures
304 : : void setupBC();
305 : :
306 : : //! Compute local contributions to domain edge integrals
307 : : void domint();
308 : :
309 : : //! Compute chare-boundary edges
310 : : void bndEdges();
311 : :
312 : : //! Compute boundary point normals
313 : : void bndint();
314 : :
315 : : //! Combine own and communicated portions of the boundary point normals
316 : : void bnorm();
317 : :
318 : : //! Convert integrals into streamable data structures
319 : : void streamable();
320 : :
321 : : //! Generate superedge-groups for domain-edge loops
322 : : void domsuped();
323 : :
324 : : //! Output mesh and particle fields to files
325 : : void out();
326 : :
327 : : //! Output mesh-based fields to file
328 : : void writeFields( CkCallback cb );
329 : :
330 : : //! Combine own and communicated portions of the integrals
331 : : void merge();
332 : :
333 : : //! Compute righ-hand side vector of transport equations
334 : : void rhs();
335 : :
336 : : //! Advance systems of equations
337 : : void solve();
338 : :
339 : : //! Optionally refine/derefine mesh
340 : : void refine();
341 : :
342 : : //! Compute time step size
343 : : void dt();
344 : :
345 : : //! Evaluate whether to save checkpoint/restart
346 : : void evalRestart();
347 : :
348 : : //! Apply boundary conditions
349 : : void BC( tk::real t );
350 : :
351 : : //! Compute gradients for next time step
352 : : void grad();
353 : :
354 : : //! Apply scalar source to solution
355 : : void src();
356 : : };
357 : :
358 : : } // inciter::
|
