Branch data Line data Source code
1 : : // *****************************************************************************
2 : : /*!
3 : : \file src/Inciter/Discretization.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 : : \details Data and functionality common to all discretization schemes
10 : : The Discretization class contains data and functionality common to all
11 : : discretization schemes.
12 : : */
13 : : // *****************************************************************************
14 : : #pragma once
15 : :
16 : : #include "Types.hpp"
17 : : #include "Timer.hpp"
18 : : #include "Fields.hpp"
19 : : #include "PUPUtil.hpp"
20 : : #include "PDFReducer.hpp"
21 : : #include "UnsMesh.hpp"
22 : : #include "History.hpp"
23 : :
24 : : #include "NoWarning/discretization.decl.h"
25 : : #include "NoWarning/refiner.decl.h"
26 : :
27 : : namespace inciter {
28 : :
29 : : //! \brief Discretization Charm++ chare array holding common functinoality to
30 : : //! all discretization schemes
31 : : class Discretization : public CBase_Discretization {
32 : :
33 : : public:
34 : : #if defined(__clang__)
35 : : #pragma clang diagnostic push
36 : : #pragma clang diagnostic ignored "-Wunused-parameter"
37 : : #pragma clang diagnostic ignored "-Wdeprecated-declarations"
38 : : #elif defined(STRICT_GNUC)
39 : : #pragma GCC diagnostic push
40 : : #pragma GCC diagnostic ignored "-Wunused-parameter"
41 : : #pragma GCC diagnostic ignored "-Wdeprecated-declarations"
42 : : #elif defined(__INTEL_COMPILER)
43 : : #pragma warning( push )
44 : : #pragma warning( disable: 1478 )
45 : : #endif
46 : : // Include Charm++ SDAG code. See http://charm.cs.illinois.edu/manuals/html/
47 : : // charm++/manual.html, Sec. "Structured Control Flow: Structured Dagger".
48 : : Discretization_SDAG_CODE
49 : : #if defined(__clang__)
50 : : #pragma clang diagnostic pop
51 : : #elif defined(STRICT_GNUC)
52 : : #pragma GCC diagnostic pop
53 : : #elif defined(__INTEL_COMPILER)
54 : : #pragma warning( pop )
55 : : #endif
56 : :
57 : : //! Constructor
58 : : explicit
59 : : Discretization(
60 : : std::size_t meshid,
61 : : const CProxy_Transporter& transporter,
62 : : const tk::CProxy_MeshWriter& meshwriter,
63 : : const tk::UnsMesh::CoordMap& coordmap,
64 : : const tk::UnsMesh::Chunk& el,
65 : : const std::map< int, std::unordered_set< std::size_t > >& nodeCommMap,
66 : : int nc );
67 : :
68 : : #if defined(__clang__)
69 : : #pragma clang diagnostic push
70 : : #pragma clang diagnostic ignored "-Wundefined-func-template"
71 : : #endif
72 : : //! Migrate constructor
73 : : // cppcheck-suppress uninitMemberVar
74 : 10748 : explicit Discretization( CkMigrateMessage* m ) : CBase_Discretization( m )
75 : 5374 : {}
76 : : #if defined(__clang__)
77 : : #pragma clang diagnostic pop
78 : : #endif
79 : :
80 : : //! Configure Charm++ reduction types
81 : : static void registerReducers();
82 : :
83 : : //! Resize mesh data structures after mesh refinement
84 : : void resizePostAMR(
85 : : const tk::UnsMesh::Chunk& chunk,
86 : : const tk::UnsMesh::Coords& coord,
87 : : const std::unordered_map< int, std::unordered_set< std::size_t > >&
88 : : nodeCommMap,
89 : : const std::set< std::size_t >& removedNodes );
90 : :
91 : : //! Get ready for (re-)computing/communicating nodal volumes
92 : : void startvol();
93 : :
94 : : //! Sum mesh volumes to nodes, start communicating them on chare-boundaries
95 : : void vol();
96 : :
97 : : //! Set Refiner Charm++ proxy
98 : : void setRefiner( const CProxy_Refiner& ref );
99 : :
100 : : //! Collect nodal volumes across chare boundaries
101 : : void comvol( int c,
102 : : const std::vector< std::size_t >& gid,
103 : : const std::vector< tk::real >& nodevol );
104 : :
105 : : //! Sum mesh volumes and contribute own mesh volume to total volume
106 : : void totalvol();
107 : :
108 : : //! Compute mesh cell statistics
109 : : void stat( tk::real mesh_volume );
110 : :
111 : : //! Decide whether to start the deactivation procedure in this time step
112 : : bool deastart();
113 : :
114 : : //! Decide whether to run deactivation procedure in this time step
115 : : bool deactivate() const;
116 : :
117 : : //! Receive deactivation report
118 : : void deactivated( int n );
119 : :
120 : : //! Decide whether to do load balancing this time step
121 : : bool lb() const;
122 : :
123 : : //! Compute total box IC volume
124 : : void boxvol();
125 : :
126 : : /** @name Accessors */
127 : : ///@{
128 : : //! Coordinates accessor as const-ref
129 : : const tk::UnsMesh::Coords& Coord() const { return m_coord; }
130 : : //! Coordinates accessor as reference
131 [ + - ][ + - ]: 126183 : tk::UnsMesh::Coords& Coord() { return m_coord; }
132 : :
133 : : //! Global ids accessors as const-ref
134 [ + + ][ + - ]: 9396 : const std::vector< std::size_t >& Gid() const { return m_gid; }
[ - - ]
135 : :
136 : : //! Local ids accessors as const-ref
137 : : const std::unordered_map< std::size_t, std::size_t >& Lid() const
138 [ + - ][ + + ]: 762256 : { return m_lid; }
[ + + ][ + + ]
[ + - ][ - - ]
[ - - ]
139 : :
140 : : //! Tetrahedron element connectivity (with local ids) accessors as const-ref
141 [ + - ][ + - ]: 24142 : const std::vector< std::size_t >& Inpoel() const { return m_inpoel; }
[ + - ][ + - ]
[ + - ][ - - ]
[ - - ][ - - ]
[ - - ]
142 : :
143 : : //! Mesh chunk accessor as const-ref
144 : : const tk::UnsMesh::Chunk& Chunk() const { return m_el; }
145 : :
146 : : //! Total mesh volume accessor
147 : : tk::real meshvol() const { return m_meshvol; }
148 : :
149 : : //! Nodal mesh volume accessors const-ref
150 : 48240 : const std::vector< tk::real >& V() const { return m_v; }
151 : : //! Nodal mesh volumes at current time step accessors as const-ref
152 [ + - ][ - - ]: 729 : const std::vector< tk::real >& Vol() const { return m_vol; }
153 : : //! Access chare-volume contributions per chare per chare-boundary node
154 : : const std::unordered_map< int,
155 : : std::unordered_map< std::size_t, tk::real > >& Cvolc() const {
156 : : return m_cvolc;
157 : : }
158 : :
159 : : //! Access number of chares
160 [ - + ]: 190 : int nchare() const { return m_nchare; }
161 : :
162 : : //! Set 'initial' flag
163 : : //! \param[in] i Value to put in 'initial'
164 : 2427 : void Initial( std::size_t i ) { m_initial = i; }
165 : : //! Query 'initial' flag
166 : : //! \return True during setup, false durign time stepping
167 [ + + ][ + + ]: 13473 : bool Initial() const { return m_initial; }
[ + + ]
168 : :
169 : : //! History points data accessor as const-ref
170 : : const std::vector< HistData >& Hist() const { return m_histdata; }
171 : :
172 : : //! Box volume accessor
173 : : tk::real& Boxvol() { return m_boxvol; }
174 : :
175 : : //! Mesh ID accessor
176 [ + - ][ + - ]: 39966 : std::size_t MeshId() const { return m_meshid; }
[ + - ]
177 : :
178 : : //! Time step size accessor
179 [ + - ][ + - ]: 516370 : tk::real Dt() const { return m_dt; }
[ + - ]
180 : : //! Time step size at previous time step accessor
181 : : tk::real Dtn() const { return m_dtn; }
182 : : //! Physical time accessor
183 [ + - ][ + + ]: 1250184 : tk::real T() const { return m_t; }
[ + + ][ + + ]
184 : : //! Iteration count accessor
185 [ + + ][ + + ]: 144796 : uint64_t It() const { return m_it; }
[ + - ][ + + ]
[ + + ][ + - ]
186 : :
187 : : //! Non-const-ref refinement iteration count accessor
188 : : uint64_t& Itr() { return m_itr; }
189 : : //! Non-const-ref field-output iteration count accessor
190 : : uint64_t& Itf() { return m_itf; }
191 : :
192 : : //! Non-const-ref number of restarts accessor
193 : : int& Nrestart() { return m_nrestart; }
194 : :
195 : : //! Timer accessor as const-ref
196 : : const tk::Timer& Timer() const { return m_timer; }
197 : : //! Timer accessor as non-const-ref
198 : 2427 : tk::Timer& Timer() { return m_timer; }
199 : :
200 : : //! Accessor to flag indicating if the mesh was refined as a value
201 : : int refined() const { return m_refined; }
202 : : //! Accessor to flag indicating if the mesh was refined as non-const-ref
203 : : int& refined() { return m_refined; }
204 : :
205 : : //! Transporter proxy accessor as const-ref
206 : : const CProxy_Transporter& Tr() const { return m_transporter; }
207 : : //! Transporter proxy accessor as non-const-ref
208 : : CProxy_Transporter& Tr() { return m_transporter; }
209 : :
210 : : //! Access bound Refiner class pointer
211 : 0 : Refiner* Ref() const {
212 : : Assert( m_refiner[ thisIndex ].ckLocal() != nullptr,
213 : : "Refiner ckLocal() null" );
214 : 0 : return m_refiner[ thisIndex ].ckLocal();
215 : : }
216 : :
217 : : //! Node communication map accessor as const-ref
218 : : const std::unordered_map< int, std::unordered_set< std::size_t > >&
219 [ + - ][ + - ]: 1055 : NodeCommMap() const { return m_nodeCommMap; }
[ - - ][ - - ]
220 : :
221 : : //! IC boxnodes accessor
222 : : const std::vector< std::unordered_set< std::size_t > >& BoxNodes() const
223 [ + - ]: 80259 : { return m_boxnodes; }
224 : : //@}
225 : :
226 : : //! Set time step size
227 : : void setdt( tk::real newdt );
228 : :
229 : : //! Prepare for next step
230 : : void next();
231 : :
232 : : //! Otput one-liner status report
233 : : void status();
234 : :
235 : : //! Construct history output filename
236 : : std::string histfilename( const std::string& id,
237 : : std::streamsize precision );
238 : :
239 : : //! Output headers for time history files (one for each point)
240 : : void histheader( std::vector< std::string >&& names );
241 : :
242 : : //! Output time history for a time step
243 : : void history( std::vector< std::vector< tk::real > >&& data );
244 : :
245 : : //! Output mesh and fields data (solution dump) to file(s)
246 : : void write( const std::vector< std::size_t >& inpoel,
247 : : const tk::UnsMesh::Coords& coord,
248 : : const std::map< int, std::vector< std::size_t > >& bface,
249 : : const std::map< int, std::vector< std::size_t > >& bnode,
250 : : const std::vector< std::size_t >& triinpoel,
251 : : const std::vector< std::string>& elemfieldnames,
252 : : const std::vector< std::string>& nodefieldnames,
253 : : const std::vector< std::string>& elemsurfnames,
254 : : const std::vector< std::string>& nodesurfnames,
255 : : const std::vector< std::vector< tk::real > >& elemfields,
256 : : const std::vector< std::vector< tk::real > >& nodefields,
257 : : const std::vector< std::vector< tk::real > >& elemsurfs,
258 : : const std::vector< std::vector< tk::real > >& nodesurfs,
259 : : CkCallback c );
260 : :
261 : : //! Zero grind-timer
262 : : void grindZero();
263 : :
264 : : //! Detect if just returned from a checkpoint and if so, zero timers
265 : : bool restarted( int nrestart );
266 : :
267 : : //! Remap mesh data due to new local ids
268 : : void remap( const std::unordered_map< std::size_t, std::size_t >& map );
269 : :
270 : : //! Decide if field output iteration count interval is hit
271 : : bool fielditer() const;
272 : : //! Decide if field output physics time interval is hit
273 : : bool fieldtime() const;
274 : : //! Decide if physics time falls into a field output time range
275 : : bool fieldrange() const;
276 : :
277 : : //! Decide if history output iteration count interval is hit
278 : : bool histiter() const;
279 : : //! Decide if history output physics time interval is hit
280 : : bool histtime() const;
281 : : //! Decide if physics time falls into a history output time range
282 : : bool histrange() const;
283 : :
284 : : //! Decide if integral output iteration count interval is hit
285 : : bool integiter() const;
286 : : //! Decide if integral output physics time interval is hit
287 : : bool integtime() const;
288 : : //! Decide if physics time falls into a integral output time range
289 : : bool integrange() const;
290 : :
291 : : //! Decide if this is the last time step
292 : : bool finished() const;
293 : :
294 : : //! Update residual (during convergence to steady state)
295 : : void residual( tk::real r );
296 : :
297 : : //! Update number of pressure linear solve iterations taken
298 : : void pit( std::size_t it );
299 : :
300 : : //! Update number of momentum/transport linear solve iterations taken
301 : : void mit( std::size_t it );
302 : :
303 : : //! Set number of mesh points (across all meshes)
304 : : void npoin( std::size_t n );
305 : :
306 : : /** @name Charm++ pack/unpack serializer member functions */
307 : : ///@{
308 : : //! \brief Pack/Unpack serialize member function
309 : : //! \param[in,out] p Charm++'s PUP::er serializer object reference
310 : 16804 : void pup( PUP::er &p ) override {
311 : 16804 : p | m_npoin;
312 : 16804 : p | m_meshid;
313 : 16804 : p | m_nchare;
314 : 16804 : p | m_it;
315 : 16804 : p | m_itr;
316 : 16804 : p | m_itf;
317 : 16804 : p | m_initial;
318 : 16804 : p | m_t;
319 : 16804 : p | m_lastDumpTime;
320 : 16804 : p | m_physFieldFloor;
321 : 16804 : p | m_physHistFloor;
322 : 16804 : p | m_physIntegFloor;
323 : 16804 : p | m_rangeFieldFloor;
324 : 16804 : p | m_rangeHistFloor;
325 : 16804 : p | m_rangeIntegFloor;
326 : 16804 : p | m_dt;
327 : 16804 : p | m_dtn;
328 : 16804 : p | m_nvol;
329 : 16804 : p | m_boxnodes;
330 : : p | m_transporter;
331 : : p | m_meshwriter;
332 : : p | m_refiner;
333 : 16804 : p | m_el;
334 [ + + ]: 16804 : if (p.isUnpacking()) {
335 : 5374 : m_inpoel = std::get< 0 >( m_el );
336 : 5374 : m_gid = std::get< 1 >( m_el );
337 [ - + ]: 5374 : m_lid = std::get< 2 >( m_el );
338 : : }
339 : : p | m_coord;
340 : 16804 : p | m_nodeCommMap;
341 : 16804 : p | m_meshvol;
342 : 16804 : p | m_v;
343 : 16804 : p | m_vol;
344 : 16804 : p | m_volc;
345 : 16804 : p | m_cvolc;
346 : 16804 : p | m_boxvol;
347 : 16804 : p | m_timer;
348 : 16804 : p | m_refined;
349 : 16804 : p( reinterpret_cast<char*>(&m_prevstatus), sizeof(Clock::time_point) );
350 : 16804 : p | m_nrestart;
351 : 16804 : p | m_histdata;
352 : 16804 : p | m_res;
353 : 16804 : p | m_res0;
354 : 16804 : p | m_res1;
355 : 16804 : p | m_dea;
356 : 16804 : p | m_deastarted;
357 : 16804 : p | m_pit;
358 : 16804 : p | m_mit;
359 : 16804 : }
360 : : //! \brief Pack/Unpack serialize operator|
361 : : //! \param[in,out] p Charm++'s PUP::er serializer object reference
362 : : //! \param[in,out] i Discretization object reference
363 : : friend void operator|( PUP::er& p, Discretization& i ) { i.pup(p); }
364 : : //@}
365 : :
366 : : private:
367 : : // Shorthand for clock, setting an internal clock type
368 : : using Clock = std::chrono::high_resolution_clock;
369 : :
370 : : //! Total number of mesh points (across all meshes)
371 : : std::size_t m_npoin;
372 : : //! Mesh ID
373 : : std::size_t m_meshid;
374 : : //! Total number of Discretization chares
375 : : int m_nchare;
376 : : //! Iteration count
377 : : uint64_t m_it;
378 : : //! Iteration count with mesh refinement
379 : : //! \details Used as the restart sequence number {RS} in saving output in
380 : : //! an ExodusII sequence
381 : : //! \see https://www.paraview.org/Wiki/Restarted_Simulation_Readers
382 : : uint64_t m_itr;
383 : : //! Field output iteration count without mesh refinement
384 : : //! \details Counts the number of field outputs to file during two
385 : : //! time steps with mesh efinement
386 : : uint64_t m_itf;
387 : : //! Flag that is nonzero during setup and zero during time stepping
388 : : std::size_t m_initial;
389 : : //! Physical time
390 : : tk::real m_t;
391 : : //! Physics time at last field output
392 : : tk::real m_lastDumpTime;
393 : : //! Recent floor of physics time divided by field output interval time
394 : : tk::real m_physFieldFloor;
395 : : //! Recent floor of physics time divided by history output interval time
396 : : tk::real m_physHistFloor;
397 : : //! Recent floor of physics time divided by integral output interval time
398 : : tk::real m_physIntegFloor;
399 : : //! Recent floors of physics time divided by field output time for ranges
400 : : std::vector< tk::real > m_rangeFieldFloor;
401 : : //! Recent floors of physics time divided by history output time for ranges
402 : : std::vector< tk::real > m_rangeHistFloor;
403 : : //! Recent floors of physics time divided by integral output time for ranges
404 : : std::vector< tk::real > m_rangeIntegFloor;
405 : : //! Physical time step size
406 : : tk::real m_dt;
407 : : //! Physical time step size at the previous time step
408 : : tk::real m_dtn;
409 : : //! \brief Number of chares from which we received nodal volume
410 : : //! contributions on chare boundaries
411 : : std::size_t m_nvol;
412 : : //! List of nodes at which box user ICs are set for each IC box
413 : : std::vector< std::unordered_set< std::size_t > > m_boxnodes;
414 : : //! Transporter proxy
415 : : CProxy_Transporter m_transporter;
416 : : //! Mesh writer proxy
417 : : tk::CProxy_MeshWriter m_meshwriter;
418 : : //! Mesh refiner proxy
419 : : CProxy_Refiner m_refiner;
420 : : //! \brief Elements of the mesh chunk we operate on
421 : : //! \details Initialized by the constructor. The first vector is the element
422 : : //! connectivity (local IDs), the second vector is the global node IDs of
423 : : //! owned elements, while the third one is a map of global->local node
424 : : //! IDs.
425 : : tk::UnsMesh::Chunk m_el;
426 : : //! Alias to element connectivity
427 : : std::vector< std::size_t >& m_inpoel = std::get<0>( m_el );
428 : : //! Alias to global node IDs of owned elements
429 : : std::vector< std::size_t >& m_gid = std::get<1>( m_el );
430 : : //! \brief Alias to local node ids associated to the global ones of owned
431 : : //! elements
432 : : std::unordered_map< std::size_t, std::size_t >& m_lid = std::get<2>( m_el );
433 : : //! Mesh point coordinates
434 : : tk::UnsMesh::Coords m_coord;
435 : : //! \brief Global mesh node IDs bordering the mesh chunk held by fellow
436 : : //! Discretization chares associated to their chare IDs
437 : : std::unordered_map< int, std::unordered_set< std::size_t > > m_nodeCommMap;
438 : : //! Total mesh volume
439 : : tk::real m_meshvol;
440 : : //! Nodal mesh volumes
441 : : //! \details This is the volume of the mesh associated to nodes of owned
442 : : //! elements (sum of surrounding cell volumes / 4) without contributions
443 : : //! from other chares on chare-boundaries
444 : : std::vector< tk::real > m_v;
445 : : //! Volume of nodes
446 : : //! \details This is the volume of the mesh associated to nodes of owned
447 : : //! elements (sum of surrounding cell volumes / 4) with contributions from
448 : : //! other chares on chare-boundaries
449 : : std::vector< tk::real > m_vol;
450 : : //! Receive buffer for volume of nodes (with global node id as key)
451 : : //! \details This is a communication buffer used to compute the volume of
452 : : //! the mesh associated to nodes of owned elements (sum of surrounding
453 : : //! cell volumes / 4) with contributions from other chares on
454 : : //! chare-boundaries.
455 : : std::unordered_map< std::size_t, tk::real > m_volc;
456 : : //! Chare-volume contributions per chare per chare-boundary node
457 : : //! \details Outer key: chare id, value: local node id and node volume
458 : : std::unordered_map< int, std::unordered_map< std::size_t, tk::real > >
459 : : m_cvolc;
460 : : //! Volume of user-defined box IC
461 : : tk::real m_boxvol;
462 : : //! Timer measuring a time step
463 : : tk::Timer m_timer;
464 : : //! 1 if mesh was refined in a time step, 0 if it was not
465 : : int m_refined;
466 : : //! Time point storing clock state at status()
467 : : Clock::time_point m_prevstatus;
468 : : //! Number of times restarted
469 : : int m_nrestart;
470 : : //! Data at history point locations
471 : : std::vector< HistData > m_histdata;
472 : : //! Current residual (during convergence to steady state)
473 : : tk::real m_res;
474 : : //! Residual at previous ETA calcuation (during convergence to steady state)
475 : : tk::real m_res0;
476 : : //! Residual at next ETA calcuation (during convergence to steady state)
477 : : tk::real m_res1;
478 : : //! Numberf of deactived chares
479 : : int m_dea;
480 : : //! Flag: 1 if deactivation procedure has started
481 : : int m_deastarted;
482 : : //! Numberf of recent pressure solve iterations taken
483 : : std::size_t m_pit;
484 : : //! Numberf of recent momentum/transport solve iterations taken
485 : : std::size_t m_mit;
486 : :
487 : : //! Set mesh coordinates based on coordinates map
488 : : tk::UnsMesh::Coords setCoord( const tk::UnsMesh::CoordMap& coordmap );
489 : : };
490 : :
491 : : } // inciter::
|
