All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
regionSplit.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "regionSplit.H"
27 #include "FaceCellWave.H"
28 #include "cyclicPolyPatch.H"
29 #include "processorPolyPatch.H"
30 #include "globalIndex.H"
31 #include "syncTools.H"
32 #include "minData.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38  defineTypeNameAndDebug(regionSplit, 0);
39 }
40 
41 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
42 
43 void Foam::regionSplit::calcNonCompactRegionSplit
44 (
45  const globalIndex& globalFaces,
46  const boolList& blockedFace,
47  const List<labelPair>& explicitConnections,
48 
49  labelList& cellRegion
50 ) const
51 {
52  // Field on cells and faces.
53  List<minData> cellData(mesh().nCells());
54  List<minData> faceData(mesh().nFaces());
55 
56  // Take over blockedFaces by seeding a negative number
57  // (so is always less than the decomposition)
58  label nUnblocked = 0;
59  forAll(faceData, facei)
60  {
61  if (blockedFace.size() && blockedFace[facei])
62  {
63  faceData[facei] = minData(-2);
64  }
65  else
66  {
67  nUnblocked++;
68  }
69  }
70 
71  // Seed unblocked faces
72  labelList seedFaces(nUnblocked);
73  List<minData> seedData(nUnblocked);
74  nUnblocked = 0;
75 
76 
77  forAll(faceData, facei)
78  {
79  if (blockedFace.empty() || !blockedFace[facei])
80  {
81  seedFaces[nUnblocked] = facei;
82  // Seed face with globally unique number
83  seedData[nUnblocked] = minData(globalFaces.toGlobal(facei));
84  nUnblocked++;
85  }
86  }
87 
88 
89  // Propagate information inwards
90  FaceCellWave<minData> deltaCalc
91  (
92  mesh(),
93  explicitConnections,
94  false, // disable walking through cyclicAMI for backwards compatibility
95  seedFaces,
96  seedData,
97  faceData,
98  cellData,
99  mesh().globalData().nTotalCells()+1
100  );
101 
102 
103  // And extract
104  cellRegion.setSize(mesh().nCells());
105  forAll(cellRegion, celli)
106  {
107  if (cellData[celli].valid(deltaCalc.data()))
108  {
109  cellRegion[celli] = cellData[celli].data();
110  }
111  else
112  {
113  // Unvisited cell -> only possible if surrounded by blocked faces.
114  // If so make up region from any of the faces
115  const cell& cFaces = mesh().cells()[celli];
116  label facei = cFaces[0];
117 
118  if (blockedFace.size() && !blockedFace[facei])
119  {
120  FatalErrorIn("regionSplit::calcNonCompactRegionSplit(..)")
121  << "Problem: unblocked face " << facei
122  << " at " << mesh().faceCentres()[facei]
123  << " on unassigned cell " << celli
124  << mesh().cellCentres()[celli]
125  << exit(FatalError);
126  }
127  cellRegion[celli] = globalFaces.toGlobal(facei);
128  }
129  }
130 }
131 
132 
133 Foam::autoPtr<Foam::globalIndex> Foam::regionSplit::calcRegionSplit
134 (
135  const bool doGlobalRegions,
136  const boolList& blockedFace,
137  const List<labelPair>& explicitConnections,
138 
139  labelList& cellRegion
140 ) const
141 {
142  // See header in regionSplit.H
143 
144 
145  if (!doGlobalRegions)
146  {
147  // Block all parallel faces to avoid comms across
148  boolList coupledOrBlockedFace(blockedFace);
149  const polyBoundaryMesh& pbm = mesh().boundaryMesh();
150 
151  if (coupledOrBlockedFace.size())
152  {
153  forAll(pbm, patchi)
154  {
155  const polyPatch& pp = pbm[patchi];
156  if (isA<processorPolyPatch>(pp))
157  {
158  label facei = pp.start();
159  forAll(pp, i)
160  {
161  coupledOrBlockedFace[facei++] = true;
162  }
163  }
164  }
165  }
166 
167  // Create dummy (local only) globalIndex
168  labelList offsets(Pstream::nProcs()+1, 0);
169  for (label i = Pstream::myProcNo()+1; i < offsets.size(); i++)
170  {
171  offsets[i] = mesh().nFaces();
172  }
173  const globalIndex globalRegions(move(offsets));
174 
175  // Minimize regions across connected cells
176  // Note: still uses global decisions so all processors are running
177  // in lock-step, i.e. slowest determines overall time.
178  // To avoid this we could switch off Pstream::parRun.
179  calcNonCompactRegionSplit
180  (
181  globalRegions,
182  coupledOrBlockedFace,
183  explicitConnections,
184  cellRegion
185  );
186 
187  // Compact
188  Map<label> globalToCompact(mesh().nCells()/8);
189  forAll(cellRegion, celli)
190  {
191  label region = cellRegion[celli];
192 
193  label globalRegion;
194 
195  Map<label>::const_iterator fnd = globalToCompact.find(region);
196  if (fnd == globalToCompact.end())
197  {
198  globalRegion = globalRegions.toGlobal(globalToCompact.size());
199  globalToCompact.insert(region, globalRegion);
200  }
201  else
202  {
203  globalRegion = fnd();
204  }
205  cellRegion[celli] = globalRegion;
206  }
207 
208 
209  // Return globalIndex with size = localSize and all regions local
210  labelList compactOffsets(Pstream::nProcs()+1, 0);
211  for (label i = Pstream::myProcNo()+1; i < compactOffsets.size(); i++)
212  {
213  compactOffsets[i] = globalToCompact.size();
214  }
215 
216  return autoPtr<globalIndex>(new globalIndex(move(compactOffsets)));
217  }
218 
219 
220 
221  // Initial global region numbers
222  const globalIndex globalRegions(mesh().nFaces());
223 
224  // Minimize regions across connected cells (including parallel)
225  calcNonCompactRegionSplit
226  (
227  globalRegions,
228  blockedFace,
229  explicitConnections,
230  cellRegion
231  );
232 
233 
234  // Now our cellRegion will have
235  // - non-local regions (i.e. originating from other processors)
236  // - non-compact locally originating regions
237  // so we'll need to compact
238 
239  // 4a: count per originating processor the number of regions
240  labelList nOriginating(Pstream::nProcs(), 0);
241  {
242  labelHashSet haveRegion(mesh().nCells()/8);
243 
244  forAll(cellRegion, celli)
245  {
246  label region = cellRegion[celli];
247 
248  // Count originating processor. Use isLocal as efficiency since
249  // most cells are locally originating.
250  if (globalRegions.isLocal(region))
251  {
252  if (haveRegion.insert(region))
253  {
254  nOriginating[Pstream::myProcNo()]++;
255  }
256  }
257  else
258  {
259  label proci = globalRegions.whichProcID(region);
260  if (haveRegion.insert(region))
261  {
262  nOriginating[proci]++;
263  }
264  }
265  }
266  }
267 
268  if (debug)
269  {
270  Pout<< "Counted " << nOriginating[Pstream::myProcNo()]
271  << " local regions." << endl;
272  }
273 
274 
275  // Global numbering for compacted local regions
276  autoPtr<globalIndex> globalCompactPtr
277  (
278  new globalIndex(nOriginating[Pstream::myProcNo()])
279  );
280  const globalIndex& globalCompact = globalCompactPtr();
281 
282 
283  // 4b: renumber
284  // Renumber into compact indices. Note that since we've already made
285  // all regions global we now need a Map to store the compacting information
286  // instead of a labelList - otherwise we could have used a straight
287  // labelList.
288 
289  // Local compaction map
290  Map<label> globalToCompact(2*nOriginating[Pstream::myProcNo()]);
291  // Remote regions we want the compact number for
292  List<labelHashSet> nonLocal(Pstream::nProcs());
293  forAll(nonLocal, proci)
294  {
295  if (proci != Pstream::myProcNo())
296  {
297  nonLocal[proci].resize(2*nOriginating[proci]);
298  }
299  }
300 
301  forAll(cellRegion, celli)
302  {
303  label region = cellRegion[celli];
304  if (globalRegions.isLocal(region))
305  {
306  // Insert new compact region (if not yet present)
307  globalToCompact.insert
308  (
309  region,
310  globalCompact.toGlobal(globalToCompact.size())
311  );
312  }
313  else
314  {
315  nonLocal[globalRegions.whichProcID(region)].insert(region);
316  }
317  }
318 
319 
320  // Now we have all the local regions compacted. Now we need to get the
321  // non-local ones from the processors to whom they are local.
322  // Convert the nonLocal (labelHashSets) to labelLists.
323 
324  labelListList sendNonLocal(Pstream::nProcs());
325  forAll(sendNonLocal, proci)
326  {
327  sendNonLocal[proci] = nonLocal[proci].toc();
328  }
329 
330  if (debug)
331  {
332  forAll(sendNonLocal, proci)
333  {
334  Pout<< " from processor " << proci
335  << " want " << sendNonLocal[proci].size()
336  << " region numbers."
337  << endl;
338  }
339  Pout<< endl;
340  }
341 
342 
343  // Get the wanted region labels into recvNonLocal
344  labelListList recvNonLocal;
345  Pstream::exchange<labelList, label>(sendNonLocal, recvNonLocal);
346 
347  // Now we have the wanted compact region labels that proci wants in
348  // recvNonLocal[proci]. Construct corresponding list of compact
349  // region labels to send back.
350 
351  labelListList sendWantedLocal(Pstream::nProcs());
352  forAll(recvNonLocal, proci)
353  {
354  const labelList& nonLocal = recvNonLocal[proci];
355  sendWantedLocal[proci].setSize(nonLocal.size());
356 
357  forAll(nonLocal, i)
358  {
359  sendWantedLocal[proci][i] = globalToCompact[nonLocal[i]];
360  }
361  }
362 
363 
364  // Send back (into recvNonLocal)
365  recvNonLocal.clear();
366  Pstream::exchange<labelList, label>(sendWantedLocal, recvNonLocal);
367  sendWantedLocal.clear();
368 
369  // Now recvNonLocal contains for every element in setNonLocal the
370  // corresponding compact number. Insert these into the local compaction
371  // map.
372 
373  forAll(recvNonLocal, proci)
374  {
375  const labelList& wantedRegions = sendNonLocal[proci];
376  const labelList& compactRegions = recvNonLocal[proci];
377 
378  forAll(wantedRegions, i)
379  {
380  globalToCompact.insert(wantedRegions[i], compactRegions[i]);
381  }
382  }
383 
384  // Finally renumber the regions
385  forAll(cellRegion, celli)
386  {
387  cellRegion[celli] = globalToCompact[cellRegion[celli]];
388  }
389 
390  return globalCompactPtr;
391 }
392 
393 
394 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
395 
396 Foam::regionSplit::regionSplit(const polyMesh& mesh, const bool doGlobalRegions)
397 :
399  labelList(mesh.nCells(), -1)
400 {
401  globalNumberingPtr_ = calcRegionSplit
402  (
403  doGlobalRegions, // do global regions
404  boolList(0, false), // blockedFaces
405  List<labelPair>(0), // explicitConnections,
406  *this
407  );
408 }
409 
410 
412 (
413  const polyMesh& mesh,
414  const boolList& blockedFace,
415  const bool doGlobalRegions
416 )
417 :
419  labelList(mesh.nCells(), -1)
420 {
421  globalNumberingPtr_ = calcRegionSplit
422  (
423  doGlobalRegions,
424  blockedFace, // blockedFaces
425  List<labelPair>(0), // explicitConnections,
426  *this
427  );
428 }
429 
430 
432 (
433  const polyMesh& mesh,
434  const boolList& blockedFace,
435  const List<labelPair>& explicitConnections,
436  const bool doGlobalRegions
437 )
438 :
440  labelList(mesh.nCells(), -1)
441 {
442  globalNumberingPtr_ = calcRegionSplit
443  (
444  doGlobalRegions,
445  blockedFace, // blockedFaces
446  explicitConnections, // explicitConnections,
447  *this
448  );
449 }
450 
451 
452 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:434
This class separates the mesh into distinct unconnected regions, each of which is then given a label ...
Definition: regionSplit.H:119
List< labelList > labelListList
A List of labelList.
Definition: labelList.H:57
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
label nFaces() const
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:429
label nCells() const
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
HashTable< label, label, Hash< label > >::const_iterator const_iterator
Definition: Map.H:59
const cellList & cells() const
List< bool > boolList
Bool container classes.
Definition: boolList.H:50
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys.
Definition: HashSet.H:211
Templated abstract base-class for optional mesh objects used to automate their allocation to the mesh...
Definition: MeshObject.H:86
dynamicFvMesh & mesh
List< label > labelList
A List of labels.
Definition: labelList.H:56
const vectorField & cellCentres() const
defineTypeNameAndDebug(combustionModel, 0)
regionSplit(const polyMesh &, const bool doGlobalRegions=Pstream::parRun())
Construct from mesh.
Definition: regionSplit.C:396
const vectorField & faceCentres() const
void setSize(const label)
Reset size of List.
Definition: List.C:281
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:411
label patchi
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Namespace for OpenFOAM.