All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
domainDecomposition.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 "domainDecomposition.H"
27 #include "dictionary.H"
28 #include "labelIOList.H"
29 #include "processorPolyPatch.H"
31 #include "fvMesh.H"
32 #include "OSspecific.H"
33 #include "Map.H"
34 #include "DynamicList.H"
35 #include "fvFieldDecomposer.H"
36 #include "IOobjectList.H"
37 #include "cellSet.H"
38 #include "faceSet.H"
39 #include "pointSet.H"
40 #include "decompositionModel.H"
41 #include "hexRef8Data.H"
42 
43 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
44 
45 void Foam::domainDecomposition::mark
46 (
47  const labelList& zoneElems,
48  const label zoneI,
49  labelList& elementToZone
50 )
51 {
52  forAll(zoneElems, i)
53  {
54  label pointi = zoneElems[i];
55 
56  if (elementToZone[pointi] == -1)
57  {
58  // First occurrence
59  elementToZone[pointi] = zoneI;
60  }
61  else if (elementToZone[pointi] >= 0)
62  {
63  // Multiple zones
64  elementToZone[pointi] = -2;
65  }
66  }
67 }
68 
69 
70 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
71 
73 (
74  const IOobject& io,
75  const fileName& dictFile
76 )
77 :
78  fvMesh(io),
79  facesInstancePointsPtr_
80  (
82  ? new pointIOField
83  (
84  IOobject
85  (
86  "points",
87  facesInstance(),
89  *this,
92  false
93  )
94  )
95  : nullptr
96  ),
97  nProcs_
98  (
99  decompositionModel::New
100  (
101  *this,
102  dictFile
103  ).lookup<int>("numberOfSubdomains")
104  ),
105  distributed_(false),
106  cellToProc_(nCells()),
107  procPointAddressing_(nProcs_),
108  procFaceAddressing_(nProcs_),
109  procCellAddressing_(nProcs_),
110  procPatchSize_(nProcs_),
111  procPatchStartIndex_(nProcs_),
112  procNeighbourProcessors_(nProcs_),
113  procProcessorPatchSize_(nProcs_),
114  procProcessorPatchStartIndex_(nProcs_),
115  procProcessorPatchSubPatchIDs_(nProcs_),
116  procProcessorPatchSubPatchStarts_(nProcs_)
117 {
119  (
120  *this,
121  dictFile
122  ).readIfPresent("distributed", distributed_);
123 }
124 
125 
126 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
127 
129 {}
130 
131 
132 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
133 
134 bool Foam::domainDecomposition::writeDecomposition(const bool decomposeSets)
135 {
136  Info<< "\nConstructing processor meshes" << endl;
137 
138  // Mark point/faces/cells that are in zones.
139  // -1 : not in zone
140  // -2 : in multiple zones
141  // >= 0 : in single given zone
142  // This will give direct lookup of elements that are in a single zone
143  // and we'll only have to revert back to searching through all zones
144  // for the duplicate elements
145 
146  // Point zones
147  labelList pointToZone(points().size(), -1);
148 
149  forAll(pointZones(), zoneI)
150  {
151  mark(pointZones()[zoneI], zoneI, pointToZone);
152  }
153 
154  // Face zones
155  labelList faceToZone(faces().size(), -1);
156 
157  forAll(faceZones(), zoneI)
158  {
159  mark(faceZones()[zoneI], zoneI, faceToZone);
160  }
161 
162  // Cell zones
163  labelList cellToZone(nCells(), -1);
164 
165  forAll(cellZones(), zoneI)
166  {
167  mark(cellZones()[zoneI], zoneI, cellToZone);
168  }
169 
170 
171  PtrList<const cellSet> cellSets;
172  PtrList<const faceSet> faceSets;
173  PtrList<const pointSet> pointSets;
174  if (decomposeSets)
175  {
176  // Read sets
177  IOobjectList objects(*this, facesInstance(), "polyMesh/sets");
178  {
179  IOobjectList cSets(objects.lookupClass(cellSet::typeName));
180  forAllConstIter(IOobjectList, cSets, iter)
181  {
182  cellSets.append(new cellSet(*iter()));
183  }
184  }
185  {
186  IOobjectList fSets(objects.lookupClass(faceSet::typeName));
187  forAllConstIter(IOobjectList, fSets, iter)
188  {
189  faceSets.append(new faceSet(*iter()));
190  }
191  }
192  {
193  IOobjectList pSets(objects.lookupClass(pointSet::typeName));
194  forAllConstIter(IOobjectList, pSets, iter)
195  {
196  pointSets.append(new pointSet(*iter()));
197  }
198  }
199  }
200 
201 
202  // Load refinement data (if any)
203  hexRef8Data baseMeshData
204  (
205  IOobject
206  (
207  "dummy",
208  facesInstance(),
210  *this,
213  false
214  )
215  );
216 
217 
218 
219  label maxProcCells = 0;
220  label totProcFaces = 0;
221  label maxProcPatches = 0;
222  label totProcPatches = 0;
223  label maxProcFaces = 0;
224 
225 
226  // Write out the meshes
227  for (label proci = 0; proci < nProcs_; proci++)
228  {
229  // Create processor points
230  const labelList& curPointLabels = procPointAddressing_[proci];
231 
232  const pointField& meshPoints = points();
233 
234  labelList pointLookup(nPoints(), -1);
235 
236  pointField procPoints(curPointLabels.size());
237 
238  forAll(curPointLabels, pointi)
239  {
240  procPoints[pointi] = meshPoints[curPointLabels[pointi]];
241 
242  pointLookup[curPointLabels[pointi]] = pointi;
243  }
244 
245  // Create processor faces
246  const labelList& curFaceLabels = procFaceAddressing_[proci];
247 
248  const faceList& meshFaces = faces();
249 
250  labelList faceLookup(nFaces(), -1);
251 
252  faceList procFaces(curFaceLabels.size());
253 
254  forAll(curFaceLabels, facei)
255  {
256  // Mark the original face as used
257  // Remember to decrement the index by one (turning index)
258  //
259  label curF = mag(curFaceLabels[facei]) - 1;
260 
261  faceLookup[curF] = facei;
262 
263  // get the original face
264  labelList origFaceLabels;
265 
266  if (curFaceLabels[facei] >= 0)
267  {
268  // face not turned
269  origFaceLabels = meshFaces[curF];
270  }
271  else
272  {
273  origFaceLabels = meshFaces[curF].reverseFace();
274  }
275 
276  // translate face labels into local point list
277  face& procFaceLabels = procFaces[facei];
278 
279  procFaceLabels.setSize(origFaceLabels.size());
280 
281  forAll(origFaceLabels, pointi)
282  {
283  procFaceLabels[pointi] = pointLookup[origFaceLabels[pointi]];
284  }
285  }
286 
287  // Create processor cells
288  const labelList& curCellLabels = procCellAddressing_[proci];
289 
290  const cellList& meshCells = cells();
291 
292  cellList procCells(curCellLabels.size());
293 
294  forAll(curCellLabels, celli)
295  {
296  const labelList& origCellLabels = meshCells[curCellLabels[celli]];
297 
298  cell& curCell = procCells[celli];
299 
300  curCell.setSize(origCellLabels.size());
301 
302  forAll(origCellLabels, cellFacei)
303  {
304  curCell[cellFacei] = faceLookup[origCellLabels[cellFacei]];
305  }
306  }
307 
308  // Create processor mesh without a boundary
309 
310  fileName processorCasePath
311  (
312  time().caseName()/fileName(word("processor") + Foam::name(proci))
313  );
314 
315  // create a database
316  Time processorDb
317  (
319  time().rootPath(),
320  processorCasePath,
321  word("system"),
322  word("constant")
323  );
324  processorDb.setTime(time());
325 
326  // create the mesh. Two situations:
327  // - points and faces come from the same time ('instance'). The mesh
328  // will get constructed in the same instance.
329  // - points come from a different time (moving mesh cases).
330  // It will read the points belonging to the faces instance and
331  // construct the procMesh with it which then gets handled as above.
332  // (so with 'old' geometry).
333  // Only at writing time will it additionally write the current
334  // points.
335 
336  autoPtr<polyMesh> procMeshPtr;
337 
338  if (facesInstancePointsPtr_.valid())
339  {
340  // Construct mesh from facesInstance.
341  pointField facesInstancePoints
342  (
343  facesInstancePointsPtr_(),
344  curPointLabels
345  );
346 
347  procMeshPtr.reset
348  (
349  new polyMesh
350  (
351  IOobject
352  (
353  this->polyMesh::name(), // region of undecomposed mesh
354  facesInstance(),
355  processorDb
356  ),
357  move(facesInstancePoints),
358  move(procFaces),
359  move(procCells)
360  )
361  );
362  }
363  else
364  {
365  procMeshPtr.reset
366  (
367  new polyMesh
368  (
369  IOobject
370  (
371  this->polyMesh::name(), // region of undecomposed mesh
372  facesInstance(),
373  processorDb
374  ),
375  move(procPoints),
376  move(procFaces),
377  move(procCells)
378  )
379  );
380  }
381  polyMesh& procMesh = procMeshPtr();
382 
383 
384  // Create processor boundary patches
385  const labelList& curPatchSizes = procPatchSize_[proci];
386 
387  const labelList& curPatchStarts = procPatchStartIndex_[proci];
388 
389  const labelList& curNeighbourProcessors =
390  procNeighbourProcessors_[proci];
391 
392  const labelList& curProcessorPatchSizes =
393  procProcessorPatchSize_[proci];
394 
395  const labelList& curProcessorPatchStarts =
396  procProcessorPatchStartIndex_[proci];
397 
398  const labelListList& curSubPatchIDs =
399  procProcessorPatchSubPatchIDs_[proci];
400 
401  const labelListList& curSubStarts =
402  procProcessorPatchSubPatchStarts_[proci];
403 
404  const polyPatchList& meshPatches = boundaryMesh();
405 
406 
407  // Count the number of inter-proc patches
408  label nInterProcPatches = 0;
409  forAll(curSubPatchIDs, procPatchi)
410  {
411  nInterProcPatches += curSubPatchIDs[procPatchi].size();
412  }
413 
414  List<polyPatch*> procPatches
415  (
416  curPatchSizes.size() + nInterProcPatches,
417  reinterpret_cast<polyPatch*>(0)
418  );
419 
420  label nPatches = 0;
421 
422  forAll(curPatchSizes, patchi)
423  {
424  // Get the face labels consistent with the field mapping
425  // (reuse the patch field mappers)
426  const polyPatch& meshPatch = meshPatches[patchi];
427 
428  fvFieldDecomposer::patchFieldDecomposer patchMapper
429  (
430  SubList<label>
431  (
432  curFaceLabels,
433  curPatchSizes[patchi],
434  curPatchStarts[patchi]
435  ),
436  meshPatch.start()
437  );
438 
439  // Map existing patches
440  procPatches[nPatches] = meshPatch.clone
441  (
442  procMesh.boundaryMesh(),
443  nPatches,
444  patchMapper.addressing(),
445  curPatchStarts[patchi]
446  ).ptr();
447 
448  nPatches++;
449  }
450 
451  forAll(curProcessorPatchSizes, procPatchi)
452  {
453  const labelList& subPatchID = curSubPatchIDs[procPatchi];
454  const labelList& subStarts = curSubStarts[procPatchi];
455 
456  label curStart = curProcessorPatchStarts[procPatchi];
457 
458  forAll(subPatchID, i)
459  {
460  label size =
461  (
462  i < subPatchID.size()-1
463  ? subStarts[i+1] - subStarts[i]
464  : curProcessorPatchSizes[procPatchi] - subStarts[i]
465  );
466 
467  if (subPatchID[i] == -1)
468  {
469  // From internal faces
470  procPatches[nPatches] =
471  new processorPolyPatch
472  (
473  size,
474  curStart,
475  nPatches,
476  procMesh.boundaryMesh(),
477  proci,
478  curNeighbourProcessors[procPatchi]
479  );
480  }
481  else
482  {
483  const coupledPolyPatch& pcPatch
484  = refCast<const coupledPolyPatch>
485  (
486  boundaryMesh()[subPatchID[i]]
487  );
488 
489  procPatches[nPatches] =
490  new processorCyclicPolyPatch
491  (
492  size,
493  curStart,
494  nPatches,
495  procMesh.boundaryMesh(),
496  proci,
497  curNeighbourProcessors[procPatchi],
498  pcPatch.name()
499  );
500  }
501 
502  curStart += size;
503 
504  nPatches++;
505  }
506  }
507 
508  // Add boundary patches
509  procMesh.addPatches(procPatches);
510 
511  // Create and add zones
512 
513  // Point zones
514  {
515  const pointZoneMesh& pz = pointZones();
516 
517  // Go through all the zoned points and find out if they
518  // belong to a zone. If so, add it to the zone as
519  // necessary
520  List<DynamicList<label>> zonePoints(pz.size());
521 
522  // Estimate size
523  forAll(zonePoints, zoneI)
524  {
525  zonePoints[zoneI].setCapacity(pz[zoneI].size() / nProcs_);
526  }
527 
528  // Use the pointToZone map to find out the single zone (if any),
529  // use slow search only for shared points.
530  forAll(curPointLabels, pointi)
531  {
532  label curPoint = curPointLabels[pointi];
533 
534  label zoneI = pointToZone[curPoint];
535 
536  if (zoneI >= 0)
537  {
538  // Single zone.
539  zonePoints[zoneI].append(pointi);
540  }
541  else if (zoneI == -2)
542  {
543  // Multiple zones. Lookup.
544  forAll(pz, zoneI)
545  {
546  label index = pz[zoneI].whichPoint(curPoint);
547 
548  if (index != -1)
549  {
550  zonePoints[zoneI].append(pointi);
551  }
552  }
553  }
554  }
555 
556  procMesh.pointZones().clearAddressing();
557  procMesh.pointZones().setSize(zonePoints.size());
558  forAll(zonePoints, zoneI)
559  {
560  procMesh.pointZones().set
561  (
562  zoneI,
563  pz[zoneI].clone
564  (
565  procMesh.pointZones(),
566  zoneI,
567  zonePoints[zoneI].shrink()
568  )
569  );
570  }
571 
572  if (pz.size())
573  {
574  // Force writing on all processors
575  procMesh.pointZones().writeOpt() = IOobject::AUTO_WRITE;
576  }
577  }
578 
579  // Face zones
580  {
581  const faceZoneMesh& fz = faceZones();
582 
583  // Go through all the zoned face and find out if they
584  // belong to a zone. If so, add it to the zone as
585  // necessary
586  List<DynamicList<label>> zoneFaces(fz.size());
587  List<DynamicList<bool>> zoneFaceFlips(fz.size());
588 
589  // Estimate size
590  forAll(zoneFaces, zoneI)
591  {
592  label procSize = fz[zoneI].size() / nProcs_;
593 
594  zoneFaces[zoneI].setCapacity(procSize);
595  zoneFaceFlips[zoneI].setCapacity(procSize);
596  }
597 
598  // Go through all the zoned faces and find out if they
599  // belong to a zone. If so, add it to the zone as
600  // necessary
601  forAll(curFaceLabels, facei)
602  {
603  // Remember to decrement the index by one (turning index)
604  //
605  label curF = mag(curFaceLabels[facei]) - 1;
606 
607  label zoneI = faceToZone[curF];
608 
609  if (zoneI >= 0)
610  {
611  // Single zone. Add the face
612  zoneFaces[zoneI].append(facei);
613 
614  label index = fz[zoneI].whichFace(curF);
615 
616  bool flip = fz[zoneI].flipMap()[index];
617 
618  if (curFaceLabels[facei] < 0)
619  {
620  flip = !flip;
621  }
622 
623  zoneFaceFlips[zoneI].append(flip);
624  }
625  else if (zoneI == -2)
626  {
627  // Multiple zones. Lookup.
628  forAll(fz, zoneI)
629  {
630  label index = fz[zoneI].whichFace(curF);
631 
632  if (index != -1)
633  {
634  zoneFaces[zoneI].append(facei);
635 
636  bool flip = fz[zoneI].flipMap()[index];
637 
638  if (curFaceLabels[facei] < 0)
639  {
640  flip = !flip;
641  }
642 
643  zoneFaceFlips[zoneI].append(flip);
644  }
645  }
646  }
647  }
648 
649  procMesh.faceZones().clearAddressing();
650  procMesh.faceZones().setSize(zoneFaces.size());
651  forAll(zoneFaces, zoneI)
652  {
653  procMesh.faceZones().set
654  (
655  zoneI,
656  fz[zoneI].clone
657  (
658  zoneFaces[zoneI].shrink(), // addressing
659  zoneFaceFlips[zoneI].shrink(), // flipmap
660  zoneI,
661  procMesh.faceZones()
662  )
663  );
664  }
665 
666  if (fz.size())
667  {
668  // Force writing on all processors
669  procMesh.faceZones().writeOpt() = IOobject::AUTO_WRITE;
670  }
671  }
672 
673  // Cell zones
674  {
675  const cellZoneMesh& cz = cellZones();
676 
677  // Go through all the zoned cells and find out if they
678  // belong to a zone. If so, add it to the zone as
679  // necessary
680  List<DynamicList<label>> zoneCells(cz.size());
681 
682  // Estimate size
683  forAll(zoneCells, zoneI)
684  {
685  zoneCells[zoneI].setCapacity(cz[zoneI].size() / nProcs_);
686  }
687 
688  forAll(curCellLabels, celli)
689  {
690  label curCelli = curCellLabels[celli];
691 
692  label zoneI = cellToZone[curCelli];
693 
694  if (zoneI >= 0)
695  {
696  // Single zone.
697  zoneCells[zoneI].append(celli);
698  }
699  else if (zoneI == -2)
700  {
701  // Multiple zones. Lookup.
702  forAll(cz, zoneI)
703  {
704  label index = cz[zoneI].whichCell(curCelli);
705 
706  if (index != -1)
707  {
708  zoneCells[zoneI].append(celli);
709  }
710  }
711  }
712  }
713 
714  procMesh.cellZones().clearAddressing();
715  procMesh.cellZones().setSize(zoneCells.size());
716  forAll(zoneCells, zoneI)
717  {
718  procMesh.cellZones().set
719  (
720  zoneI,
721  cz[zoneI].clone
722  (
723  zoneCells[zoneI].shrink(),
724  zoneI,
725  procMesh.cellZones()
726  )
727  );
728  }
729 
730  if (cz.size())
731  {
732  // Force writing on all processors
733  procMesh.cellZones().writeOpt() = IOobject::AUTO_WRITE;
734  }
735  }
736 
737  // Set the precision of the points data to be min 10
739 
740  procMesh.write();
741 
742  // Write points if pointsInstance differing from facesInstance
743  if (facesInstancePointsPtr_.valid())
744  {
745  pointIOField pointsInstancePoints
746  (
747  IOobject
748  (
749  "points",
750  pointsInstance(),
752  procMesh,
755  false
756  ),
757  move(procPoints)
758  );
759  pointsInstancePoints.write();
760  }
761 
762 
763  // Decompose any sets
764  if (decomposeSets)
765  {
766  forAll(cellSets, i)
767  {
768  const cellSet& cs = cellSets[i];
769  cellSet set(procMesh, cs.name(), cs.size()/nProcs_);
770  forAll(curCellLabels, i)
771  {
772  if (cs.found(curCellLabels[i]))
773  {
774  set.insert(i);
775  }
776  }
777  set.write();
778  }
779  forAll(faceSets, i)
780  {
781  const faceSet& cs = faceSets[i];
782  faceSet set(procMesh, cs.name(), cs.size()/nProcs_);
783  forAll(curFaceLabels, i)
784  {
785  if (cs.found(mag(curFaceLabels[i])-1))
786  {
787  set.insert(i);
788  }
789  }
790  set.write();
791  }
792  forAll(pointSets, i)
793  {
794  const pointSet& cs = pointSets[i];
795  pointSet set(procMesh, cs.name(), cs.size()/nProcs_);
796  forAll(curPointLabels, i)
797  {
798  if (cs.found(curPointLabels[i]))
799  {
800  set.insert(i);
801  }
802  }
803  set.write();
804  }
805  }
806 
807 
808  // Optional hexRef8 data
809  hexRef8Data
810  (
811  IOobject
812  (
813  "dummy",
814  facesInstance(),
816  procMesh,
819  false
820  ),
821  baseMeshData,
822  procCellAddressing_[proci],
823  procPointAddressing_[proci]
824  ).write();
825 
826 
827  // Statistics
828 
829  Info<< endl
830  << "Processor " << proci << nl
831  << " Number of cells = " << procMesh.nCells()
832  << endl;
833 
834  maxProcCells = max(maxProcCells, procMesh.nCells());
835 
836  label nBoundaryFaces = 0;
837  label nProcPatches = 0;
838  label nProcFaces = 0;
839 
840  forAll(procMesh.boundaryMesh(), patchi)
841  {
842  if (isA<processorPolyPatch>(procMesh.boundaryMesh()[patchi]))
843  {
844  const processorPolyPatch& ppp =
845  refCast<const processorPolyPatch>
846  (
847  procMesh.boundaryMesh()[patchi]
848  );
849 
850  Info<< " Number of faces shared with processor "
851  << ppp.neighbProcNo() << " = " << ppp.size() << endl;
852 
853  nProcPatches++;
854  nProcFaces += ppp.size();
855  }
856  else
857  {
858  nBoundaryFaces += procMesh.boundaryMesh()[patchi].size();
859  }
860  }
861 
862  Info<< " Number of processor patches = " << nProcPatches << nl
863  << " Number of processor faces = " << nProcFaces << nl
864  << " Number of boundary faces = " << nBoundaryFaces << endl;
865 
866  totProcFaces += nProcFaces;
867  totProcPatches += nProcPatches;
868  maxProcPatches = max(maxProcPatches, nProcPatches);
869  maxProcFaces = max(maxProcFaces, nProcFaces);
870 
871  // create and write the addressing information
872  labelIOList pointProcAddressing
873  (
874  IOobject
875  (
876  "pointProcAddressing",
877  procMesh.facesInstance(),
878  procMesh.meshSubDir,
879  procMesh,
882  ),
883  procPointAddressing_[proci]
884  );
885  pointProcAddressing.write();
886 
888  (
889  IOobject
890  (
891  "faceProcAddressing",
892  procMesh.facesInstance(),
893  procMesh.meshSubDir,
894  procMesh,
897  ),
898  procFaceAddressing_[proci]
899  );
900  faceProcAddressing.write();
901 
902  labelIOList cellProcAddressing
903  (
904  IOobject
905  (
906  "cellProcAddressing",
907  procMesh.facesInstance(),
908  procMesh.meshSubDir,
909  procMesh,
912  ),
913  procCellAddressing_[proci]
914  );
915  cellProcAddressing.write();
916 
917  // Write patch map for backwards compatibility.
918  // (= identity map for original patches, -1 for processor patches)
919  label nMeshPatches = curPatchSizes.size();
920  labelList procBoundaryAddressing(identity(nMeshPatches));
921  procBoundaryAddressing.setSize(nMeshPatches+nProcPatches, -1);
922 
923  labelIOList boundaryProcAddressing
924  (
925  IOobject
926  (
927  "boundaryProcAddressing",
928  procMesh.facesInstance(),
929  procMesh.meshSubDir,
930  procMesh,
933  ),
934  procBoundaryAddressing
935  );
936  boundaryProcAddressing.write();
937  }
938 
939  scalar avgProcCells = scalar(nCells())/nProcs_;
940  scalar avgProcPatches = scalar(totProcPatches)/nProcs_;
941  scalar avgProcFaces = scalar(totProcFaces)/nProcs_;
942 
943  // In case of all faces on one processor. Just to avoid division by 0.
944  if (totProcPatches == 0)
945  {
946  avgProcPatches = 1;
947  }
948  if (totProcFaces == 0)
949  {
950  avgProcFaces = 1;
951  }
952 
953  Info<< nl
954  << "Number of processor faces = " << totProcFaces/2 << nl
955  << "Max number of cells = " << maxProcCells
956  << " (" << 100.0*(maxProcCells-avgProcCells)/avgProcCells
957  << "% above average " << avgProcCells << ")" << nl
958  << "Max number of processor patches = " << maxProcPatches
959  << " (" << 100.0*(maxProcPatches-avgProcPatches)/avgProcPatches
960  << "% above average " << avgProcPatches << ")" << nl
961  << "Max number of faces between processors = " << maxProcFaces
962  << " (" << 100.0*(maxProcFaces-avgProcFaces)/avgProcFaces
963  << "% above average " << avgProcFaces << ")" << nl
964  << endl;
965 
966  return true;
967 }
968 
969 
970 // ************************************************************************* //
label nPatches
Definition: readKivaGrid.H:402
PtrList< polyPatch > polyPatchList
container classes for polyPatch
Definition: polyPatchList.H:45
List< labelList > labelListList
A List of labelList.
Definition: labelList.H:57
bool writeDecomposition(const bool decomposeSets)
Write decomposition.
ZoneMesh< faceZone, polyMesh > faceZoneMesh
A ZoneMesh with the type faceZone.
#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
const word & name() const
Return name.
Definition: IOobject.H:303
static const decompositionModel & New(const polyMesh &mesh, const fileName &decompDictFile="")
Read (optionally from absolute path) & register on mesh.
domainDecomposition(const IOobject &io, const fileName &dictFile)
Construct from IOobject and decomposition dictionary name.
const fileName & facesInstance() const
Return the current instance directory for faces.
Definition: polyMesh.C:808
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
vectorIOField pointIOField
pointIOField is a vectorIOField.
Definition: pointIOField.H:42
PtrList< labelIOList > & faceProcAddressing
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
Definition: polyMesh.H:312
static unsigned int defaultPrecision()
Return the default precision.
Definition: IOstream.H:461
IOobject(const word &name, const fileName &instance, const objectRegistry &registry, readOption r=NO_READ, writeOption w=NO_WRITE, bool registerObject=true)
Construct from name, instance, registry, io options.
Definition: IOobject.C:177
List< face > faceList
Definition: faceListFwd.H:43
label nCells() const
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
Definition: ListOps.C:104
~domainDecomposition()
Destructor.
Functions used by OpenFOAM that are specific to POSIX compliant operating systems and need to be repl...
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
T clone(const T &t)
Definition: List.H:54
const fileName & pointsInstance() const
Return the current instance directory for points.
Definition: polyMesh.C:802
const cellShapeList & cells
const pointField & points
label nPoints
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
List< label > labelList
A List of labels.
Definition: labelList.H:56
ZoneMesh< pointZone, polyMesh > pointZoneMesh
A ZoneMesh with the type pointZone.
static word controlDictName
The default control dictionary name (normally "controlDict")
Definition: Time.H:197
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
static const char nl
Definition: Ostream.H:260
objects
static autoPtr< dictionary > New(Istream &)
Construct top-level dictionary on freestore from Istream.
Definition: dictionaryIO.C:69
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
fvMesh(const IOobject &io)
Construct from IOobject.
Definition: fvMesh.C:256
void setSize(const label)
Reset size of List.
Definition: List.C:281
label patchi
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
List< cell > cellList
list of cells
Definition: cellList.H:42
ZoneMesh< cellZone, polyMesh > cellZoneMesh
A ZoneMesh with the type cellZone.
IOList< label > labelIOList
Label container classes.
Definition: labelIOList.H:42
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:812
polyMesh(const IOobject &io)
Construct from IOobject.
Definition: polyMesh.C:163