Branch data Line data Source code
1 : : // *****************************************************************************
2 : : /*!
3 : : \file src/Inciter/OversetFE.hpp
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 OversetFE for a PDE system with continuous Galerkin FE + RK
9 : : \details OversetFE advances a system of partial differential equations
10 : : using a continuous Galerkin (CG) finite element (FE) spatial discretization
11 : : (using linear shapefunctions on tetrahedron elements) combined with a
12 : : Runge-Kutta (RK) time stepping scheme and overset grids.
13 : :
14 : : There are a potentially large number of OversetFE Charm++ chares created by
15 : : Transporter. Each OversetFE gets a chunk of the full load (part of the mesh)
16 : : and does the same: initializes and advances a number of PDE systems in time.
17 : :
18 : : The implementation uses the Charm++ runtime system and is fully
19 : : asynchronous, overlapping computation and communication. The algorithm
20 : : utilizes the structured dagger (SDAG) Charm++ functionality.
21 : : */
22 : : // *****************************************************************************
23 : : #ifndef OversetFE_h
24 : : #define OversetFE_h
25 : :
26 : : #include <vector>
27 : : #include <map>
28 : :
29 : : #include "Types.hpp"
30 : : #include "Fields.hpp"
31 : : #include "Table.hpp"
32 : : #include "DerivedData.hpp"
33 : : #include "NodeDiagnostics.hpp"
34 : : #include "Ghosts.hpp"
35 : :
36 : : #include "NoWarning/oversetfe.decl.h"
37 : :
38 : : namespace inciter {
39 : :
40 : : //! OversetFE Charm++ chare array used to advance PDEs in time with OversetFE+RK
41 : : class OversetFE : public CBase_OversetFE {
42 : :
43 : : public:
44 : : #if defined(__clang__)
45 : : #pragma clang diagnostic push
46 : : #pragma clang diagnostic ignored "-Wunused-parameter"
47 : : #pragma clang diagnostic ignored "-Wdeprecated-declarations"
48 : : #elif defined(STRICT_GNUC)
49 : : #pragma GCC diagnostic push
50 : : #pragma GCC diagnostic ignored "-Wunused-parameter"
51 : : #pragma GCC diagnostic ignored "-Wdeprecated-declarations"
52 : : #elif defined(__INTEL_COMPILER)
53 : : #pragma warning( push )
54 : : #pragma warning( disable: 1478 )
55 : : #endif
56 : : // Include Charm++ SDAG code. See http://charm.cs.illinois.edu/manuals/html/
57 : : // charm++/manual.html, Sec. "Structured Control Flow: Structured Dagger".
58 : : OversetFE_SDAG_CODE
59 : : #if defined(__clang__)
60 : : #pragma clang diagnostic pop
61 : : #elif defined(STRICT_GNUC)
62 : : #pragma GCC diagnostic pop
63 : : #elif defined(__INTEL_COMPILER)
64 : : #pragma warning( pop )
65 : : #endif
66 : :
67 : : //! Constructor
68 : : explicit OversetFE( const CProxy_Discretization& disc,
69 : : const CProxy_Ghosts& ghostsproxy,
70 : : const std::map< int, std::vector< std::size_t > >& bface,
71 : : const std::map< int, std::vector< std::size_t > >& bnode,
72 : : const std::vector< std::size_t >& triinpoel );
73 : :
74 : : #if defined(__clang__)
75 : : #pragma clang diagnostic push
76 : : #pragma clang diagnostic ignored "-Wundefined-func-template"
77 : : #endif
78 : : //! Migrate constructor
79 : : // cppcheck-suppress uninitMemberVar
80 : 984 : explicit OversetFE( CkMigrateMessage* msg ) : CBase_OversetFE( msg ) {}
81 : : #if defined(__clang__)
82 : : #pragma clang diagnostic pop
83 : : #endif
84 : :
85 : : //! Configure Charm++ custom reduction types initiated from this chare array
86 : : static void registerReducers();
87 : :
88 : : //! Return from migration
89 : : void ResumeFromSync() override;
90 : :
91 : : //! Start setup for solution
92 : : void setup();
93 : :
94 : : //! Receive total box IC volume and set conditions in box
95 : : void box( tk::real v, const std::vector< tk::real >& blkvols );
96 : :
97 : : //! Transfer solution to other solver and mesh if coupled
98 : : void transferSol();
99 : :
100 : : // Start time stepping
101 : : void start();
102 : :
103 : : //! Advance equations to next time step
104 : : void advance( tk::real newdt, std::array< tk::real, 6 > F );
105 : :
106 : : //! Compute left-hand side of transport equations
107 : : void lhs();
108 : :
109 : : //! Receive contributions to duual-face normals on chare boundaries
110 : : void comdfnorm(
111 : : const std::unordered_map< tk::UnsMesh::Edge,
112 : : std::array< tk::real, 3 >,
113 : : tk::UnsMesh::Hash<2>, tk::UnsMesh::Eq<2> >& dfnorm );
114 : :
115 : : //! Receive boundary point normals on chare-boundaries
116 : : void comnorm( const std::unordered_map< int,
117 : : std::unordered_map< std::size_t, std::array< tk::real, 4 > > >& innorm );
118 : :
119 : : //! Receive mesh block information for nodes on chare-boundaries
120 : : void comblk( const std::vector< std::size_t >& gid,
121 : : const std::vector< std::set< std::size_t > >& mb );
122 : :
123 : : //! Receive contributions to gradients on chare-boundaries
124 : : void comChBndGrad( const std::vector< std::size_t >& gid,
125 : : const std::vector< std::vector< tk::real > >& G );
126 : :
127 : : //! Receive contributions to right-hand side vector on chare-boundaries
128 : : void comrhs( const std::vector< std::size_t >& gid,
129 : : const std::vector< std::vector< tk::real > >& R );
130 : :
131 : : //! Update solution at the end of time step
132 : : void update( const tk::Fields& a );
133 : :
134 : : //! Optionally refine/derefine mesh
135 : : void refine( const std::vector< tk::real >& l2res );
136 : :
137 : : //! Extract field output to file
138 : : void extractFieldOutput(
139 : : const std::vector< std::size_t >& /* ginpoel */,
140 : : const tk::UnsMesh::Chunk& /*chunk*/,
141 : : const tk::UnsMesh::Coords& /*coord*/,
142 : : const std::unordered_map< std::size_t, tk::UnsMesh::Edge >& /* addedNodes */,
143 : : const std::unordered_map< std::size_t, std::size_t >& /*addedTets*/,
144 : : const tk::NodeCommMap& /*nodeCommMap*/,
145 : : const std::map< int, std::vector< std::size_t > >& /*bface*/,
146 : : const std::map< int, std::vector< std::size_t > >& /* bnode */,
147 : : const std::vector< std::size_t >& /*triinpoel*/,
148 : : CkCallback /*c*/ ) {}
149 : :
150 : : //! Const-ref access to current solution
151 : : //! \return Const-ref to current solution
152 : : const tk::Fields& solution() const { return m_u; }
153 : :
154 : : //! Evaluate whether to continue with next time step
155 : : void step();
156 : :
157 : : // Evaluate whether to do load balancing
158 : : void evalLB( int nrestart );
159 : :
160 : : //! Evaluate whether to continue with next time step stage
161 : : void stage();
162 : :
163 : : //! Continue to next time step
164 : : void next();
165 : :
166 : : //! Size communication buffers (no-op)
167 : : void resizeComm() {}
168 : :
169 : : //! Setup node-neighborhood (no-op)
170 : : void nodeNeighSetup() {}
171 : :
172 : : //! Receive new mesh from Refiner (no-op)
173 : : void resizePostAMR(
174 : : const std::vector< std::size_t >&,
175 : : const tk::UnsMesh::Chunk&,
176 : : const tk::UnsMesh::Coords&,
177 : : const std::unordered_map< std::size_t, tk::UnsMesh::Edge >&,
178 : : const std::unordered_map< std::size_t, std::size_t >&,
179 : : const std::set< std::size_t >&,
180 : : const std::unordered_map< std::size_t, std::size_t >&,
181 : : const tk::NodeCommMap&,
182 : : const std::map< int, std::vector< std::size_t > >&,
183 : : const std::map< int, std::vector< std::size_t > >&,
184 : : const std::vector< std::size_t >&,
185 : : const std::unordered_map< std::size_t, std::set< std::size_t > >& ) {}
186 : :
187 : : //! Resizing data structures after mesh refinement has completed (no-op)
188 : : void resized() {}
189 : :
190 : : /** @name Charm++ pack/unpack serializer member functions */
191 : : ///@{
192 : : //! \brief Pack/Unpack serialize member function
193 : : //! \param[in,out] p Charm++'s PUP::er serializer object reference
194 : 2959 : void pup( PUP::er &p ) override {
195 : 2959 : p | m_disc;
196 : 2959 : p | m_nsol;
197 : 2959 : p | m_ngrad;
198 : 2959 : p | m_nrhs;
199 : 2959 : p | m_nbnorm;
200 : 2959 : p | m_ndfnorm;
201 : 2959 : p | m_nmblk;
202 : 2959 : p | m_bnode;
203 : 2959 : p | m_bface;
204 : 2959 : p | m_triinpoel;
205 : 2959 : p | m_bndel;
206 : 2959 : p | m_dfnorm;
207 : 2959 : p | m_dfnormc;
208 : 2959 : p | m_dfn;
209 : 2959 : p | m_esup;
210 : 2959 : p | m_psup;
211 : 2959 : p | m_u;
212 : 2959 : p | m_uc;
213 : 2959 : p | m_un;
214 : 2959 : p | m_rhs;
215 : 2959 : p | m_rhsc;
216 : 2959 : p | m_chBndGrad;
217 : 2959 : p | m_dirbc;
218 : 2959 : p | m_chBndGradc;
219 : 2959 : p | m_blank;
220 : : p | m_diag;
221 : 2959 : p | m_bnorm;
222 : 2959 : p | m_bnormc;
223 : 2959 : p | m_symbcnodes;
224 : 2959 : p | m_farfieldbcnodes;
225 : 2959 : p | m_symbctri;
226 : 2959 : p | m_timedepbcnodes;
227 : 2959 : p | m_timedepbcFn;
228 : 2959 : p | m_stage;
229 : 2959 : p | m_boxnodes;
230 : 2959 : p | m_edgenode;
231 : 2959 : p | m_edgeid;
232 : 2959 : p | m_dtp;
233 : 2959 : p | m_tp;
234 : 2959 : p | m_finished;
235 : 2959 : p | m_movedmesh;
236 : 2959 : p | m_nusermeshblk;
237 : 2959 : p | m_nodeblockid;
238 : 2959 : p | m_nodeblockidc;
239 : 2959 : p | m_ixfer;
240 : : p | m_surfForce;
241 : : p | m_surfTorque;
242 : : p | m_centMass;
243 : : p | m_centMassVel;
244 : 2959 : p | m_angVelMesh;
245 : 2959 : }
246 : : //! \brief Pack/Unpack serialize operator|
247 : : //! \param[in,out] p Charm++'s PUP::er serializer object reference
248 : : //! \param[in,out] i OversetFE object reference
249 : : friend void operator|( PUP::er& p, OversetFE& i ) { i.pup(p); }
250 : : //@}
251 : :
252 : : private:
253 : : using ncomp_t = tk::ncomp_t;
254 : :
255 : : //! Discretization proxy
256 : : CProxy_Discretization m_disc;
257 : : //! Counter for high order solution vector nodes updated
258 : : std::size_t m_nsol;
259 : : //! Counter for nodal gradients updated
260 : : std::size_t m_ngrad;
261 : : //! Counter for right-hand side vector nodes updated
262 : : std::size_t m_nrhs;
263 : : //! Counter for receiving boundary point normals
264 : : std::size_t m_nbnorm;
265 : : //! Counter for receiving dual-face normals on chare-boundary edges
266 : : std::size_t m_ndfnorm;
267 : : //! Counter for mesh block nodes updated
268 : : std::size_t m_nmblk;
269 : : //! Boundary node lists mapped to side set ids used in the input file
270 : : std::map< int, std::vector< std::size_t > > m_bnode;
271 : : //! Boundary face lists mapped to side set ids used in the input file
272 : : std::map< int, std::vector< std::size_t > > m_bface;
273 : : //! Boundary triangle face connecitivity where BCs are set by user
274 : : std::vector< std::size_t > m_triinpoel;
275 : : //! Elements along mesh boundary
276 : : std::vector< std::size_t > m_bndel;
277 : : //! Dual-face normals along edges
278 : : std::unordered_map< tk::UnsMesh::Edge, std::array< tk::real, 3 >,
279 : : tk::UnsMesh::Hash<2>, tk::UnsMesh::Eq<2> > m_dfnorm;
280 : : //! Receive buffer for dual-face normals along chare-boundary edges
281 : : std::unordered_map< tk::UnsMesh::Edge, std::array< tk::real, 3 >,
282 : : tk::UnsMesh::Hash<2>, tk::UnsMesh::Eq<2> > m_dfnormc;
283 : : //! Streamable dual-face normals
284 : : std::vector< tk::real > m_dfn;
285 : : //! El;ements surrounding points
286 : : std::pair< std::vector< std::size_t >, std::vector< std::size_t > > m_esup;
287 : : //! Points surrounding points
288 : : std::pair< std::vector< std::size_t >, std::vector< std::size_t > > m_psup;
289 : : //! Unknown/solution vector at mesh nodes
290 : : tk::Fields m_u;
291 : : //! \brief Copy of unknown/solution vector at mesh nodes for m2m transfer,
292 : : //! appended with a solution-transfer-flag. This flag indicates
293 : : //! appropriate solution transfers for the overset procedure.
294 : : //! Value 0: Do not transfer solution, 1: transfer sol, 2: blank nodes
295 : : //! TODO: avoid creating this copy
296 : : tk::Fields m_uc;
297 : : //! Unknown/solution vector at mesh nodes at previous time
298 : : tk::Fields m_un;
299 : : //! Right-hand side vector (for the high order system)
300 : : tk::Fields m_rhs;
301 : : //! Receive buffer for communication of the right hand side
302 : : //! \details Key: global node id, value: rhs for all scalar components per
303 : : //! node.
304 : : std::unordered_map< std::size_t, std::vector< tk::real > > m_rhsc;
305 : : //! Nodal gradients at chare-boundary nodes
306 : : tk::Fields m_chBndGrad;
307 : : //! Boundary conditions evaluated and assigned to local mesh node IDs
308 : : //! \details Vector of pairs of bool and boundary condition value associated
309 : : //! to local mesh node IDs at which the user has set Dirichlet boundary
310 : : //! conditions for all PDEs integrated. The bool indicates whether the BC
311 : : //! is set at the node for that component the if true, the real value is
312 : : //! the increment (from t to dt) in the BC specified for a component.
313 : : std::unordered_map< std::size_t,
314 : : std::vector< std::pair< bool, tk::real > > > m_dirbc;
315 : : //! Receive buffer for communication of the nodal gradients
316 : : //! \details Key: chare id, value: gradients for all scalar components per
317 : : //! node
318 : : std::unordered_map< std::size_t, std::vector< tk::real > > m_chBndGradc;
319 : : //! Blanking coefficient for overset, indicating hole in the background mesh
320 : : std::vector< tk::real > m_blank;
321 : : //! Diagnostics object
322 : : NodeDiagnostics m_diag;
323 : : //! Face normals in boundary points associated to side sets
324 : : //! \details Key: local node id, value: unit normal and inverse distance
325 : : //! square between face centroids and points, outer key: side set id
326 : : std::unordered_map< int,
327 : : std::unordered_map< std::size_t, std::array< tk::real, 4 > > > m_bnorm;
328 : : //! \brief Receive buffer for communication of the boundary point normals
329 : : //! associated to side sets
330 : : //! \details Key: global node id, value: normals (first 3 components),
331 : : //! inverse distance squared (4th component), outer key, side set id
332 : : std::unordered_map< int,
333 : : std::unordered_map< std::size_t, std::array< tk::real, 4 > > > m_bnormc;
334 : : //! Unique set of nodes at which symmetry BCs are set
335 : : std::unordered_set< std::size_t > m_symbcnodes;
336 : : //! Unique set of nodes at which farfield BCs are set
337 : : std::unordered_set< std::size_t > m_farfieldbcnodes;
338 : : //! Vector with 1 at symmetry BC boundary triangles
339 : : std::vector< int > m_symbctri;
340 : : //! \brief Unique set of nodes at which time dependent BCs are set
341 : : // for each time dependent BC
342 : : std::vector< std::unordered_set< std::size_t > > m_timedepbcnodes;
343 : : //! \brief User defined discrete function of time used in the time dependent
344 : : // BCs associated with (index in vector) the number of distinct time
345 : : // dependent BCs specified. This index is the same as the index in
346 : : // m_timedepbcnodes.
347 : : std::vector< tk::Table<5> > m_timedepbcFn;
348 : : //! Runge-Kutta stage counter
349 : : std::size_t m_stage;
350 : : //! Mesh node ids at which user-defined box ICs are defined (multiple boxes)
351 : : std::vector< std::unordered_set< std::size_t > > m_boxnodes;
352 : : //! Local node IDs of edges
353 : : std::vector< std::size_t > m_edgenode;
354 : : //! Edge ids in the order of access
355 : : std::vector< std::size_t > m_edgeid;
356 : : //! Time step size for each mesh node
357 : : std::vector< tk::real > m_dtp;
358 : : //! Physical time for each mesh node
359 : : std::vector< tk::real > m_tp;
360 : : //! True in the last time step
361 : : int m_finished;
362 : : //! True if overset mesh moved
363 : : int m_movedmesh;
364 : : //! Number of mesh-blocks with user-defined ICs
365 : : std::size_t m_nusermeshblk;
366 : : //! Local node ids associated with mesh block ids
367 : : std::unordered_map< std::size_t, std::set< std::size_t > > m_nodeblockid;
368 : : //! Receive buffer for communication of the mesh block ids
369 : : //! \details Key: mesh block id, value: set of global node ids for nodes
370 : : //! in this mesh block.
371 : : std::unordered_map< std::size_t, std::set< std::size_t > > m_nodeblockidc;
372 : : //! Counter for two-way transfer
373 : : std::size_t m_ixfer;
374 : : //! Surface force on mesh for rigid body motion
375 : : std::array< tk::real, 3 > m_surfForce;
376 : : //! Surface torque on mesh for rigid body motion
377 : : std::array< tk::real, 3 > m_surfTorque;
378 : : //! Center of mass of rigid body
379 : : std::array< tk::real, 3 > m_centMass;
380 : : //! Velocity of the center of mass of rigid body
381 : : std::array< tk::real, 3 > m_centMassVel;
382 : : //! Angular velocity of the rigid body
383 : : tk::real m_angVelMesh;
384 : :
385 : : //! Access bound Discretization class pointer
386 : 391959 : Discretization* Disc() const {
387 : : Assert( m_disc[ thisIndex ].ckLocal() != nullptr, "ckLocal() null" );
388 : 783918 : return m_disc[ thisIndex ].ckLocal();
389 : : }
390 : :
391 : : //! Compute normal of dual-mesh associated to edge
392 : : std::array< tk::real, 3 >
393 : : edfnorm( const tk::UnsMesh::Edge& edge,
394 : : const std::unordered_map< tk::UnsMesh::Edge,
395 : : std::vector< std::size_t >,
396 : : tk::UnsMesh::Hash<2>, tk::UnsMesh::Eq<2> >& esued ) const;
397 : :
398 : : //! Compute chare-boundary edges
399 : : void bndEdges();
400 : :
401 : : //! Start (re-)computing boundare point-, and dual-face normals
402 : : void norm();
403 : :
404 : : //! Compute dual-face normals associated to edges
405 : : void dfnorm();
406 : :
407 : : //! Compute boundary point normals
408 : : void
409 : : bnorm( const std::unordered_map< int,
410 : : std::unordered_set< std::size_t > >& bcnodes );
411 : :
412 : : //! Set flags informing solution transfer decisions
413 : : void setTransferFlags(
414 : : std::size_t dirn );
415 : :
416 : : //! \brief Finish computing dual-face and boundary point normals and apply
417 : : //! boundary conditions on the initial conditions
418 : : void normfinal();
419 : :
420 : : //! Continue setup for solution, after communication for mesh blocks
421 : : void continueSetup();
422 : :
423 : : //! Output mesh field data and continue to next time step
424 : : void out();
425 : :
426 : : //! Output mesh-based fields to file
427 : : void writeFields( CkCallback c );
428 : :
429 : : //! Combine own and communicated contributions to normals
430 : : void mergelhs();
431 : :
432 : : //! Compute gradients
433 : : void chBndGrad();
434 : :
435 : : //! Compute righ-hand side vector of transport equations
436 : : void rhs();
437 : :
438 : : //! Advance systems of equations
439 : : void solve();
440 : :
441 : : //! Compute time step size
442 : : void dt();
443 : :
444 : : //! Evaluate whether to save checkpoint/restart
445 : : void evalRestart();
446 : :
447 : : //! Query/update boundary-conditions-related data structures from user input
448 : : void getBCNodes();
449 : :
450 : : // \brief Apply the transferred solution to the solution vector based on
451 : : // transfer flags previously set up
452 : : void applySolTransfer( std::size_t dirn );
453 : :
454 : : //! Apply boundary conditions
455 : : void BC();
456 : : };
457 : :
458 : : } // inciter::
459 : :
460 : : #endif // OversetFE_h
|
