metisLikeDecomp.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2017-2022 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28#include "metisLikeDecomp.H"
29#include "Time.H"
30#include "globalIndex.H"
31
32// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
33
35(
36 const labelList& adjncy,
37 const labelList& xadj,
38 const List<scalar>& cWeights,
39 labelList& decomp
40) const
41{
42 if (!Pstream::parRun())
43 {
44 return decomposeSerial
45 (
46 adjncy,
47 xadj,
48 cWeights,
49 decomp
50 );
51 }
52
53 if (debug)
54 {
55 Info<< type() << "Decomp : running in parallel."
56 << " Decomposing all of graph on master processor." << endl;
57 }
58
59 // Protect against zero-sized offset list
60 const label numCells = max(0, (xadj.size()-1));
61
62 const globalIndex globalAdjncy(adjncy.size());
63 const globalIndex globalCells(numCells);
64
65 List<label> allAdjncy(globalAdjncy.gather(adjncy));
66
67 // Gathering xadj to master is similar to globalIndex gather()
68 // except for the following:
69 //
70 // - gathered list is size+1
71 // - apply local to global renumbering
72
74 const label startOfRequests = UPstream::nRequests();
75
76
77 List<label> allXadj;
78 if (Pstream::master())
79 {
80 allXadj.resize(globalCells.totalSize()+1);
81 allXadj.last() = globalAdjncy.totalSize(); // Final end offset
82
83 // My values - no renumbering required
84 SubList<label>(allXadj, globalCells.localSize(0)) =
85 SubList<label>(xadj, globalCells.localSize(0));
86
87 for (const int proci : globalCells.subProcs())
88 {
89 SubList<label> procSlot(allXadj, globalCells.range(proci));
90
91 if (procSlot.empty())
92 {
93 // Nothing to do
94 }
95 else
96 {
98 (
99 commsType,
100 proci,
101 procSlot.data_bytes(),
102 procSlot.size_bytes(),
105 );
106 }
107 }
108 }
109 else
110 {
111 // Send my part of the graph (local numbering)
112
113 if (!numCells)
114 {
115 // Nothing to do
116 }
117 else
118 {
119 SubList<label> procSlot(xadj, numCells);
120
122 (
123 commsType,
125 procSlot.cdata_bytes(),
126 procSlot.size_bytes(),
129 );
130 }
131 }
132
133 if (commsType == UPstream::commsTypes::nonBlocking)
134 {
135 // Wait for all to finish
136 UPstream::waitRequests(startOfRequests);
137 }
138
139 // Local to global renumbering
140 if (Pstream::master())
141 {
142 for (const int proci : globalCells.subProcs())
143 {
144 SubList<label> procSlot(allXadj, globalCells.range(proci));
145
146 globalAdjncy.inplaceToGlobal(proci, procSlot);
147 }
148 }
149
150 // Ignore zero-sized weights ... and poorly sized ones too
151 List<scalar> allWeights;
152 if
153 (
155 (
156 (cWeights.size() == globalCells.localSize()), andOp<bool>()
157 )
158 )
159 {
160 allWeights = globalCells.gather(cWeights);
161 }
162
163
164 // Global decomposition
165 labelList allDecomp;
166
167 if (Pstream::master())
168 {
170 (
171 allAdjncy,
172 allXadj,
173 allWeights,
174 allDecomp
175 );
176
177 allAdjncy.clear(); // Not needed anymore
178 allXadj.clear(); // ...
179 allWeights.clear(); // ...
180 }
181
182 // The processor-local decomposition (output)
183 decomp.resize_nocopy(globalCells.localSize());
184 globalCells.scatter(allDecomp, decomp);
185
186 return 0;
187}
188
189
190// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
191
193(
194 const word& derivedType,
195 const dictionary& decompDict,
196 const word& regionName,
197 int select
198)
199:
200 decompositionMethod(decompDict, regionName),
201 coeffsDict_(findCoeffsDict(derivedType + "Coeffs", select))
202{}
203
204
205// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
206
208(
209 const polyMesh& mesh,
210 const pointField& points,
211 const scalarField& pointWeights
212) const
213{
214 if (points.size() != mesh.nCells())
215 {
217 << "Can only use this decomposition method for entire mesh" << nl
218 << "and supply one coordinate (cellCentre) for every cell." << nl
219 << "The number of coordinates " << points.size() << nl
220 << "The number of cells in the mesh " << mesh.nCells() << nl
221 << exit(FatalError);
222 }
223
224 CompactListList<label> cellCells;
225 calcCellCells
226 (
227 mesh,
229 mesh.nCells(),
230 true,
231 cellCells
232 );
233
234 // Decompose using default weights
235 labelList decomp;
236 decomposeGeneral
237 (
238 cellCells.values(),
239 cellCells.offsets(),
240 pointWeights,
241 decomp
242 );
243
244 return decomp;
245}
246
247
249(
250 const polyMesh& mesh,
251 const labelList& agglom,
252 const pointField& agglomPoints,
253 const scalarField& agglomWeights
254) const
255{
256 if (agglom.size() != mesh.nCells())
257 {
259 << "Size of cell-to-coarse map " << agglom.size()
260 << " differs from number of cells in mesh " << mesh.nCells()
261 << exit(FatalError);
262 }
263
264 // Make Metis CSR (Compressed Storage Format) storage
265 // adjncy : contains neighbours (= edges in graph)
266 // xadj(celli) : start of information in adjncy for celli
267
268 CompactListList<label> cellCells;
269 calcCellCells
270 (
271 mesh,
272 agglom,
273 agglomPoints.size(),
274 true,
275 cellCells
276 );
277
278 // Decompose using default weights
279 labelList decomp;
280 decomposeGeneral
281 (
282 cellCells.values(),
283 cellCells.offsets(),
284 agglomWeights,
285 decomp
286 );
287
288
289 // Rework back into decomposition for original mesh
290 labelList fineDistribution(agglom.size());
291
292 forAll(fineDistribution, i)
293 {
294 fineDistribution[i] = decomp[agglom[i]];
295 }
296
297 return fineDistribution;
298}
299
300
302(
303 const labelListList& globalCellCells,
304 const pointField& cellCentres,
305 const scalarField& cellWeights
306) const
307{
308 if (cellCentres.size() != globalCellCells.size())
309 {
311 << "Inconsistent number of cells (" << globalCellCells.size()
312 << ") and number of cell centres (" << cellCentres.size()
313 << ")." << exit(FatalError);
314 }
315
316 // Make Metis CSR (Compressed Storage Format) storage
317 // adjncy : contains neighbours (= edges in graph)
318 // xadj(celli) : start of information in adjncy for celli
319
320
321 auto cellCells(CompactListList<label>::pack(globalCellCells));
322
323 // Decompose using default weights
324 labelList decomp;
325 decomposeGeneral
326 (
327 cellCells.values(),
328 cellCells.offsets(),
329 cellWeights,
330 decomp
331 );
332
333 return decomp;
334}
335
336// ************************************************************************* //
A packed storage unstructured matrix of objects of type <T> using an offset table for access.
const List< T > & values() const noexcept
Return the packed matrix of values.
const labelList & offsets() const noexcept
Return the offset table (= size()+1)
void resize_nocopy(const label len)
Adjust allocated size of list without necessarily.
Definition: ListI.H:146
void resize(const label len)
Adjust allocated size of list.
Definition: ListI.H:139
void clear()
Clear the list, i.e. set size to zero.
Definition: ListI.H:116
virtual bool read()
Re-read model coefficients if they have changed.
A List obtained as a section of another List.
Definition: SubList.H:70
char * data_bytes() noexcept
Return pointer to the underlying array serving as data storage,.
Definition: UListI.H:251
bool empty() const noexcept
True if the UList is empty (ie, size() is zero)
Definition: UListI.H:427
const char * cdata_bytes() const noexcept
Return pointer to the underlying array serving as data storage,.
Definition: UListI.H:244
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
std::streamsize size_bytes() const noexcept
Number of contiguous bytes for the List data.
Definition: UListI.H:258
T & last()
Return the last element of the list.
Definition: UListI.H:216
commsTypes
Types of communications.
Definition: UPstream.H:67
@ nonBlocking
"nonBlocking"
static int & msgType() noexcept
Message tag of standard messages.
Definition: UPstream.H:556
static constexpr int masterNo() noexcept
Process index of the master (always 0)
Definition: UPstream.H:451
static label worldComm
Default communicator (all processors)
Definition: UPstream.H:293
static label nRequests()
Get number of outstanding requests.
Definition: UPstream.C:90
static void waitRequests(const label start=0)
Wait until all requests (from start onwards) have finished.
Definition: UPstream.C:100
static bool & parRun() noexcept
Test if this a parallel run.
Definition: UPstream.H:433
Abstract base class for domain decomposition.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
virtual bool write()
Write the output fields.
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:68
label localSize() const
My local size.
Definition: globalIndexI.H:207
labelRange range() const
Return start/size range of local processor data.
Definition: globalIndexI.H:232
static void gather(const labelUList &offsets, const label comm, const ProcIDsContainer &procIDs, const UList< Type > &fld, List< Type > &allFld, const int tag=UPstream::msgType(), const UPstream::commsTypes=UPstream::commsTypes::nonBlocking)
Collect data in processor order on master (== procIDs[0]).
void inplaceToGlobal(labelUList &labels) const
From local to global index (inplace)
Definition: globalIndexI.H:303
static void scatter(const labelUList &offsets, const label comm, const ProcIDsContainer &procIDs, const UList< Type > &allFld, UList< Type > &fld, const int tag=UPstream::msgType(), const UPstream::commsTypes=UPstream::commsTypes::nonBlocking)
Distribute data in processor order.
labelRange subProcs() const noexcept
Range of process indices for addressed sub-offsets (processes)
Definition: globalIndexI.H:159
label totalSize() const
Global sum of localSizes.
Definition: globalIndexI.H:125
Domain decomposition using METIS-like data structures.
virtual label decomposeSerial(const labelList &adjncy, const labelList &xadj, const List< scalar > &cellWeights, labelList &decomp) const =0
Decomposition with metis-like parameters.
virtual label decomposeGeneral(const labelList &adjncy, const labelList &xadj, const List< scalar > &cellWeights, labelList &decomp) const
Serial and/or collect/distribute for parallel operation.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
label nCells() const noexcept
Number of mesh cells.
splitCell * master() const
Definition: splitCell.H:113
bool decompose() const noexcept
Query the decompose flag (normally off)
A class for handling words, derived from Foam::string.
Definition: word.H:68
dynamicFvMesh & mesh
Foam::word regionName(Foam::polyMesh::defaultRegion)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
const pointField & points
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i)
Definition: labelList.C:38
messageStream Info
Information stream (stdout output on master, null elsewhere)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
error FatalError
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce (copy) and return value.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333