Branch data Line data Source code
1 : : // *****************************************************************************
2 : : /*!
3 : : \file src/Inciter/Partitioner.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 Charm++ chare partitioner nodegroup used to perform mesh
9 : : partitioning
10 : : \details Charm++ chare partitioner nodegroup used to perform mesh read and
11 : : partitioning, one worker per compute node.
12 : : */
13 : : // *****************************************************************************
14 : :
15 : : #include <numeric>
16 : :
17 : : #include "NoWarning/kokkos_cr.hpp"
18 : :
19 : : #include "Partitioner.hpp"
20 : : #include "DerivedData.hpp"
21 : : #include "Reorder.hpp"
22 : : #include "MeshReader.hpp"
23 : : #include "CGPDE.hpp"
24 : : #include "DGPDE.hpp"
25 : : #include "Inciter/Options/Scheme.hpp"
26 : : #include "UnsMesh.hpp"
27 : : #include "ContainerUtil.hpp"
28 : : #include "Callback.hpp"
29 : :
30 : : namespace inciter {
31 : :
32 : : extern ctr::InputDeck g_inputdeck;
33 : :
34 : : } // inciter::
35 : :
36 : : using inciter::Partitioner;
37 : :
38 : 574 : Partitioner::Partitioner(
39 : : std::size_t meshid,
40 : : const std::string& filename,
41 : : const tk::PartitionerCallback& cbp,
42 : : const tk::RefinerCallback& cbr,
43 : : const tk::SorterCallback& cbs,
44 : : const CProxy_Transporter& host,
45 : : const CProxy_Refiner& refiner,
46 : : const CProxy_Sorter& sorter,
47 : : const tk::CProxy_MeshWriter& meshwriter,
48 : : const std::vector< Scheme >& scheme,
49 : : const std::map< int, std::vector< std::size_t > >& bface,
50 : : const std::map< int, std::vector< std::size_t > >& faces,
51 : 574 : const std::map< int, std::vector< std::size_t > >& bnode ) :
52 : : m_meshid( meshid ),
53 : : m_cbp( cbp ),
54 : : m_cbr( cbr ),
55 : : m_cbs( cbs ),
56 : : m_host( host ),
57 : : m_refiner( refiner ),
58 : : m_sorter( sorter ),
59 : : m_meshwriter( meshwriter ),
60 : : m_scheme( scheme ),
61 : : m_ginpoel(),
62 : : m_coord(),
63 : : m_inpoel(),
64 : : m_lid(),
65 : : m_elemBlockId(),
66 : : m_ndist( 0 ),
67 : : m_nchare( 0 ),
68 : : m_nface(),
69 : : m_chinpoel(),
70 : : m_chcoordmap(),
71 : : m_chbface(),
72 : : m_chtriinpoel(),
73 : : m_chbnode(),
74 : : m_chelemblockid(),
75 : : m_bface( bface ),
76 [ + - ][ + - ]: 574 : m_bnode( bnode )
[ + - ][ + - ]
[ + - ][ + - ]
[ + - ]
77 : : // *****************************************************************************
78 : : // Constructor
79 : : //! \param[in] meshid Mesh ID
80 : : //! \param[in] filename Input mesh filename to read from
81 : : //! \param[in] cbp Charm++ callbacks for Partitioner
82 : : //! \param[in] cbr Charm++ callbacks for Refiner
83 : : //! \param[in] cbs Charm++ callbacks for Sorter
84 : : //! \param[in] host Host Charm++ proxy we are being called from
85 : : //! \param[in] refiner Mesh refiner proxy
86 : : //! \param[in] sorter Mesh reordering (sorter) proxy
87 : : //! \param[in] meshwriter Mesh writer proxy
88 : : //! \param[in] scheme Discretization scheme
89 : : //! \param[in] bface File-internal elem ids of side sets (whole mesh)
90 : : //! \param[in] faces Elem-relative face ids of side sets (whole mesh)
91 : : //! \param[in] bnode Node lists of side sets (whole mesh)
92 : : // *****************************************************************************
93 : : {
94 : : // The following call has to be made on all MPI ranks immediately after
95 : : // initializing MPI. This has to be done from a Charm++ nodegroup, since there
96 : : // exists only one rank per node on a nodegroup, in SMP mode.
97 : : // see https://github.com/trilinos/Trilinos/issues/11197#issuecomment-1301325163
98 [ + - ]: 574 : Kokkos::initialize();
99 : :
100 : : // Create mesh reader
101 [ + - ]: 1148 : tk::MeshReader mr( filename );
102 : :
103 : : // Read this compute node's chunk of the mesh (graph and coords) from file
104 : 1148 : std::vector< std::size_t > triinpoel;
105 : 574 : mr.readMeshPart( m_ginpoel, m_inpoel, triinpoel, m_lid, m_coord,
106 [ + - ]: 574 : m_elemBlockId, CkNumNodes(), CkMyNode() );
107 : :
108 : : // Compute triangle connectivity for side sets, reduce boundary face for side
109 : : // sets to this compute node only and to compute-node-local face ids
110 [ + - ]: 574 : m_triinpoel = mr.triinpoel( m_bface, faces, m_ginpoel, triinpoel );
111 : :
112 : : // Reduce boundary node lists (global ids) for side sets to this compute node
113 : : // only
114 [ + - ]: 574 : ownBndNodes( m_lid, m_bnode );
115 : :
116 : : // Sum number of cells across distributed mesh
117 [ + - ]: 1148 : std::vector< std::size_t > meshdata{ meshid, m_ginpoel.size()/4 };
118 [ + - ]: 574 : contribute( meshdata, CkReduction::sum_ulong, m_cbp.get< tag::load >() );
119 : 574 : }
120 : :
121 : : void
122 : 574 : Partitioner::ownBndNodes(
123 : : const std::unordered_map< std::size_t, std::size_t >& lid,
124 : : std::map< int, std::vector< std::size_t > >& bnode )
125 : : // *****************************************************************************
126 : : // Keep only those nodes for side sets that reside on this compute node
127 : : //! \param[in] lid Global->local node IDs of elements of this compute node's
128 : : //! mesh chunk
129 : : //! \param[in,out] bnode Global ids of nodes for side sets for whole mesh
130 : : //! \details This function overwrites the input boundary node lists map with the
131 : : //! nodes that reside on the caller compute node.
132 : : // *****************************************************************************
133 : : {
134 : 1148 : std::map< int, std::vector< std::size_t > > bnode_own;
135 : :
136 [ + + ]: 1007 : for (const auto& [ setid, nodes ] : bnode) {
137 [ + - ]: 433 : auto& b = bnode_own[ setid ];
138 [ + + ]: 95037 : for (auto n : nodes) {
139 [ + - ]: 94604 : auto i = lid.find( n );
140 [ + + ][ + - ]: 94604 : if (i != end(lid)) b.push_back( n );
141 : : }
142 [ + + ][ + - ]: 433 : if (b.empty()) bnode_own.erase( setid );
143 : : }
144 : :
145 : 574 : bnode = std::move(bnode_own);
146 : 574 : }
147 : :
148 : : void
149 : 574 : Partitioner::partition( int nchare )
150 : : // *****************************************************************************
151 : : // Partition the computational mesh into a number of chares
152 : : //! \param[in] nchare Number of parts the mesh will be partitioned into
153 : : //! \details This function calls the mesh partitioner to partition the mesh. The
154 : : //! number of partitions equals the number nchare argument which must be no
155 : : //! lower than the number of compute nodes.
156 : : // *****************************************************************************
157 : : {
158 [ - + ][ - - ]: 574 : Assert( nchare >= CkNumNodes(), "Number of chares must not be lower than the "
[ - - ][ - - ]
159 : : "number of compute nodes" );
160 : :
161 : : // Generate element IDs for Zoltan
162 [ + - ]: 1148 : std::vector< long long > gelemid( m_ginpoel.size()/4 );
163 : 574 : std::iota( begin(gelemid), end(gelemid), 0 );
164 : :
165 : 574 : m_nchare = nchare;
166 : 574 : const auto alg = g_inputdeck.get< tag::partitioning >();
167 : : const auto che = tk::zoltan::geomPartMesh( alg,
168 [ + - ]: 574 : centroids( m_inpoel, m_coord ),
169 : : gelemid,
170 [ + - ]: 574 : nchare );
171 : :
172 [ - + ][ - - ]: 574 : if ( g_inputdeck.get< tag::cmd, tag::feedback >() ) m_host.pepartitioned();
173 : :
174 [ + - ]: 574 : contribute( sizeof(std::size_t), &m_meshid, CkReduction::nop,
175 : 574 : m_cbp.get< tag::partitioned >() );
176 : :
177 [ - + ][ - - ]: 574 : Assert( che.size() == gelemid.size(), "Size of ownership array (chare ID "
[ - - ][ - - ]
178 : : "of elements) after mesh partitioning does not equal the number of "
179 : : "mesh graph elements" );
180 : :
181 : : // Categorize mesh elements (given by their gobal node IDs) by target chare
182 : : // and distribute to their compute nodes based on mesh partitioning.
183 [ + - ][ + - ]: 574 : distribute( categorize( che ) );
184 : 574 : }
185 : :
186 : : void
187 : 2128 : Partitioner::addMesh(
188 : : int fromnode,
189 : : const std::unordered_map< int, // chare id
190 : : std::tuple<
191 : : std::vector< std::size_t >, // tet connectivity
192 : : tk::UnsMesh::CoordMap, // node coords
193 : : std::unordered_map< int, std::vector< std::size_t > >, // bface conn
194 : : std::unordered_map< int, std::vector< std::size_t > >, // bnodes
195 : : std::vector< std::size_t > // elem-blocks
196 : : > >& chmesh )
197 : : // *****************************************************************************
198 : : // Receive mesh associated to chares we own after refinement
199 : : //! \param[in] fromnode Compute node call coming from
200 : : //! \param[in] chmesh Map associating mesh connectivities to global node ids
201 : : //! and node coordinates for mesh chunks we are assigned by the partitioner
202 : : // *****************************************************************************
203 : : {
204 : : // Store mesh connectivity and global node coordinates categorized by chares.
205 : : // The send side also writes to the data written here, so concat.
206 [ + + ]: 8462 : for (const auto& [ chareid, chunk ] : chmesh) {
207 [ + - ][ - + ]: 6334 : Assert( node(chareid) == CkMyNode(), "Compute node "
[ - - ][ - - ]
[ - - ][ - - ]
[ - - ]
208 : : + std::to_string(CkMyNode()) +
209 : : " received a mesh whose chare it does not own" );
210 : : // Store domain element (tetrahedron) connectivity
211 : 6334 : const auto& inpoel = std::get< 0 >( chunk );
212 [ + - ]: 6334 : auto& inp = m_chinpoel[ chareid ]; // will store tetrahedron connectivity
213 [ + - ]: 6334 : inp.insert( end(inp), begin(inpoel), end(inpoel) );
214 : : // Store own mesh block id and associated tet ids
215 : 6334 : const auto& elblock = std::get<4>( chunk );
216 [ + - ]: 6334 : auto& cheb = m_chelemblockid[ chareid ];
217 [ + - ]: 6334 : cheb.insert( end(cheb), begin(elblock), end(elblock) );
218 : : // Store mesh node coordinates associated to global node IDs
219 : 6334 : const auto& coord = std::get< 1 >( chunk );
220 [ + - ][ - + ]: 6334 : Assert( tk::uniquecopy(inpoel).size() == coord.size(), "Size mismatch" );
[ - - ][ - - ]
[ - - ]
221 [ + - ]: 6334 : auto& chcm = m_chcoordmap[ chareid ]; // will store node coordinates
222 [ + - ]: 6334 : chcm.insert( begin(coord), end(coord) ); // concatenate node coords
223 : : // Store boundary side set id + face ids + face connectivities
224 : 6334 : const auto& bconn = std::get< 2 >( chunk );
225 [ + - ]: 6334 : auto& bface = m_chbface[ chareid ]; // for side set id + boundary face ids
226 [ + - ]: 6334 : auto& t = m_chtriinpoel[ chareid ]; // for boundary face connectivity
227 [ + - ]: 6334 : auto& f = m_nface[ chareid ]; // use counter for chare
228 [ + + ]: 12447 : for (const auto& [ setid, faceids ] : bconn) {
229 [ + - ]: 6113 : auto& b = bface[ setid ];
230 [ + + ]: 78014 : for (std::size_t i=0; i<faceids.size()/3; ++i) {
231 [ + - ]: 71901 : b.push_back( f++ );
232 [ + - ]: 71901 : t.push_back( faceids[i*3+0] );
233 [ + - ]: 71901 : t.push_back( faceids[i*3+1] );
234 [ + - ]: 71901 : t.push_back( faceids[i*3+2] );
235 : : }
236 : : }
237 : : // Store boundary side set id + node lists
238 : 6334 : const auto& bnode = std::get< 3 >( chunk );
239 [ + - ]: 6334 : auto& nodes = m_chbnode[ chareid ]; // for side set id + boundary nodes
240 [ + + ]: 7601 : for (const auto& [ setid, bnodes ] : bnode) {
241 [ + - ]: 1267 : auto& b = nodes[ setid ];
242 [ + - ]: 1267 : b.insert( end(b), begin(bnodes), end(bnodes) );
243 : : }
244 : : }
245 : :
246 [ + - ]: 2128 : thisProxy[ fromnode ].recvMesh();
247 : 2128 : }
248 : :
249 : : int
250 : 12668 : Partitioner::node( int id ) const
251 : : // *****************************************************************************
252 : : // Return nodegroup id for chare id
253 : : //! \param[in] id Chare id
254 : : //! \return Nodegroup that creates the chare
255 : : //! \details This is computed based on a simple contiguous linear
256 : : //! distribution of chare ids to compute nodes.
257 : : // *****************************************************************************
258 : : {
259 [ - + ][ - - ]: 12668 : Assert( m_nchare > 0, "Number of chares must be a positive number" );
[ - - ][ - - ]
260 : 12668 : auto p = id / (m_nchare / CkNumNodes());
261 [ + + ]: 12668 : if (p >= CkNumNodes()) p = CkNumNodes()-1;
262 [ - + ][ - - ]: 12668 : Assert( p < CkNumNodes(), "Assigning to nonexistent node" );
[ - - ][ - - ]
263 : 12668 : return p;
264 : : }
265 : :
266 : : void
267 : 2128 : Partitioner::recvMesh()
268 : : // *****************************************************************************
269 : : // Acknowledge received mesh chunk and its nodes after mesh refinement
270 : : // *****************************************************************************
271 : : {
272 [ + + ]: 2128 : if (--m_ndist == 0) {
273 [ - + ]: 494 : if (g_inputdeck.get< tag::cmd, tag::feedback >()) m_host.pedistributed();
274 : 494 : contribute( sizeof(std::size_t), &m_meshid, CkReduction::nop,
275 : 494 : m_cbp.get< tag::distributed >() );
276 : : }
277 : 2128 : }
278 : :
279 : : void
280 : 574 : Partitioner::refine()
281 : : // *****************************************************************************
282 : : // Optionally start refining the mesh
283 : : // *****************************************************************************
284 : : {
285 [ + - ]: 574 : auto dist = distribution( m_nchare );
286 : :
287 : 574 : std::size_t error = 0;
288 [ - + ]: 574 : if (m_chinpoel.size() < static_cast<std::size_t>(dist[1])) {
289 : :
290 : 0 : error = 1;
291 : :
292 : : } else {
293 : :
294 [ + + ]: 2517 : for (int c=0; c<dist[1]; ++c) {
295 : : // compute chare ID
296 : 1943 : auto cid = CkMyNode() * dist[0] + c;
297 : :
298 [ + - ][ + - ]: 1943 : Assert(tk::cref_find(m_chinpoel,cid).size()/4 ==
[ - + ][ - - ]
[ - - ][ - - ]
299 : : tk::cref_find(m_chelemblockid,cid).size(),
300 : : "incorrect number of elements in elem-block id vector");
301 : : // create refiner Charm++ chare array element using dynamic insertion
302 [ + - ]: 3886 : m_refiner[ cid ].insert( m_meshid,
303 : 1943 : m_host,
304 : 1943 : m_sorter,
305 : 1943 : m_meshwriter,
306 : 1943 : m_scheme,
307 : 1943 : m_cbr,
308 : 1943 : m_cbs,
309 [ + - ]: 1943 : tk::cref_find(m_chinpoel,cid),
310 [ + - ]: 1943 : tk::cref_find(m_chcoordmap,cid),
311 [ + - ]: 1943 : tk::cref_find(m_chbface,cid),
312 [ + - ]: 1943 : tk::cref_find(m_chtriinpoel,cid),
313 [ + - ]: 1943 : tk::cref_find(m_chbnode,cid),
314 [ + - ][ + - ]: 1943 : tk::cref_find(m_chelemblockid,cid),
315 : : m_nchare );
316 : : }
317 : :
318 : : }
319 : :
320 : 574 : tk::destroy( m_ginpoel );
321 : 574 : tk::destroy( m_coord );
322 : 574 : tk::destroy( m_inpoel );
323 : 574 : tk::destroy( m_lid );
324 : 574 : tk::destroy( m_elemBlockId );
325 : 574 : tk::destroy( m_nface );
326 : 574 : tk::destroy( m_nodech );
327 : 574 : tk::destroy( m_linnodes );
328 : 574 : tk::destroy( m_chinpoel );
329 : 574 : tk::destroy( m_chcoordmap );
330 : 574 : tk::destroy( m_chbface );
331 : 574 : tk::destroy( m_chtriinpoel );
332 : 574 : tk::destroy( m_chbnode );
333 : 574 : tk::destroy( m_chelemblockid );
334 : 574 : tk::destroy( m_bnodechares );
335 : 574 : tk::destroy( m_bface );
336 : 574 : tk::destroy( m_triinpoel );
337 : 574 : tk::destroy( m_bnode );
338 : :
339 [ + - ]: 1148 : std::vector< std::size_t > meshdata{ m_meshid, error };
340 [ + - ]: 574 : contribute( meshdata, CkReduction::max_ulong, m_cbp.get<tag::refinserted>() );
341 : 574 : }
342 : :
343 : : std::array< std::vector< tk::real >, 3 >
344 : 574 : Partitioner::centroids( const std::vector< std::size_t >& inpoel,
345 : : const tk::UnsMesh::Coords& coord )
346 : : // *****************************************************************************
347 : : // Compute element centroid coordinates
348 : : //! \param[in] inpoel Mesh connectivity with local ids
349 : : //! \param[in] coord Node coordinates
350 : : //! \return Centroids for all cells on this compute node
351 : : // *****************************************************************************
352 : : {
353 [ - + ][ - - ]: 574 : Assert( tk::uniquecopy(inpoel).size() == coord[0].size(), "Size mismatch" );
[ - - ][ - - ]
354 : :
355 : 574 : const auto& x = coord[0];
356 : 574 : const auto& y = coord[1];
357 : 574 : const auto& z = coord[2];
358 : :
359 : : // Make room for element centroid coordinates
360 : 574 : std::array< std::vector< tk::real >, 3 > cent;
361 : 574 : auto& cx = cent[0];
362 : 574 : auto& cy = cent[1];
363 : 574 : auto& cz = cent[2];
364 : 574 : auto num = inpoel.size()/4;
365 [ + - ]: 574 : cx.resize( num );
366 [ + - ]: 574 : cy.resize( num );
367 [ + - ]: 574 : cz.resize( num );
368 : :
369 : : // Compute element centroids for mesh passed in
370 [ + + ]: 447131 : for (std::size_t e=0; e<num; ++e) {
371 : 446557 : auto A = inpoel[e*4+0];
372 : 446557 : auto B = inpoel[e*4+1];
373 : 446557 : auto C = inpoel[e*4+2];
374 : 446557 : auto D = inpoel[e*4+3];
375 : 446557 : cx[e] = (x[A] + x[B] + x[C] + x[D]) / 4.0;
376 : 446557 : cy[e] = (y[A] + y[B] + y[C] + y[D]) / 4.0;
377 : 446557 : cz[e] = (z[A] + z[B] + z[C] + z[D]) / 4.0;
378 : : }
379 : :
380 : 574 : return cent;
381 : : }
382 : :
383 : : std::unordered_map< int, Partitioner::MeshData >
384 : 574 : Partitioner::categorize( const std::vector< std::size_t >& target ) const
385 : : // *****************************************************************************
386 : : // Categorize mesh data by target
387 : : //! \param[in] target Target chares of mesh elements, size: number of
388 : : //! elements in the chunk of the mesh graph on this compute node.
389 : : //! \return Vector of global mesh node ids connecting elements owned by each
390 : : //! target chare.
391 : : // *****************************************************************************
392 : : {
393 [ - + ][ - - ]: 574 : Assert( target.size() == m_ginpoel.size()/4, "Size mismatch");
[ - - ][ - - ]
394 : :
395 : : using Face = tk::UnsMesh::Face;
396 : :
397 : : // Build hash map associating side set id to boundary faces
398 : : std::unordered_map< Face, int,
399 : 1148 : tk::UnsMesh::Hash<3>, tk::UnsMesh::Eq<3> > faceside;
400 [ + + ]: 2508 : for (const auto& [ setid, faceids ] : m_bface)
401 [ + + ]: 164242 : for (auto f : faceids)
402 : 162308 : faceside[ {{ m_triinpoel[f*3+0],
403 : 162308 : m_triinpoel[f*3+1],
404 [ + - ]: 162308 : m_triinpoel[f*3+2] }} ] = setid;
405 : :
406 : : // Build hash map associating side set ids to boundary nodes
407 : 1148 : std::unordered_map< std::size_t, std::unordered_set< int > > nodeside;
408 [ + + ]: 1001 : for (const auto& [ setid, nodes ] : m_bnode)
409 [ + + ]: 74070 : for (auto n : nodes)
410 [ + - ][ + - ]: 73643 : nodeside[ n ].insert( setid );
411 : :
412 : : // Categorize mesh data (tets, node coordinates, and boundary data) by target
413 : : // chare based on which chare the partitioner assigned elements (tets) to
414 : 574 : std::unordered_map< int, MeshData > chmesh;
415 [ + + ]: 447131 : for (std::size_t e=0; e<target.size(); ++e) {
416 : : // Construct a tetrahedron with global node ids
417 : 446557 : tk::UnsMesh::Tet t{{ m_ginpoel[e*4+0], m_ginpoel[e*4+1],
418 : 446557 : m_ginpoel[e*4+2], m_ginpoel[e*4+3] }};
419 : : // Categorize tetrahedron (domain element) connectivity
420 [ + - ]: 446557 : auto& mesh = chmesh[ static_cast<int>(target[e]) ];
421 : 446557 : auto& inpoel = std::get< 0 >( mesh );
422 [ + - ]: 446557 : inpoel.insert( end(inpoel), begin(t), end(t) );
423 : : // Categorize mesh block ids and associated local tet ids (local tet ids
424 : : // basically correspond to the inpoel entries updated in line above).
425 : 446557 : [[maybe_unused]] bool added(false);
426 : 446557 : auto& elblock = std::get< 3 >( mesh );
427 [ + + ]: 906981 : for (const auto& [blid, elset] : m_elemBlockId) {
428 [ + - ][ + + ]: 460424 : if (elset.find(e) != elset.end()) {
429 [ + - ]: 446557 : elblock.push_back(blid);
430 : 446557 : added = true;
431 : : }
432 : : }
433 [ - + ][ - - ]: 446557 : Assert(added, "Tet element not found in any block after partitioning");
[ - - ][ - - ]
434 : : // Categorize boundary face connectivity
435 : 446557 : auto& bconn = std::get< 1 >( mesh );
436 : 446557 : std::array<Face,4> face{{ {{t[0],t[2],t[1]}}, {{t[0],t[1],t[3]}},
437 : 446557 : {{t[0],t[3],t[2]}}, {{t[1],t[2],t[3]}} }};
438 [ + + ]: 2232785 : for (const auto& f : face) {
439 [ + - ]: 1786228 : auto it = faceside.find( f );
440 [ + + ]: 1786228 : if (it != end(faceside)) {
441 [ + - ]: 162308 : auto& s = bconn[ it->second ];
442 [ + - ]: 162308 : s.insert( end(s), begin(f), end(f) );
443 : : }
444 : : }
445 : : // Categorize boundary node lists
446 : 446557 : auto& bnode = std::get< 2 >( mesh );
447 [ + + ]: 2232785 : for (const auto& n : t) {
448 [ + - ]: 1786228 : auto it = nodeside.find( n );
449 [ + + ]: 1786228 : if (it != end(nodeside))
450 [ + + ]: 649128 : for (auto s : it->second)
451 [ + - ][ + - ]: 334712 : bnode[ s ].push_back( n );
452 : : }
453 : : }
454 : :
455 : : // Make boundary node lists unique per side set
456 [ + + ]: 8639 : for (auto& c : chmesh)
457 [ + + ]: 9843 : for (auto& n : std::get<2>(c.second))
458 [ + - ]: 1778 : tk::unique( n.second );
459 : :
460 : : // Make sure all compute nodes have target chares assigned
461 [ - + ][ - - ]: 574 : Assert( !chmesh.empty(), "No elements have been assigned to a chare" );
[ - - ][ - - ]
462 : :
463 : : // This check should always be done, hence ErrChk and not Assert, as it
464 : : // can result from particular pathological combinations of (1) too large
465 : : // degree of virtualization, (2) too many compute nodes, and/or (3) too small
466 : : // of a mesh and not due to programmer error.
467 [ + + ]: 8639 : for(const auto& c : chmesh)
468 [ - + ][ - - ]: 8065 : ErrChk( !std::get<0>(c.second).empty(),
[ - - ][ - - ]
469 : : "Overdecomposition of the mesh is too large compared to the "
470 : : "number of work units computed based on the degree of "
471 : : "virtualization desired. As a result, there would be at least "
472 : : "one work unit with no mesh elements to work on, i.e., nothing "
473 : : "to do. Solution 1: decrease the virtualization to a lower "
474 : : "value using the command-line argument '-u'. Solution 2: "
475 : : "decrease the number processing elements (PEs and/or compute "
476 : : "nodes) using the charmrun command-line argument '+pN' where N is "
477 : : "the number of PEs (or in SMP-mode in combination with +ppn to "
478 : : "reduce the number of compute nodes), which implicitly increases "
479 : : "the size (and thus decreases the number) of work units.)" );
480 : :
481 : 1148 : return chmesh;
482 : : }
483 : :
484 : : tk::UnsMesh::CoordMap
485 : 8065 : Partitioner::coordmap( const std::vector< std::size_t >& inpoel )
486 : : // *****************************************************************************
487 : : // Extract coordinates associated to global nodes of a mesh chunk
488 : : //! \param[in] inpoel Mesh connectivity
489 : : //! \return Map storing the coordinates of unique nodes associated to global
490 : : //! node IDs in mesh given by inpoel
491 : : // *****************************************************************************
492 : : {
493 [ - + ][ - - ]: 8065 : Assert( inpoel.size() % 4 == 0, "Incomplete mesh connectivity" );
[ - - ][ - - ]
494 : :
495 : 8065 : tk::UnsMesh::CoordMap map;
496 : :
497 [ + - ][ + + ]: 362419 : for (auto g : tk::uniquecopy(inpoel)) {
498 [ + - ]: 354354 : auto i = tk::cref_find( m_lid, g );
499 [ + - ]: 354354 : auto& c = map[g];
500 : 354354 : c[0] = m_coord[0][i];
501 : 354354 : c[1] = m_coord[1][i];
502 : 354354 : c[2] = m_coord[2][i];
503 : : }
504 : :
505 [ + - ][ - + ]: 8065 : Assert( tk::uniquecopy(inpoel).size() == map.size(), "Size mismatch" );
[ - - ][ - - ]
[ - - ]
506 : :
507 : 8065 : return map;
508 : : }
509 : :
510 : : void
511 : 574 : Partitioner::distribute( std::unordered_map< int, MeshData >&& mesh )
512 : : // *****************************************************************************
513 : : // Distribute mesh to target compute nodes after mesh partitioning
514 : : //! \param[in] mesh Mesh data categorized by target by target chares
515 : : // *****************************************************************************
516 : : {
517 [ + - ]: 574 : auto dist = distribution( m_nchare );
518 : :
519 : : // Extract mesh data whose chares are on ("owned by") this compute node
520 [ + + ]: 2517 : for (int c=0; c<dist[1]; ++c) {
521 : 1943 : auto chid = CkMyNode() * dist[0] + c; // compute owned chare ID
522 [ + - ]: 1943 : const auto it = mesh.find( chid ); // attempt to find its mesh data
523 [ + + ]: 1943 : if (it != end(mesh)) { // if found
524 : : // Store own tetrahedron connectivity
525 : 1731 : const auto& inpoel = std::get<0>( it->second );
526 [ + - ]: 1731 : auto& inp = m_chinpoel[ chid ]; // will store own mesh connectivity
527 [ + - ]: 1731 : inp.insert( end(inp), begin(inpoel), end(inpoel) );
528 : : // Store own mesh block id and associated tet ids
529 : 1731 : const auto& elblock = std::get<3>( it->second );
530 [ + - ]: 1731 : auto& cheb = m_chelemblockid[ chid ];
531 [ + - ]: 1731 : cheb.insert( end(cheb), begin(elblock), end(elblock) );
532 : : // Store own node coordinates
533 [ + - ]: 1731 : auto& chcm = m_chcoordmap[ chid ]; // will store own node coordinates
534 [ + - ]: 3462 : auto cm = coordmap( inpoel ); // extract node coordinates
535 [ + - ]: 1731 : chcm.insert( begin(cm), end(cm) ); // concatenate node coords
536 : : // Store own boundary face connectivity
537 : 1731 : const auto& bconn = std::get<1>( it->second );
538 [ + - ]: 1731 : auto& bface = m_chbface[ chid ]; // will store own boundary faces
539 [ + - ]: 1731 : auto& t = m_chtriinpoel[ chid ]; // wil store own boundary face conn
540 [ + - ]: 1731 : auto& f = m_nface[ chid ]; // use counter for chare
541 [ + + ]: 4004 : for (const auto& [ setid, faceids ] : bconn) {
542 [ + - ]: 2273 : auto& b = bface[ setid ];
543 [ + + ]: 92680 : for (std::size_t i=0; i<faceids.size()/3; ++i) {
544 [ + - ]: 90407 : b.push_back( f++ );
545 [ + - ]: 90407 : t.push_back( faceids[i*3+0] );
546 [ + - ]: 90407 : t.push_back( faceids[i*3+1] );
547 [ + - ]: 90407 : t.push_back( faceids[i*3+2] );
548 : : }
549 : : }
550 : : // Store own boundary node lists
551 : 1731 : const auto& bnode = std::get<2>( it->second );
552 [ + - ]: 1731 : auto& nodes = m_chbnode[ chid ]; // will store own boundary nodes
553 [ + + ]: 2242 : for (const auto& [ setid, nodeids ] : bnode) {
554 [ + - ]: 511 : auto& b = nodes[ setid ];
555 [ + - ]: 511 : b.insert( end(b), begin(nodeids), end(nodeids) );
556 : : }
557 : : // Remove chare ID and mesh data
558 [ + - ]: 1731 : mesh.erase( it );
559 : : }
560 [ + - ][ - + ]: 1943 : Assert( mesh.find(chid) == end(mesh), "Not all owned mesh data stored" );
[ - - ][ - - ]
[ - - ]
561 : : }
562 : :
563 : : // Construct export map (associating mesh connectivities with global node
564 : : // indices and node coordinates) for mesh chunks associated to chare IDs
565 : : // owned by chares we do not own.
566 : : std::unordered_map< int, // target compute node
567 : : std::unordered_map< int, // chare ID
568 : : std::tuple<
569 : : // (domain-element) tetrahedron connectivity
570 : : std::vector< std::size_t >,
571 : : // (domain) node IDs & coordinates
572 : : tk::UnsMesh::CoordMap,
573 : : // boundary side set + face connectivity
574 : : std::unordered_map< int, std::vector< std::size_t > >,
575 : : // boundary side set + node list
576 : : std::unordered_map< int, std::vector< std::size_t > >,
577 : : // Mesh block ids + local tet ids
578 : : std::vector< std::size_t >
579 : 1148 : > > > exp;
580 : :
581 [ + + ]: 6908 : for (const auto& c : mesh)
582 [ + - ][ + - ]: 6334 : exp[ node(c.first) ][ c.first ] =
[ + - ]
583 [ + - ]: 12668 : std::make_tuple( std::get<0>(c.second),
584 [ + - ]: 12668 : coordmap(std::get<0>(c.second)),
585 : 6334 : std::get<1>(c.second),
586 : 6334 : std::get<2>(c.second),
587 : 12668 : std::get<3>(c.second) );
588 : :
589 : : // Export chare IDs and mesh we do not own to fellow compute nodes
590 [ + + ]: 574 : if (exp.empty()) {
591 [ - + ][ - - ]: 80 : if ( g_inputdeck.get< tag::cmd, tag::feedback >() ) m_host.pedistributed();
592 [ + - ]: 80 : contribute( sizeof(std::size_t), &m_meshid, CkReduction::nop,
593 : 80 : m_cbp.get< tag::distributed >() );
594 : : } else {
595 : 494 : m_ndist += exp.size();
596 [ + + ]: 2622 : for (const auto& [ targetnode, chunk ] : exp)
597 [ + - ][ + - ]: 2128 : thisProxy[ targetnode ].addMesh( CkMyNode(), chunk );
598 : : }
599 : 574 : }
600 : :
601 : : std::array< int, 2 >
602 : 1148 : Partitioner::distribution( int npart ) const
603 : : // *****************************************************************************
604 : : // Compute chare (partition) distribution
605 : : //! \param[in] npart Total number of chares (partitions) to distribute
606 : : //! \return Chunksize, i.e., number of chares per all compute nodes except the
607 : : //! last one, and the number of chares for this compute node.
608 : : //! \details Chare ids are distributed to compute nodes in a linear continguous
609 : : //! order with the last compute node taking the remainder if the number of
610 : : //! compute nodes is not divisible by the number chares. For example, if
611 : : //! nchare=7 and nnode=3, the chare distribution is node0: 0 1, node1: 2 3,
612 : : //! and node2: 4 5 6. As a result of this distribution, all compute nodes will
613 : : //! have their chare-categorized element connectivity filled with the global
614 : : //! mesh node IDs associated to the Charm++ chare IDs each compute node owns.
615 : : // *****************************************************************************
616 : : {
617 : 1148 : auto chunksize = npart / CkNumNodes();
618 : 1148 : auto mynchare = chunksize;
619 [ + + ]: 1148 : if (CkMyNode() == CkNumNodes()-1) mynchare += npart % CkNumNodes();
620 : 1148 : return {{ chunksize, mynchare }};
621 : : }
622 : :
623 : 0 : Partitioner::~Partitioner()
624 : : // *****************************************************************************
625 : : // Destructor
626 : : // *****************************************************************************
627 : : {
628 : : // The following call has to be made on all MPI ranks to free all resources.
629 : : // see https://github.com/trilinos/Trilinos/issues/11197#issuecomment-1301325163
630 : 0 : Kokkos::finalize();
631 : 0 : }
632 : :
633 : : #include "NoWarning/partitioner.def.h"
|
