domainDecompositionMesh.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) 2011-2016 OpenFOAM Foundation
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
26InClass
27 domainDecomposition
28
29Description
30 Private member of domainDecomposition.
31 Decomposes the mesh into bits
32
33\*---------------------------------------------------------------------------*/
34
35#include "domainDecomposition.H"
36#include "IOstreams.H"
37#include "bitSet.H"
38#include "cyclicPolyPatch.H"
39
40// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
41
42void Foam::domainDecomposition::append(labelList& lst, const label elem)
43{
44 label sz = lst.size();
45 lst.setSize(sz+1);
46 lst[sz] = elem;
47}
48
49
50void Foam::domainDecomposition::addInterProcFace
51(
52 const label facei,
53 const label ownerProc,
54 const label nbrProc,
55
56 List<Map<label>>& nbrToInterPatch,
57 List<DynamicList<DynamicList<label>>>& interPatchFaces
58) const
59{
60 // Introduce turning index only for internal faces (are duplicated).
61 const label ownerIndex = facei+1;
62 const label nbrIndex = -(facei+1);
63
64 const auto patchiter = nbrToInterPatch[ownerProc].cfind(nbrProc);
65
66 if (patchiter.found())
67 {
68 // Existing interproc patch. Add to both sides.
69 const label toNbrProcPatchi = *patchiter;
70 interPatchFaces[ownerProc][toNbrProcPatchi].append(ownerIndex);
71
72 if (isInternalFace(facei))
73 {
74 label toOwnerProcPatchi = nbrToInterPatch[nbrProc][ownerProc];
75 interPatchFaces[nbrProc][toOwnerProcPatchi].append(nbrIndex);
76 }
77 }
78 else
79 {
80 // Create new interproc patches.
81 const label toNbrProcPatchi = nbrToInterPatch[ownerProc].size();
82 nbrToInterPatch[ownerProc].insert(nbrProc, toNbrProcPatchi);
83
84 DynamicList<label> oneFace;
85 oneFace.append(ownerIndex);
86 interPatchFaces[ownerProc].append(oneFace);
87
88 if (isInternalFace(facei))
89 {
90 label toOwnerProcPatchi = nbrToInterPatch[nbrProc].size();
91 nbrToInterPatch[nbrProc].insert(ownerProc, toOwnerProcPatchi);
92 oneFace.clear();
93 oneFace.append(nbrIndex);
94 interPatchFaces[nbrProc].append(oneFace);
95 }
96 }
97}
98
99
101{
102 // Decide which cell goes to which processor
103 distributeCells();
104
105 // Distribute the cells according to the given processor label
106
107 // calculate the addressing information for the original mesh
108 Info<< "\nCalculating original mesh data" << endl;
109
110 // set references to the original mesh
111 const polyBoundaryMesh& patches = boundaryMesh();
112 const faceList& fcs = faces();
113 const labelList& owner = faceOwner();
114 const labelList& neighbour = faceNeighbour();
115
116 // loop through the list of processor labels for the cell and add the
117 // cell shape to the list of cells for the appropriate processor
118
119 Info<< "\nDistributing cells to processors" << endl;
120
121 // Cells per processor
122 procCellAddressing_ = invertOneToMany(nProcs_, cellToProc_);
123
124 Info<< "\nDistributing faces to processors" << endl;
125
126 // Loop through all internal faces and decide which processor they belong to
127 // First visit all internal faces. If cells at both sides belong to the
128 // same processor, the face is an internal face. If they are different,
129 // it belongs to both processors.
130
131 procFaceAddressing_.setSize(nProcs_);
132
133 // Internal faces
134 forAll(neighbour, facei)
135 {
136 if (cellToProc_[owner[facei]] == cellToProc_[neighbour[facei]])
137 {
138 // Face internal to processor. Notice no turning index.
139 procFaceAddressing_[cellToProc_[owner[facei]]].append(facei+1);
140 }
141 }
142
143 // for all processors, set the size of start index and patch size
144 // lists to the number of patches in the mesh
145 forAll(procPatchSize_, proci)
146 {
147 procPatchSize_[proci].setSize(patches.size());
148 procPatchStartIndex_[proci].setSize(patches.size());
149 }
150
151 forAll(patches, patchi)
152 {
153 // Reset size and start index for all processors
154 forAll(procPatchSize_, proci)
155 {
156 procPatchSize_[proci][patchi] = 0;
157 procPatchStartIndex_[proci][patchi] =
158 procFaceAddressing_[proci].size();
159 }
160
161 const label patchStart = patches[patchi].start();
162
163 if (!isA<cyclicPolyPatch>(patches[patchi]))
164 {
165 // Normal patch. Add faces to processor where the cell
166 // next to the face lives
167
168 const labelUList& patchFaceCells =
169 patches[patchi].faceCells();
170
171 forAll(patchFaceCells, facei)
172 {
173 const label curProc = cellToProc_[patchFaceCells[facei]];
174
175 // add the face without turning index
176 procFaceAddressing_[curProc].append(patchStart+facei+1);
177
178 // increment the number of faces for this patch
179 procPatchSize_[curProc][patchi]++;
180 }
181 }
182 else
183 {
184 const cyclicPolyPatch& pp = refCast<const cyclicPolyPatch>
185 (
186 patches[patchi]
187 );
188 // cyclic: check opposite side on this processor
189 const labelUList& patchFaceCells = pp.faceCells();
190
191 const labelUList& nbrPatchFaceCells =
192 pp.neighbPatch().faceCells();
193
194 forAll(patchFaceCells, facei)
195 {
196 const label curProc = cellToProc_[patchFaceCells[facei]];
197 const label nbrProc = cellToProc_[nbrPatchFaceCells[facei]];
198 if (curProc == nbrProc)
199 {
200 // add the face without turning index
201 procFaceAddressing_[curProc].append(patchStart+facei+1);
202 // increment the number of faces for this patch
203 procPatchSize_[curProc][patchi]++;
204 }
205 }
206 }
207 }
208
209
210 // Done internal bits of the new mesh and the ordinary patches.
211
212
213 // Per processor, from neighbour processor to the inter-processor patch
214 // that communicates with that neighbour
215 List<Map<label>> procNbrToInterPatch(nProcs_);
216
217 // Per processor the faces per inter-processor patch
218 List<DynamicList<DynamicList<label>>> interPatchFaces(nProcs_);
219
220 // Processor boundaries from internal faces
221 forAll(neighbour, facei)
222 {
223 label ownerProc = cellToProc_[owner[facei]];
224 label nbrProc = cellToProc_[neighbour[facei]];
225
226 if (ownerProc != nbrProc)
227 {
228 // inter - processor patch face found.
229 addInterProcFace
230 (
231 facei,
232 ownerProc,
233 nbrProc,
234
235 procNbrToInterPatch,
236 interPatchFaces
237 );
238 }
239 }
240
241 // Add the proper processor faces to the sub information. For faces
242 // originating from internal faces this is always -1.
243 List<labelListList> subPatchIDs(nProcs_);
244 List<labelListList> subPatchStarts(nProcs_);
245 forAll(interPatchFaces, proci)
246 {
247 label nInterfaces = interPatchFaces[proci].size();
248
249 subPatchIDs[proci].setSize(nInterfaces, labelList(1, label(-1)));
250 subPatchStarts[proci].setSize(nInterfaces, labelList(1, Zero));
251 }
252
253
254 // Special handling needed for the case that multiple processor cyclic
255 // patches are created on each local processor domain, e.g. if a 3x3 case
256 // is decomposed using the decomposition:
257 //
258 // | 1 | 0 | 2 |
259 // cyclic left | 2 | 0 | 1 | cyclic right
260 // | 2 | 0 | 1 |
261 //
262 // - processors 1 and 2 will both have pieces of both cyclic left- and
263 // right sub-patches present
264 // - the interface patch faces are stored in a single list, where each
265 // sub-patch is referenced into the list using a patch start index and
266 // size
267 // - if the patches are in order (in the boundary file) of left, right
268 // - processor 1 will send: left, right
269 // - processor 1 will need to receive in reverse order: right, left
270 // - similarly for processor 2
271 // - the sub-patches are therefore generated in 4 passes of the patch lists
272 // 1. add faces from owner patch where local proc i < nbr proc i
273 // 2. add faces from nbr patch where local proc i < nbr proc i
274 // 3. add faces from owner patch where local proc i > nbr proc i
275 // 4. add faces from nbr patch where local proc i > nbr proc i
276
277 processInterCyclics
278 (
279 patches,
280 interPatchFaces,
281 procNbrToInterPatch,
282 subPatchIDs,
283 subPatchStarts,
284 true,
285 lessOp<label>()
286 );
287
288 processInterCyclics
289 (
290 patches,
291 interPatchFaces,
292 procNbrToInterPatch,
293 subPatchIDs,
294 subPatchStarts,
295 false,
296 lessOp<label>()
297 );
298
299 processInterCyclics
300 (
301 patches,
302 interPatchFaces,
303 procNbrToInterPatch,
304 subPatchIDs,
305 subPatchStarts,
306 false,
307 greaterOp<label>()
308 );
309
310 processInterCyclics
311 (
312 patches,
313 interPatchFaces,
314 procNbrToInterPatch,
315 subPatchIDs,
316 subPatchStarts,
317 true,
318 greaterOp<label>()
319 );
320
321
322 // Sort inter-proc patch by neighbour
323 labelList order;
324 forAll(procNbrToInterPatch, proci)
325 {
326 label nInterfaces = procNbrToInterPatch[proci].size();
327
328 procNeighbourProcessors_[proci].setSize(nInterfaces);
329 procProcessorPatchSize_[proci].setSize(nInterfaces);
330 procProcessorPatchStartIndex_[proci].setSize(nInterfaces);
331 procProcessorPatchSubPatchIDs_[proci].setSize(nInterfaces);
332 procProcessorPatchSubPatchStarts_[proci].setSize(nInterfaces);
333
334 //Info<< "Processor " << proci << endl;
335
336 // Get sorted neighbour processors
337 const Map<label>& curNbrToInterPatch = procNbrToInterPatch[proci];
338 labelList nbrs = curNbrToInterPatch.toc();
339
340 sortedOrder(nbrs, order);
341
342 DynamicList<DynamicList<label>>& curInterPatchFaces =
343 interPatchFaces[proci];
344
345 forAll(nbrs, i)
346 {
347 const label nbrProc = nbrs[i];
348 const label interPatch = curNbrToInterPatch[nbrProc];
349
350 procNeighbourProcessors_[proci][i] = nbrProc;
351 procProcessorPatchSize_[proci][i] =
352 curInterPatchFaces[interPatch].size();
353 procProcessorPatchStartIndex_[proci][i] =
354 procFaceAddressing_[proci].size();
355
356 // Add size as last element to substarts and transfer
357 append
358 (
359 subPatchStarts[proci][interPatch],
360 curInterPatchFaces[interPatch].size()
361 );
362 procProcessorPatchSubPatchIDs_[proci][i].transfer
363 (
364 subPatchIDs[proci][interPatch]
365 );
366 procProcessorPatchSubPatchStarts_[proci][i].transfer
367 (
368 subPatchStarts[proci][interPatch]
369 );
370
371 //Info<< " nbr:" << nbrProc << endl;
372 //Info<< " interpatch:" << interPatch << endl;
373 //Info<< " size:" << procProcessorPatchSize_[proci][i] << endl;
374 //Info<< " start:" << procProcessorPatchStartIndex_[proci][i]
375 // << endl;
376 //Info<< " subPatches:"
377 // << procProcessorPatchSubPatchIDs_[proci][i]
378 // << endl;
379 //Info<< " subStarts:"
380 // << procProcessorPatchSubPatchStarts_[proci][i] << endl;
381
382 // And add all the face labels for interPatch
383 DynamicList<label>& interPatchFaces =
384 curInterPatchFaces[interPatch];
385
386 forAll(interPatchFaces, j)
387 {
388 procFaceAddressing_[proci].append(interPatchFaces[j]);
389 }
390 interPatchFaces.clearStorage();
391 }
392 curInterPatchFaces.clearStorage();
393 procFaceAddressing_[proci].shrink();
394 }
395
396
399// forAll(procPatchStartIndex_, proci)
400// {
401// Info<< "Processor:" << proci << endl;
402//
403// Info<< " total faces:" << procFaceAddressing_[proci].size()
404// << endl;
405//
406// const labelList& curProcPatchStartIndex = procPatchStartIndex_[proci];
407//
408// forAll(curProcPatchStartIndex, patchi)
409// {
410// Info<< " patch:" << patchi
411// << "\tstart:" << curProcPatchStartIndex[patchi]
412// << "\tsize:" << procPatchSize_[proci][patchi]
413// << endl;
414// }
415// }
416// Info<< endl;
417//
418// forAll(procNeighbourProcessors_, proci)
419// {
420// Info<< "Processor " << proci << endl;
421//
422// forAll(procNeighbourProcessors_[proci], i)
423// {
424// Info<< " nbr:" << procNeighbourProcessors_[proci][i] << endl;
425// Info<< " size:" << procProcessorPatchSize_[proci][i] << endl;
426// Info<< " start:" << procProcessorPatchStartIndex_[proci][i]
427// << endl;
428// }
429// }
430// Info<< endl;
431//
432// forAll(procFaceAddressing_, proci)
433// {
434// Info<< "Processor:" << proci << endl;
435//
436// Info<< " faces:" << procFaceAddressing_[proci] << endl;
437// }
438
439
440 Info<< "\nDistributing points to processors" << endl;
441 // For every processor, loop through the list of faces for the processor.
442 // For every face, loop through the list of points and mark the point as
443 // used for the processor. Collect the list of used points for the
444 // processor.
445
446 forAll(procPointAddressing_, proci)
447 {
448 bitSet pointsInUse(nPoints(), false);
449
450 // For each of the faces used
451 for (const label facei : procFaceAddressing_[proci])
452 {
453 // Because of the turning index, some face labels may be -ve
454 const labelList& facePoints = fcs[mag(facei) - 1];
455
456 // Mark the face points as being used
457 pointsInUse.set(facePoints);
458 }
459
460 procPointAddressing_[proci] = pointsInUse.sortedToc();
461 }
462}
463
464// ************************************************************************* //
Useful combination of include files which define Sin, Sout and Serr and the use of IO streams general...
std::enable_if< std::is_same< bool, TypeT >::value, bool >::type set(const label i, bool val=true)
A bitSet::set() method for a list of bool.
Definition: List.H:330
label size() const noexcept
The number of elements in the list.
Definition: UPtrListI.H:106
void decomposeMesh()
Decompose mesh.
UPtrList< const labelUList > faceCells() const
Return a list of faceCells for each patch.
label start() const noexcept
The start label of boundary faces in the polyMesh face list.
bool append() const noexcept
True if output format uses an append mode.
const polyBoundaryMesh & patches
label nPoints
rAUs append(new volScalarField(IOobject::groupName("rAU", phase1.name()), 1.0/(U1Eqn.A()+byDt(max(phase1.residualAlpha() - alpha1, scalar(0)) *rho1))))
List< label > labelList
A List of labels.
Definition: List.H:66
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
labelList sortedOrder(const UList< T > &input)
Return the (stable) sort order for the list.
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
labelListList invertOneToMany(const label len, const labelUList &map)
Invert one-to-many map. Unmapped elements will be size 0.
Definition: ListOps.C:114
List< face > faceList
A List of faces.
Definition: faceListFwd.H:47
UList< label > labelUList
A UList of labels.
Definition: UList.H:85
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333