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 -------------------------------------------------------------------------------
10 License
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 InClass
27  domainDecomposition
28 
29 Description
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 
42 void 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 
50 void 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 // ************************************************************************* //
Foam::labelList
List< label > labelList
A List of labels.
Definition: List.H:71
IOstreams.H
Useful combination of include files which define Sin, Sout and Serr and the use of IO streams general...
cyclicPolyPatch.H
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:350
append
rAUs append(new volScalarField(IOobject::groupName("rAU", phase1.name()), 1.0/(U1Eqn.A()+byDt(max(phase1.residualAlpha() - alpha1, scalar(0)) *rho1))))
bitSet.H
Foam::polyBoundaryMesh::start
label start() const
The start label of the boundary faces in the polyMesh face list.
Definition: polyBoundaryMesh.C:638
domainDecomposition.H
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
nPoints
label nPoints
Definition: gmvOutputHeader.H:2
Foam::Info
messageStream Info
Information stream (uses stdout - output is on the master only)
Foam::domainDecomposition::decomposeMesh
void decomposeMesh()
Decompose mesh.
Foam::faceList
List< face > faceList
A List of faces.
Definition: faceListFwd.H:47
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
patches
const polyBoundaryMesh & patches
Definition: convertProcessorPatches.H:65
Foam::List::set
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:325
Foam::sortedOrder
labelList sortedOrder(const UList< T > &input)
Return the (stable) sort order for the list.
Foam::labelUList
UList< label > labelUList
A UList of labels.
Definition: UList.H:80
Foam::invertOneToMany
labelListList invertOneToMany(const label len, const labelUList &map)
Invert one-to-many map. Unmapped elements will be size 0.
Definition: ListOps.C:114