Branch data Line data Source code
1 : : // *****************************************************************************
2 : : /*!
3 : : \file src/Inciter/Sorter.cpp
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 Mesh sorter for global distributed mesh reordering
10 : : \see Sorter.h for more info.
11 : : */
12 : : // *****************************************************************************
13 : :
14 : : #include <vector>
15 : : #include <algorithm>
16 : :
17 : : #include "Sorter.hpp"
18 : : #include "Reorder.hpp"
19 : : #include "DerivedData.hpp"
20 : : #include "InciterConfig.hpp"
21 : :
22 : : namespace inciter {
23 : :
24 : : extern ctr::Config g_cfg;
25 : :
26 : : } // inciter::
27 : :
28 : : using inciter::Sorter;
29 : :
30 : 2459 : Sorter::Sorter( std::size_t meshid,
31 : : const CProxy_Transporter& transporter,
32 : : const tk::CProxy_MeshWriter& meshwriter,
33 : : const tk::SorterCallback& cbs,
34 : : const CProxy_Discretization& discretization,
35 : : const CProxy_RieCG& riecg,
36 : : const CProxy_LaxCG& laxcg,
37 : : const CProxy_ZalCG& zalcg,
38 : : const CProxy_KozCG& kozcg,
39 : : const CProxy_ChoCG& chocg,
40 : : const CProxy_LohCG& lohcg,
41 : : const tk::CProxy_ConjugateGradients& cgpre,
42 : : const tk::CProxy_ConjugateGradients& cgmom,
43 : : CkCallback reorderRefiner,
44 : : const std::vector< std::size_t >& ginpoel,
45 : : const tk::UnsMesh::CoordMap& coordmap,
46 : : const tk::UnsMesh::Chunk& el,
47 : : const std::map< int, std::vector< std::size_t > >& bface,
48 : : const std::vector< std::size_t >& triinpoel,
49 : : const std::map< int, std::vector< std::size_t > >& bnode,
50 : 2459 : int nchare ) :
51 : 2459 : m_meshid( meshid ),
52 [ + - ]: 2459 : m_host( transporter ),
53 [ + - ]: 2459 : m_meshwriter( meshwriter ),
54 : 2459 : m_cbs( cbs ),
55 [ + - ]: 2459 : m_discretization( discretization ),
56 [ + - ]: 2459 : m_riecg( riecg ),
57 [ + - ]: 2459 : m_laxcg( laxcg ),
58 [ + - ]: 2459 : m_zalcg( zalcg ),
59 [ + - ]: 2459 : m_kozcg( kozcg ),
60 [ + - ]: 2459 : m_chocg( chocg ),
61 [ + - ]: 2459 : m_lohcg( lohcg ),
62 [ + - ]: 2459 : m_cgpre( cgpre ),
63 [ + - ]: 2459 : m_cgmom( cgmom ),
64 : 2459 : m_reorderRefiner( reorderRefiner ),
65 [ + - ]: 2459 : m_ginpoel( ginpoel ),
66 [ + - ]: 2459 : m_coordmap( coordmap ),
67 [ + - ]: 2459 : m_el( el ),
68 : 2459 : m_nbnd( 0 ),
69 [ + - ]: 2459 : m_bface( bface ),
70 [ + - ]: 2459 : m_triinpoel( triinpoel ),
71 [ + - ]: 2459 : m_bnode( bnode ),
72 : 2459 : m_nchare( nchare ),
73 [ + - ]: 2459 : m_nodeset( begin(ginpoel), end(ginpoel) ),
74 : 2459 : m_noffset( 0 ),
75 : 2459 : m_nodech(),
76 : 2459 : m_chnode(),
77 : 2459 : m_nodeCommMap(),
78 : 2459 : m_reordcomm(),
79 : 2459 : m_start( 0 ),
80 : 2459 : m_newnodes(),
81 : 2459 : m_newcoordmap(),
82 : 2459 : m_reqnodes(),
83 : 2459 : m_lower( 0 ),
84 : 2459 : m_upper( 0 )
85 : : // *****************************************************************************
86 : : // Constructor: prepare owned mesh node IDs for reordering
87 : : //! \param[in] meshid Mesh ID
88 : : //! \param[in] transporter Transporter (host) Charm++ proxy
89 : : //! \param[in] meshwriter Mesh writer Charm++ proxy
90 : : //! \param[in] cbs Charm++ callbacks for Sorter
91 : : //! \param[in] discretization Discretization Charm++ proxy
92 : : //! \param[in] riecg RieCG Charm++ proxy
93 : : //! \param[in] laxcg RieCG Charm++ proxy
94 : : //! \param[in] zalcg ZalCG Charm++ proxy
95 : : //! \param[in] kozcg KozCG Charm++ proxy
96 : : //! \param[in] chocg ChoCG Charm++ proxy
97 : : //! \param[in] lohcg LohCG Charm++ proxy
98 : : //! \param[in] cgpre ConjugateGradients Charm++ proxy for pressure solve
99 : : //! \param[in] cgmom ConjugateGradients Charm++ proxy for momentum solve
100 : : //! \param[in] reorderRefiner Callback to use to send reordered mesh to Refiner
101 : : //! \param[in] ginpoel Mesh connectivity (this chare) using global node IDs
102 : : //! \param[in] coordmap Mesh node coordinates (this chare) for global node IDs
103 : : //! \param[in] el Elements of the mesh chunk we operate on
104 : : //! \param[in] bface Face lists mapped to side set ids
105 : : //! \param[in] triinpoel Interconnectivity of points and boundary-faces
106 : : //! \param[in] bnode Node ids mapped to side set ids
107 : : //! \param[in] nchare Total number of Charm++ worker chares
108 : : // *****************************************************************************
109 : : {
110 : : // Ensure boundary face ids will not index out of face connectivity
111 [ + - ][ - + ]: 247647 : Assert( std::all_of( begin(m_bface), end(m_bface),
[ - - ][ - - ]
[ - - ]
112 : : [&](const auto& s)
113 : : { return std::all_of( begin(s.second), end(s.second),
114 : : [&](auto f){ return f*3+2 < m_triinpoel.size(); } ); } ),
115 : : "Boundary face data structures inconsistent" );
116 : 2459 : }
117 : :
118 : : void
119 : 2459 : Sorter::setup( std::size_t npoin )
120 : : // *****************************************************************************
121 : : // Setup chare mesh boundary node communication map
122 : : //! \param[in] npoin Total number of mesh points in mesh. Note that the number
123 : : //! of mesh points does not have to be exactly the total number of points in
124 : : //! the mesh. It can be a larger number, but not less. This is only used here
125 : : //! to assign nodes to workers that will assign ids to mesh nodes during node
126 : : //! reordering.
127 : : // *****************************************************************************
128 : : {
129 : : // Compute the number of nodes (chunksize) a chare will build a node
130 : : // communication map for. We compute two values of chunksize: one for when
131 : : // the global node ids are abounded between [0...npoin-1], inclusive, and
132 : : // another one for when the global node ids are assigned by a hash algorithm
133 : : // during initial mesh refinement. In the latter case, the maximum
134 : : // representable value of a std::size_t is assumed to be the large global node
135 : : // id and is used to compute the chunksize. To compute the bin id, we attempt
136 : : // to use the first chunksize first: if it gives a chare id that is
137 : : // (strictly) lower than the number of chares, that's good. If not, we compute
138 : : // the bin id based on the second chunksize, which almost always will give a
139 : : // bin id strictly lower than the number of chares, except if the global node
140 : : // id assigned by the hash algorithm in Refiner hits the maximum
141 : : // representable number in std::size_t. If that is the case, we just assign
142 : : // that node to the last chare.
143 : 2459 : auto N = static_cast< std::size_t >( m_nchare );
144 : : std::array< std::size_t, 2 > chunksize{{
145 : 2459 : npoin / N, std::numeric_limits< std::size_t >::max() / N }};
146 : :
147 : : // Find chare-boundary nodes of our mesh chunk. This algorithm collects the
148 : : // global mesh node ids on the chare boundary. A node is on a chare boundary
149 : : // if it belongs to a face of a tetrahedron that has no neighbor tet at a
150 : : // face. The nodes are categorized to bins that will be sent to different
151 : : // chares to build point-to-point communication maps across all chares. The
152 : : // binning is determined by the global node id divided by the chunksizes. See
153 : : // discussion above on how we use two chunksizes for global node ids assigned
154 : : // by the hash algorithm in Refiner (if initial mesh refinement has been
155 : : // done).
156 : 2459 : std::map< int, std::unordered_set< std::size_t > > chbnd;
157 [ + - ]: 2459 : auto el = tk::global2local( m_ginpoel ); // generate local mesh data
158 : 2459 : const auto& inpoel = std::get< 0 >( el ); // local connectivity
159 [ + - ]: 2459 : auto esup = tk::genEsup( inpoel, 4 ); // elements surrounding points
160 [ + - ]: 2459 : auto esuel = tk::genEsuelTet( inpoel, esup ); // elems surrounding elements
161 [ + + ]: 974902 : for (std::size_t e=0; e<esuel.size()/4; ++e) {
162 : 972443 : auto mark = e*4;
163 [ + + ]: 4862215 : for (std::size_t f=0; f<4; ++f)
164 [ + + ]: 3889772 : if (esuel[mark+f] == -1)
165 [ + + ]: 1446064 : for (std::size_t n=0; n<3; ++n) {
166 : 1084548 : auto g = m_ginpoel[ mark+tk::lpofa[f][n] ];
167 : 1084548 : auto bin = g / chunksize[0];
168 [ + + ]: 1084548 : if (bin >= N) bin = g / chunksize[1];
169 [ - + ]: 1084548 : if (bin >= N) bin = N - 1;
170 [ - + ][ - - ]: 1084548 : Assert( bin < N, "Will index out of number of chares" );
[ - - ][ - - ]
171 [ + - ][ + - ]: 1084548 : chbnd[ static_cast< int >( bin ) ].insert( g );
172 : : }
173 : : }
174 : :
175 : : // Send boundary data in bins to chares that will compute communication maps
176 : : // for the data in the bin. These bins form a distributed table. Note that
177 : : // we only send data to those chares that have data to work on. The receiving
178 : : // sides do not know in advance if they receive messages or not. Completion
179 : : // is detected by having the receiver respond back and counting the responses
180 : : // on the sender side, i.e., this chare.
181 : 2459 : m_nbnd = chbnd.size();
182 [ - + ]: 2459 : if (m_nbnd == 0) {
183 [ - - ]: 0 : contribute( sizeof(std::size_t), &m_meshid, CkReduction::nop,
184 : 0 : m_cbs.get< tag::queried >() );
185 : : } else {
186 [ + + ]: 16840 : for (const auto& [ targetchare, bnd ] : chbnd) {
187 [ + - ][ + - ]: 14381 : thisProxy[ targetchare ].query( thisIndex, bnd );
188 : : }
189 : : }
190 : 2459 : }
191 : :
192 : : void
193 : 14381 : Sorter::query( int fromch, const std::unordered_set< std::size_t >& bnd )
194 : : // *****************************************************************************
195 : : // Incoming query for a list of mesh nodes for which this chare compiles node
196 : : // communication maps
197 : : //! \param[in] fromch Sender chare ID
198 : : //! \param[in] bnd Chare-boundary data from another chare
199 : : // *****************************************************************************
200 : : {
201 : : // Store incoming nodes in node->chare and its inverse, chare->node, maps
202 [ + - ][ + - ]: 199214 : for (auto n : bnd) m_nodech[ n ].push_back( fromch );
[ + + ]
203 : 14381 : m_chnode[ fromch ].insert( begin(bnd), end(bnd) );
204 : :
205 : : // Report back to chare message received from
206 [ + - ][ + - ]: 14381 : thisProxy[ fromch ].recvquery();
207 : 14381 : }
208 : :
209 : : void
210 : 14381 : Sorter::recvquery()
211 : : // *****************************************************************************
212 : : // Receive receipt of boundary node lists to query
213 : : // *****************************************************************************
214 : : {
215 [ + + ]: 14381 : if (--m_nbnd == 0)
216 : 2459 : contribute( sizeof(std::size_t), &m_meshid, CkReduction::nop,
217 : 2459 : m_cbs.get< tag::queried >() );
218 : 14381 : }
219 : :
220 : : void
221 : 2459 : Sorter::response()
222 : : // *****************************************************************************
223 : : // Respond to boundary node list queries
224 : : // *****************************************************************************
225 : : {
226 : : std::unordered_map< int, std::map< int, std::unordered_set< std::size_t > > >
227 : 2459 : exp;
228 : :
229 : : // Compute node communication map to be sent back to chares
230 [ + + ]: 16840 : for (const auto& [ neighborchare, bndnodes ] : m_chnode) {
231 [ + - ]: 14381 : auto& nc = exp[ neighborchare ];
232 [ + + ]: 199214 : for (auto n : bndnodes)
233 [ + - ][ + + ]: 512066 : for (auto d : tk::cref_find(m_nodech,n))
234 [ + + ]: 327233 : if (d != neighborchare)
235 [ + - ][ + - ]: 142400 : nc[d].insert( n );
236 : : }
237 : :
238 : : // Send communication maps to chares that issued a query to us. Communication
239 : : // maps were computed above for those chares that queried this map from us.
240 : : // This data form a distributed table and we only work on a chunk of it. Note
241 : : // that we only send data back to those chares that have queried us. The
242 : : // receiving sides do not know in advance if the receive messages or not.
243 : : // Completion is detected by having the receiver respond back and counting
244 : : // the responses on the sender side, i.e., this chare.
245 : 2459 : m_nbnd = exp.size();
246 [ + + ]: 2459 : if (m_nbnd == 0)
247 [ + - ]: 1139 : contribute( sizeof(std::size_t), &m_meshid, CkReduction::nop,
248 : 1139 : m_cbs.get< tag::responded >() );
249 : : else
250 [ + + ]: 15701 : for (const auto& [ targetchare, maps ] : exp)
251 [ + - ][ + - ]: 14381 : thisProxy[ targetchare ].bnd( thisIndex, maps );
252 : 2459 : }
253 : :
254 : : void
255 : 14381 : Sorter::bnd( int fromch,
256 : : const std::map< int, std::unordered_set< std::size_t > >& nodeCommMap )
257 : : // *****************************************************************************
258 : : // Receive boundary node communication maps for our mesh chunk
259 : : //! \param[in] fromch Sender chare ID
260 : : //! \param[in] nodeCommMap Communication map assembled by chare fromch
261 : : // *****************************************************************************
262 : : {
263 [ + + ]: 76535 : for (const auto& [ neighborchare, map ] : nodeCommMap) {
264 [ + - ][ + - ]: 62154 : m_nodeCommMap[ neighborchare ].insert( begin(map), end(map) );
265 : : }
266 : :
267 : : // Report back to chare message received from
268 [ + - ][ + - ]: 14381 : thisProxy[ fromch ].recvbnd();
269 : 14381 : }
270 : :
271 : : void
272 : 14381 : Sorter::recvbnd()
273 : : // *****************************************************************************
274 : : // Receive receipt of boundary node communication map
275 : : // *****************************************************************************
276 : : {
277 [ + + ]: 14381 : if (--m_nbnd == 0)
278 : 1320 : contribute( sizeof(std::size_t), &m_meshid, CkReduction::nop,
279 : 1320 : m_cbs.get< tag::responded >() );
280 : 14381 : }
281 : :
282 : : void
283 : 2459 : Sorter::start()
284 : : // *****************************************************************************
285 : : // Start reordering (if enabled)
286 : : // *****************************************************************************
287 : : {
288 [ + + ]: 2459 : if (g_cfg.get< tag::feedback >()) m_host.chcomm();
289 : :
290 : 2459 : tk::destroy( m_nodech );
291 : 2459 : tk::destroy( m_chnode );
292 : :
293 [ + + ]: 2459 : if (g_cfg.get< tag::reorder >())
294 : 9 : mask(); // continue with mesh node reordering if requested (or required)
295 : : else
296 : 2450 : createDiscWorkers(); // skip mesh node reordering
297 : 2459 : }
298 : :
299 : : void
300 : 9 : Sorter::mask()
301 : : // *****************************************************************************
302 : : // Start preparing for mesh node reordering in parallel
303 : : // *****************************************************************************
304 : : {
305 : : // Compute asymmetric communcation map that will be used for reordering. This
306 : : // communication map is asymmetric because it associates global mesh node IDs
307 : : // to chares only with lower IDs than thisIndex. That is because this chare
308 : : // will need to receive new (reorderd) node IDs only from chares with lower
309 : : // IDs than thisIndex during node reordering. Since it only stores data for
310 : : // lower chare IDs, it is asymmetric. Note that because of this algorithm the
311 : : // type of m_nodeCommMap is an ordered map, because of the std::none_of()
312 : : // algorithm needs to look at ALL chares this chare potentially communicates
313 : : // nodes with that have lower chare IDs that thisIndex. Since the map is
314 : : // ordered, it can walk through from the beginning of m_nodeCommMap until the
315 : : // outer loop variable c, which is the chare ID the outer loop works on in a
316 : : // given cycle.
317 [ + + ]: 65 : for (auto c=m_nodeCommMap.cbegin(); c!=m_nodeCommMap.cend(); ++c) {
318 [ + + ]: 56 : if (thisIndex > c->first) {
319 [ + - ]: 28 : auto& n = m_reordcomm[ c->first ];
320 [ + + ]: 215 : for (auto j : c->second) {
321 [ + - ][ + + ]: 187 : if (std::none_of( m_nodeCommMap.cbegin(), c,
322 : 249 : [j]( const auto& s ) {
323 [ + - ]: 249 : return s.second.find(j) != end(s.second); } ))
324 : : {
325 [ + - ]: 126 : n.insert(j);
326 : : }
327 : : }
328 [ + + ][ + - ]: 28 : if (n.empty()) m_reordcomm.erase( c->first );
329 : : }
330 : : }
331 : :
332 : : // Count up total number of nodes this chare will need to receive
333 : 9 : auto nrecv = tk::sumvalsize( m_reordcomm );
334 : :
335 [ + - ][ + - ]: 9 : if ( g_cfg.get< tag::feedback >() ) m_host.chmask();
336 : :
337 : : // Compute number of mesh node IDs we will assign IDs to
338 : 9 : auto nuniq = m_nodeset.size() - nrecv;
339 : :
340 : : // Start computing offsets for node reordering
341 [ + - ]: 9 : thisProxy.offset( thisIndex, nuniq );
342 : 9 : }
343 : :
344 : : void
345 : 65 : Sorter::offset( int c, std::size_t u )
346 : : // *****************************************************************************
347 : : // Receive number of uniquely assigned global mesh node IDs from chares with
348 : : // lower IDs than thisIndex
349 : : //! \param[in] c Chare ID
350 : : //! \param[in] u Number of mesh node IDs chare c will assign IDs to
351 : : //! \details This function computes the offset each chare will need to start
352 : : //! assigning its new node IDs from. The offset for a chare is the
353 : : //! offset for the previous chare plus the number of node IDs the previous
354 : : //! chare (uniquely) assigns new IDs for minus the number of node IDs the
355 : : //! previous chare receives from others (lower chares). This is computed here
356 : : //! in a parallel/distributed fashion by each chare sending its number of node
357 : : //! IDs (that it uniquely assigns) to all chares. Note that each chare would
358 : : //! only need to send this information to chares with higher IDs, but instead
359 : : //! this function is called in a broadcast fashion, because that is more
360 : : //! efficient than individual calls to only chares with higher IDs. Therefore
361 : : //! when computing the offsets, we only count the lower chares. When this is
362 : : //! done, we have the precise asymmetric communication map as well as the
363 : : //! start offset on all chares and so we can start the distributed global mesh
364 : : //! node ID reordering.
365 : : // *****************************************************************************
366 : : {
367 [ + + ]: 65 : if (c < thisIndex) m_start += u;
368 [ + + ]: 65 : if (++m_noffset == m_nchare) reorder();
369 : 65 : }
370 : :
371 : : void
372 : 9 : Sorter::reorder()
373 : : // *****************************************************************************
374 : : // Reorder global mesh node IDs
375 : : // *****************************************************************************
376 : : {
377 : : // Activate SDAG waits for arriving requests from other chares requesting new
378 : : // node IDs for node IDs we assign new IDs to during reordering; and for
379 : : // computing/receiving lower and upper bounds of global node IDs our chare's
380 : : // linear system will operate on after reordering.
381 [ + - ][ + - ]: 9 : thisProxy[ thisIndex ].wait4prep();
382 : :
383 : : // Send out request for new global node IDs for nodes we do not reorder
384 [ + + ]: 29 : for (const auto& [ targetchare, nodes ] : m_reordcomm)
385 [ + - ][ + - ]: 20 : thisProxy[ targetchare ].request( thisIndex, nodes );
386 : :
387 : : // Lambda to decide if node is assigned a new ID by this chare. If node is not
388 : : // found in the asymmetric communication map, it is owned, i.e., this chare
389 : : // assigns its new id.
390 : 4620 : auto ownnode = [ this ]( std::size_t p ) {
391 : 2310 : return std::all_of( m_reordcomm.cbegin(), m_reordcomm.cend(),
392 : 773 : [&](const auto& s)
393 [ + - ]: 3083 : { return s.second.find(p) == s.second.cend(); } );
394 : 9 : };
395 : :
396 : : // Reorder our chunk of the mesh node IDs. Looping through all of our node
397 : : // IDs, we test if we are to assign a new ID to a node ID, and if so, we
398 : : // assign a new ID, i.e., reorder, by constructing a map associating new to
399 : : // old IDs (m_newnodes). We also count up the reordered nodes, which serves as
400 : : // the new node id. We also store the node coordinates associated to the new
401 : : // node ID.
402 [ + + ]: 2319 : for (auto p : m_nodeset)
403 [ + - ][ + + ]: 2310 : if (ownnode(p)) {
404 [ + - ]: 2184 : m_newnodes[ p ] = m_start; // assign new node ID (reorder)
405 [ + - ][ + - ]: 2184 : m_newcoordmap.emplace( m_start, tk::cref_find(m_coordmap,p) );
406 : 2184 : ++m_start;
407 : : }
408 : :
409 : : // Trigger SDAG wait indicating that reordering our node IDs are complete
410 [ + - ]: 9 : reorderowned_complete();
411 : :
412 : : // If all our nodes have new IDs assigned, reordering complete on this chare
413 [ + + ][ + - ]: 9 : if (m_newnodes.size() == m_nodeset.size()) finish();
414 : 9 : }
415 : :
416 : : void
417 : 20 : Sorter::request( int c, const std::unordered_set< std::size_t >& nd )
418 : : // *****************************************************************************
419 : : // Request new global node IDs for old node IDs
420 : : //! \param[in] c Chare request coming from and to which we send new IDs to
421 : : //! \param[in] nd Set of old node IDs whose new IDs are requested
422 : : // *****************************************************************************
423 : : {
424 : : // Queue up requesting chare and node IDs
425 [ + - ][ + - ]: 20 : m_reqnodes.push_back( { c, nd } );
426 : : // Trigger SDAG wait signaling that node IDs have been requested from us
427 : 20 : nodes_requested_complete();
428 : 20 : }
429 : :
430 : : void
431 : : // cppcheck-suppress unusedFunction
432 : 20 : Sorter::prepare()
433 : : // *****************************************************************************
434 : : // Find new node IDs for old ones and return them to the requestor(s)
435 : : // *****************************************************************************
436 : : {
437 : : // Find and return new node IDs to sender
438 [ + + ]: 40 : for (const auto& [ requestorchare, nodes ] : m_reqnodes) {
439 : : std::unordered_map< std::size_t,
440 : 20 : std::tuple< std::size_t, tk::UnsMesh::Coord > > n;
441 [ + + ]: 146 : for (auto p : nodes) {
442 [ + - ]: 126 : auto newid = tk::cref_find( m_newnodes, p );
443 [ + - ]: 126 : n.emplace( p,
444 [ + - ][ + - ]: 252 : std::make_tuple( newid, tk::cref_find(m_newcoordmap,newid) ) );
445 : : }
446 [ + - ][ + - ]: 20 : thisProxy[ requestorchare ].neworder( n );
447 : 20 : }
448 : :
449 : 20 : tk::destroy( m_reqnodes ); // Clear queue of requests just fulfilled
450 : :
451 : : // Re-enable SDAG wait for preparing new node requests
452 [ + - ][ + - ]: 20 : thisProxy[ thisIndex ].wait4prep();
453 : :
454 : : // Re-enable trigger signaling that reordering of owned node IDs are
455 : : // complete right away
456 : 20 : reorderowned_complete();
457 : 20 : }
458 : :
459 : : void
460 : 20 : Sorter::neworder( const std::unordered_map< std::size_t,
461 : : std::tuple< std::size_t, tk::UnsMesh::Coord > >& nodes )
462 : : // *****************************************************************************
463 : : // Receive new (reordered) global node IDs
464 : : //! \param[in] nodes Map associating new to old node IDs
465 : : // *****************************************************************************
466 : : {
467 : : // Store new node IDs associated to old ones, and node coordinates associated
468 : : // to new node IDs.
469 [ + + ]: 146 : for (const auto& [ oldid, newnodes ] : nodes) {
470 : 126 : auto newid = std::get< 0 >( newnodes );
471 [ + - ]: 126 : m_newnodes[ oldid ] = newid;
472 [ + - ]: 126 : m_newcoordmap.emplace( newid, std::get< 1 >( newnodes ) );
473 : : }
474 : :
475 : : // If all our nodes have new IDs assigned, reorder complete on this PE
476 [ + + ]: 20 : if (m_newnodes.size() == m_nodeset.size()) finish();
477 : 20 : }
478 : :
479 : : void
480 : 9 : Sorter::finish()
481 : : // *****************************************************************************
482 : : // Compute final result of reordering
483 : : //! \details Reordering is now complete on this chare. We now remap all mesh
484 : : //! data to reflect the new ordering.
485 : : // *****************************************************************************
486 : : {
487 : : // Update elem connectivity with the reordered node IDs
488 : 9 : tk::remap( m_ginpoel, m_newnodes );
489 : :
490 : : // Update node coordinate map with the reordered IDs
491 : 9 : m_coordmap = m_newcoordmap;
492 : :
493 : : // Update mesh chunk data structure held in our state with new node order
494 [ + - ]: 9 : m_el = tk::global2local( m_ginpoel );
495 : :
496 : : // Update symmetric chare-node communication map with the reordered IDs
497 [ + + ]: 65 : for (auto& [ neighborchare, map ] : m_nodeCommMap) {
498 : 56 : std::unordered_set< std::size_t > n;
499 [ + - ][ + - ]: 430 : for (auto p : map) n.insert( tk::cref_find( m_newnodes, p ) );
[ + + ]
500 : 56 : map = std::move( n );
501 : 56 : }
502 : :
503 : : // Update boundary face-node connectivity with the reordered node IDs
504 : 9 : tk::remap( m_triinpoel, m_newnodes );
505 : :
506 : : // Update boundary node lists with the reordered node IDs
507 [ + - ][ + + ]: 39 : for (auto& [ setid, nodes ] : m_bnode) tk::remap( nodes, m_newnodes );
508 : :
509 : : // Update mesh in Refiner after reordering
510 : 9 : m_reorderRefiner.send();
511 : :
512 : : // Progress report to host
513 [ + - ]: 9 : if ( g_cfg.get< tag::feedback >() ) m_host.chreordered();
514 : :
515 : 9 : createDiscWorkers();
516 : 9 : }
517 : :
518 : : void
519 : 9 : Sorter::mesh( std::vector< std::size_t >& ginpoel,
520 : : tk::UnsMesh::CoordMap& coordmap,
521 : : std::vector< std::size_t >& triinpoel,
522 : : std::map< int, std::vector< std::size_t > >& bnode )
523 : : // *****************************************************************************
524 : : // Update mesh data we hold for whoever calls this function
525 : : //! \param[in,out] ginpoel Mesh connectivity using global IDs
526 : : //! \param[in,out] coordmap Map of mesh node coordinates
527 : : //! \param[in,out] triinpoel Boundary face-node connectivity
528 : : //! \param[in] bnode Node lists of side sets
529 : : // *****************************************************************************
530 : : {
531 : 9 : ginpoel = m_ginpoel;
532 : 9 : coordmap = m_coordmap;
533 : 9 : triinpoel = m_triinpoel;
534 : 9 : bnode = m_bnode;
535 : 9 : }
536 : :
537 : : void
538 : 2459 : Sorter::createDiscWorkers()
539 : : // *****************************************************************************
540 : : // Create Discretization chare array elements on this PE
541 : : //! \details We create chare array elements by calling the insert() member
542 : : //! function, which allows specifying the PE on which the array element is
543 : : //! created. and we send each chare array element the chunk of mesh it will
544 : : //! operate on.
545 : : // *****************************************************************************
546 : : {
547 : : // Create worker array element using Charm++ dynamic chare array element
548 : : // insertion.
549 : :
550 [ + - ]: 4918 : m_discretization[ thisIndex ].insert( m_meshid, m_host,
551 [ + - ]: 2459 : m_meshwriter, m_coordmap, m_el, m_nodeCommMap, m_nchare );
552 : :
553 : 2459 : contribute( sizeof(std::size_t), &m_meshid, CkReduction::nop,
554 : 2459 : m_cbs.get< tag::discinserted >() );
555 : 2459 : }
556 : :
557 : : void
558 : 2459 : Sorter::createWorkers()
559 : : // *****************************************************************************
560 : : // Create worker chare array element
561 : : // *****************************************************************************
562 : : {
563 : : // Make sure (bound) base is already created and accessible
564 [ + - ][ + - ]: 2459 : Assert( m_discretization[ thisIndex ].ckLocal() != nullptr,
[ - + ][ - - ]
[ - - ][ - - ]
565 : : "About to pass nullptr" );
566 : :
567 : : // Create worker array element using Charm++ dynamic chare array element
568 : : // insertion.
569 : 2459 : const auto& solver = g_cfg.get< tag::solver >();
570 [ + + ]: 2459 : if (solver == "riecg") {
571 [ + - ]: 1118 : m_riecg[ thisIndex ].insert( m_discretization, m_bface, m_bnode,
572 [ + - ]: 559 : m_triinpoel );
573 : : }
574 [ + + ]: 1900 : else if (solver == "laxcg") {
575 [ + - ]: 632 : m_laxcg[ thisIndex ].insert( m_discretization, m_bface, m_bnode,
576 [ + - ]: 316 : m_triinpoel );
577 : : }
578 [ + + ]: 1584 : else if (solver == "zalcg") {
579 [ + - ]: 870 : m_zalcg[ thisIndex ].insert( m_discretization, m_bface, m_bnode,
580 [ + - ]: 435 : m_triinpoel );
581 : : }
582 [ + + ]: 1149 : else if (solver == "kozcg") {
583 [ + - ]: 918 : m_kozcg[ thisIndex ].insert( m_discretization, m_bface, m_bnode,
584 [ + - ]: 459 : m_triinpoel );
585 : : }
586 [ + + ]: 690 : else if (solver == "chocg") {
587 [ + - ]: 730 : m_chocg[ thisIndex ].insert( m_discretization, m_cgpre, m_cgmom, m_bface,
588 [ + - ]: 365 : m_bnode, m_triinpoel );
589 : : }
590 [ + - ]: 325 : else if (solver == "lohcg") {
591 [ + - ]: 650 : m_lohcg[ thisIndex ].insert( m_discretization, m_cgpre, m_bface, m_bnode,
592 [ + - ]: 325 : m_triinpoel );
593 : : }
594 : : else {
595 [ - - ][ - - ]: 0 : Throw( "Unknown solver: " + solver );
[ - - ]
596 : : }
597 : :
598 [ + + ]: 2459 : if ( g_cfg.get< tag::feedback >() ) m_host.chcreated();
599 : :
600 : 2459 : contribute( sizeof(std::size_t), &m_meshid, CkReduction::nop,
601 : 2459 : m_cbs.get< tag::workinserted >() );
602 : :
603 : : // Free up some memory
604 : 2459 : tk::destroy( m_ginpoel );
605 : 2459 : tk::destroy( m_coordmap );
606 : 2459 : tk::destroy( m_bface );
607 : 2459 : tk::destroy( m_triinpoel );
608 : 2459 : tk::destroy( m_bnode );
609 : 2459 : tk::destroy( m_nodeset );
610 : 2459 : tk::destroy( m_nodech );
611 : 2459 : tk::destroy( m_chnode );
612 : 2459 : tk::destroy( m_nodeCommMap );
613 : 2459 : tk::destroy( m_reordcomm );
614 : 2459 : tk::destroy( m_newnodes );
615 : 2459 : tk::destroy( m_reqnodes );
616 : 2459 : }
617 : :
618 : : #include "NoWarning/sorter.def.h"
|
