domainDecomposition.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 | Copyright (C) 2011-2016 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  fvMesh(io),
75  facesInstancePointsPtr_
76  (
78  ? new pointIOField
79  (
80  IOobject
81  (
82  "points",
83  facesInstance(),
85  *this,
88  false
89  )
90  )
91  : NULL
92  ),
93  nProcs_
94  (
95  readInt
96  (
97  decompositionModel::New
98  (
99  *this
100  ).lookup("numberOfSubdomains")
101  )
102  ),
103  distributed_(false),
104  cellToProc_(nCells()),
105  procPointAddressing_(nProcs_),
106  procFaceAddressing_(nProcs_),
107  procCellAddressing_(nProcs_),
108  procPatchSize_(nProcs_),
109  procPatchStartIndex_(nProcs_),
110  procNeighbourProcessors_(nProcs_),
111  procProcessorPatchSize_(nProcs_),
112  procProcessorPatchStartIndex_(nProcs_),
113  procProcessorPatchSubPatchIDs_(nProcs_),
114  procProcessorPatchSubPatchStarts_(nProcs_)
115 {
117  (
118  *this
119  ).readIfPresent("distributed", distributed_);
120 }
121 
122 
123 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
124 
126 {}
127 
128 
129 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
130 
131 bool Foam::domainDecomposition::writeDecomposition(const bool decomposeSets)
132 {
133  Info<< "\nConstructing processor meshes" << endl;
134 
135  // Mark point/faces/cells that are in zones.
136  // -1 : not in zone
137  // -2 : in multiple zones
138  // >= 0 : in single given zone
139  // This will give direct lookup of elements that are in a single zone
140  // and we'll only have to revert back to searching through all zones
141  // for the duplicate elements
142 
143  // Point zones
144  labelList pointToZone(points().size(), -1);
145 
146  forAll(pointZones(), zoneI)
147  {
148  mark(pointZones()[zoneI], zoneI, pointToZone);
149  }
150 
151  // Face zones
152  labelList faceToZone(faces().size(), -1);
153 
154  forAll(faceZones(), zoneI)
155  {
156  mark(faceZones()[zoneI], zoneI, faceToZone);
157  }
158 
159  // Cell zones
160  labelList cellToZone(nCells(), -1);
161 
162  forAll(cellZones(), zoneI)
163  {
164  mark(cellZones()[zoneI], zoneI, cellToZone);
165  }
166 
167 
168  PtrList<const cellSet> cellSets;
169  PtrList<const faceSet> faceSets;
170  PtrList<const pointSet> pointSets;
171  if (decomposeSets)
172  {
173  // Read sets
174  IOobjectList objects(*this, facesInstance(), "polyMesh/sets");
175  {
176  IOobjectList cSets(objects.lookupClass(cellSet::typeName));
177  forAllConstIter(IOobjectList, cSets, iter)
178  {
179  cellSets.append(new cellSet(*iter()));
180  }
181  }
182  {
183  IOobjectList fSets(objects.lookupClass(faceSet::typeName));
184  forAllConstIter(IOobjectList, fSets, iter)
185  {
186  faceSets.append(new faceSet(*iter()));
187  }
188  }
189  {
190  IOobjectList pSets(objects.lookupClass(pointSet::typeName));
191  forAllConstIter(IOobjectList, pSets, iter)
192  {
193  pointSets.append(new pointSet(*iter()));
194  }
195  }
196  }
197 
198 
199  // Load refinement data (if any)
200  hexRef8Data baseMeshData
201  (
202  IOobject
203  (
204  "dummy",
205  facesInstance(),
207  *this,
210  false
211  )
212  );
213 
214 
215 
216  label maxProcCells = 0;
217  label totProcFaces = 0;
218  label maxProcPatches = 0;
219  label totProcPatches = 0;
220  label maxProcFaces = 0;
221 
222 
223  // Write out the meshes
224  for (label proci = 0; proci < nProcs_; proci++)
225  {
226  // Create processor points
227  const labelList& curPointLabels = procPointAddressing_[proci];
228 
229  const pointField& meshPoints = points();
230 
231  labelList pointLookup(nPoints(), -1);
232 
233  pointField procPoints(curPointLabels.size());
234 
235  forAll(curPointLabels, pointi)
236  {
237  procPoints[pointi] = meshPoints[curPointLabels[pointi]];
238 
239  pointLookup[curPointLabels[pointi]] = pointi;
240  }
241 
242  // Create processor faces
243  const labelList& curFaceLabels = procFaceAddressing_[proci];
244 
245  const faceList& meshFaces = faces();
246 
247  labelList faceLookup(nFaces(), -1);
248 
249  faceList procFaces(curFaceLabels.size());
250 
251  forAll(curFaceLabels, facei)
252  {
253  // Mark the original face as used
254  // Remember to decrement the index by one (turning index)
255  //
256  label curF = mag(curFaceLabels[facei]) - 1;
257 
258  faceLookup[curF] = facei;
259 
260  // get the original face
261  labelList origFaceLabels;
262 
263  if (curFaceLabels[facei] >= 0)
264  {
265  // face not turned
266  origFaceLabels = meshFaces[curF];
267  }
268  else
269  {
270  origFaceLabels = meshFaces[curF].reverseFace();
271  }
272 
273  // translate face labels into local point list
274  face& procFaceLabels = procFaces[facei];
275 
276  procFaceLabels.setSize(origFaceLabels.size());
277 
278  forAll(origFaceLabels, pointi)
279  {
280  procFaceLabels[pointi] = pointLookup[origFaceLabels[pointi]];
281  }
282  }
283 
284  // Create processor cells
285  const labelList& curCellLabels = procCellAddressing_[proci];
286 
287  const cellList& meshCells = cells();
288 
289  cellList procCells(curCellLabels.size());
290 
291  forAll(curCellLabels, celli)
292  {
293  const labelList& origCellLabels = meshCells[curCellLabels[celli]];
294 
295  cell& curCell = procCells[celli];
296 
297  curCell.setSize(origCellLabels.size());
298 
299  forAll(origCellLabels, cellFacei)
300  {
301  curCell[cellFacei] = faceLookup[origCellLabels[cellFacei]];
302  }
303  }
304 
305  // Create processor mesh without a boundary
306 
307  fileName processorCasePath
308  (
309  time().caseName()/fileName(word("processor") + Foam::name(proci))
310  );
311 
312  // make the processor directory
313  mkDir(time().rootPath()/processorCasePath);
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  xferMove(facesInstancePoints),
358  xferMove(procFaces),
359  xferMove(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  xferMove(procPoints),
376  xferMove(procFaces),
377  xferMove(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.directAddressing(),
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  pcPatch.transform()
500  );
501  }
502 
503  curStart += size;
504 
505  nPatches++;
506  }
507  }
508 
509  // Add boundary patches
510  procMesh.addPatches(procPatches);
511 
512  // Create and add zones
513 
514  // Point zones
515  {
516  const pointZoneMesh& pz = pointZones();
517 
518  // Go through all the zoned points and find out if they
519  // belong to a zone. If so, add it to the zone as
520  // necessary
521  List<DynamicList<label>> zonePoints(pz.size());
522 
523  // Estimate size
524  forAll(zonePoints, zoneI)
525  {
526  zonePoints[zoneI].setCapacity(pz[zoneI].size() / nProcs_);
527  }
528 
529  // Use the pointToZone map to find out the single zone (if any),
530  // use slow search only for shared points.
531  forAll(curPointLabels, pointi)
532  {
533  label curPoint = curPointLabels[pointi];
534 
535  label zoneI = pointToZone[curPoint];
536 
537  if (zoneI >= 0)
538  {
539  // Single zone.
540  zonePoints[zoneI].append(pointi);
541  }
542  else if (zoneI == -2)
543  {
544  // Multiple zones. Lookup.
545  forAll(pz, zoneI)
546  {
547  label index = pz[zoneI].whichPoint(curPoint);
548 
549  if (index != -1)
550  {
551  zonePoints[zoneI].append(pointi);
552  }
553  }
554  }
555  }
556 
557  procMesh.pointZones().clearAddressing();
558  procMesh.pointZones().setSize(zonePoints.size());
559  forAll(zonePoints, zoneI)
560  {
561  procMesh.pointZones().set
562  (
563  zoneI,
564  pz[zoneI].clone
565  (
566  procMesh.pointZones(),
567  zoneI,
568  zonePoints[zoneI].shrink()
569  )
570  );
571  }
572 
573  if (pz.size())
574  {
575  // Force writing on all processors
576  procMesh.pointZones().writeOpt() = IOobject::AUTO_WRITE;
577  }
578  }
579 
580  // Face zones
581  {
582  const faceZoneMesh& fz = faceZones();
583 
584  // Go through all the zoned face and find out if they
585  // belong to a zone. If so, add it to the zone as
586  // necessary
587  List<DynamicList<label>> zoneFaces(fz.size());
588  List<DynamicList<bool>> zoneFaceFlips(fz.size());
589 
590  // Estimate size
591  forAll(zoneFaces, zoneI)
592  {
593  label procSize = fz[zoneI].size() / nProcs_;
594 
595  zoneFaces[zoneI].setCapacity(procSize);
596  zoneFaceFlips[zoneI].setCapacity(procSize);
597  }
598 
599  // Go through all the zoned faces and find out if they
600  // belong to a zone. If so, add it to the zone as
601  // necessary
602  forAll(curFaceLabels, facei)
603  {
604  // Remember to decrement the index by one (turning index)
605  //
606  label curF = mag(curFaceLabels[facei]) - 1;
607 
608  label zoneI = faceToZone[curF];
609 
610  if (zoneI >= 0)
611  {
612  // Single zone. Add the face
613  zoneFaces[zoneI].append(facei);
614 
615  label index = fz[zoneI].whichFace(curF);
616 
617  bool flip = fz[zoneI].flipMap()[index];
618 
619  if (curFaceLabels[facei] < 0)
620  {
621  flip = !flip;
622  }
623 
624  zoneFaceFlips[zoneI].append(flip);
625  }
626  else if (zoneI == -2)
627  {
628  // Multiple zones. Lookup.
629  forAll(fz, zoneI)
630  {
631  label index = fz[zoneI].whichFace(curF);
632 
633  if (index != -1)
634  {
635  zoneFaces[zoneI].append(facei);
636 
637  bool flip = fz[zoneI].flipMap()[index];
638 
639  if (curFaceLabels[facei] < 0)
640  {
641  flip = !flip;
642  }
643 
644  zoneFaceFlips[zoneI].append(flip);
645  }
646  }
647  }
648  }
649 
650  procMesh.faceZones().clearAddressing();
651  procMesh.faceZones().setSize(zoneFaces.size());
652  forAll(zoneFaces, zoneI)
653  {
654  procMesh.faceZones().set
655  (
656  zoneI,
657  fz[zoneI].clone
658  (
659  zoneFaces[zoneI].shrink(), // addressing
660  zoneFaceFlips[zoneI].shrink(), // flipmap
661  zoneI,
662  procMesh.faceZones()
663  )
664  );
665  }
666 
667  if (fz.size())
668  {
669  // Force writing on all processors
670  procMesh.faceZones().writeOpt() = IOobject::AUTO_WRITE;
671  }
672  }
673 
674  // Cell zones
675  {
676  const cellZoneMesh& cz = cellZones();
677 
678  // Go through all the zoned cells and find out if they
679  // belong to a zone. If so, add it to the zone as
680  // necessary
681  List<DynamicList<label>> zoneCells(cz.size());
682 
683  // Estimate size
684  forAll(zoneCells, zoneI)
685  {
686  zoneCells[zoneI].setCapacity(cz[zoneI].size() / nProcs_);
687  }
688 
689  forAll(curCellLabels, celli)
690  {
691  label curCelli = curCellLabels[celli];
692 
693  label zoneI = cellToZone[curCelli];
694 
695  if (zoneI >= 0)
696  {
697  // Single zone.
698  zoneCells[zoneI].append(celli);
699  }
700  else if (zoneI == -2)
701  {
702  // Multiple zones. Lookup.
703  forAll(cz, zoneI)
704  {
705  label index = cz[zoneI].whichCell(curCelli);
706 
707  if (index != -1)
708  {
709  zoneCells[zoneI].append(celli);
710  }
711  }
712  }
713  }
714 
715  procMesh.cellZones().clearAddressing();
716  procMesh.cellZones().setSize(zoneCells.size());
717  forAll(zoneCells, zoneI)
718  {
719  procMesh.cellZones().set
720  (
721  zoneI,
722  cz[zoneI].clone
723  (
724  zoneCells[zoneI].shrink(),
725  zoneI,
726  procMesh.cellZones()
727  )
728  );
729  }
730 
731  if (cz.size())
732  {
733  // Force writing on all processors
734  procMesh.cellZones().writeOpt() = IOobject::AUTO_WRITE;
735  }
736  }
737 
738  // Set the precision of the points data to be min 10
740 
741  procMesh.write();
742 
743  // Write points if pointsInstance differing from facesInstance
744  if (facesInstancePointsPtr_.valid())
745  {
746  pointIOField pointsInstancePoints
747  (
748  IOobject
749  (
750  "points",
751  pointsInstance(),
753  procMesh,
756  false
757  ),
758  xferMove(procPoints)
759  );
760  pointsInstancePoints.write();
761  }
762 
763 
764  // Decompose any sets
765  if (decomposeSets)
766  {
767  forAll(cellSets, i)
768  {
769  const cellSet& cs = cellSets[i];
770  cellSet set(procMesh, cs.name(), cs.size()/nProcs_);
771  forAll(curCellLabels, i)
772  {
773  if (cs.found(curCellLabels[i]))
774  {
775  set.insert(i);
776  }
777  }
778  set.write();
779  }
780  forAll(faceSets, i)
781  {
782  const faceSet& cs = faceSets[i];
783  faceSet set(procMesh, cs.name(), cs.size()/nProcs_);
784  forAll(curFaceLabels, i)
785  {
786  if (cs.found(mag(curFaceLabels[i])-1))
787  {
788  set.insert(i);
789  }
790  }
791  set.write();
792  }
793  forAll(pointSets, i)
794  {
795  const pointSet& cs = pointSets[i];
796  pointSet set(procMesh, cs.name(), cs.size()/nProcs_);
797  forAll(curPointLabels, i)
798  {
799  if (cs.found(curPointLabels[i]))
800  {
801  set.insert(i);
802  }
803  }
804  set.write();
805  }
806  }
807 
808 
809  // Optional hexRef8 data
810  hexRef8Data
811  (
812  IOobject
813  (
814  "dummy",
815  facesInstance(),
817  procMesh,
820  false
821  ),
822  baseMeshData,
823  procCellAddressing_[proci],
824  procPointAddressing_[proci]
825  ).write();
826 
827 
828  // Statistics
829 
830  Info<< endl
831  << "Processor " << proci << nl
832  << " Number of cells = " << procMesh.nCells()
833  << endl;
834 
835  maxProcCells = max(maxProcCells, procMesh.nCells());
836 
837  label nBoundaryFaces = 0;
838  label nProcPatches = 0;
839  label nProcFaces = 0;
840 
841  forAll(procMesh.boundaryMesh(), patchi)
842  {
843  if (isA<processorPolyPatch>(procMesh.boundaryMesh()[patchi]))
844  {
845  const processorPolyPatch& ppp =
846  refCast<const processorPolyPatch>
847  (
848  procMesh.boundaryMesh()[patchi]
849  );
850 
851  Info<< " Number of faces shared with processor "
852  << ppp.neighbProcNo() << " = " << ppp.size() << endl;
853 
854  nProcPatches++;
855  nProcFaces += ppp.size();
856  }
857  else
858  {
859  nBoundaryFaces += procMesh.boundaryMesh()[patchi].size();
860  }
861  }
862 
863  Info<< " Number of processor patches = " << nProcPatches << nl
864  << " Number of processor faces = " << nProcFaces << nl
865  << " Number of boundary faces = " << nBoundaryFaces << endl;
866 
867  totProcFaces += nProcFaces;
868  totProcPatches += nProcPatches;
869  maxProcPatches = max(maxProcPatches, nProcPatches);
870  maxProcFaces = max(maxProcFaces, nProcFaces);
871 
872  // create and write the addressing information
873  labelIOList pointProcAddressing
874  (
875  IOobject
876  (
877  "pointProcAddressing",
878  procMesh.facesInstance(),
879  procMesh.meshSubDir,
880  procMesh,
883  ),
884  procPointAddressing_[proci]
885  );
886  pointProcAddressing.write();
887 
889  (
890  IOobject
891  (
892  "faceProcAddressing",
893  procMesh.facesInstance(),
894  procMesh.meshSubDir,
895  procMesh,
898  ),
899  procFaceAddressing_[proci]
900  );
901  faceProcAddressing.write();
902 
903  labelIOList cellProcAddressing
904  (
905  IOobject
906  (
907  "cellProcAddressing",
908  procMesh.facesInstance(),
909  procMesh.meshSubDir,
910  procMesh,
913  ),
914  procCellAddressing_[proci]
915  );
916  cellProcAddressing.write();
917 
918  // Write patch map for backwards compatibility.
919  // (= identity map for original patches, -1 for processor patches)
920  label nMeshPatches = curPatchSizes.size();
921  labelList procBoundaryAddressing(identity(nMeshPatches));
922  procBoundaryAddressing.setSize(nMeshPatches+nProcPatches, -1);
923 
924  labelIOList boundaryProcAddressing
925  (
926  IOobject
927  (
928  "boundaryProcAddressing",
929  procMesh.facesInstance(),
930  procMesh.meshSubDir,
931  procMesh,
934  ),
935  procBoundaryAddressing
936  );
937  boundaryProcAddressing.write();
938  }
939 
940  scalar avgProcCells = scalar(nCells())/nProcs_;
941  scalar avgProcPatches = scalar(totProcPatches)/nProcs_;
942  scalar avgProcFaces = scalar(totProcFaces)/nProcs_;
943 
944  // In case of all faces on one processor. Just to avoid division by 0.
945  if (totProcPatches == 0)
946  {
947  avgProcPatches = 1;
948  }
949  if (totProcFaces == 0)
950  {
951  avgProcFaces = 1;
952  }
953 
954  Info<< nl
955  << "Number of processor faces = " << totProcFaces/2 << nl
956  << "Max number of cells = " << maxProcCells
957  << " (" << 100.0*(maxProcCells-avgProcCells)/avgProcCells
958  << "% above average " << avgProcCells << ")" << nl
959  << "Max number of processor patches = " << maxProcPatches
960  << " (" << 100.0*(maxProcPatches-avgProcPatches)/avgProcPatches
961  << "% above average " << avgProcPatches << ")" << nl
962  << "Max number of faces between processors = " << maxProcFaces
963  << " (" << 100.0*(maxProcFaces-avgProcFaces)/avgProcFaces
964  << "% above average " << avgProcFaces << ")" << nl
965  << endl;
966 
967  return true;
968 }
969 
970 
971 // ************************************************************************* //
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:428
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
static const decompositionModel & New(const polyMesh &mesh, const fileName &decompDictFile="")
Read (optionallly from absolute path) & register on mesh.
const fileName & facesInstance() const
Return the current instance directory for faces.
Definition: polyMesh.C:778
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
int readInt(Istream &)
Definition: intIO.C:31
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:309
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:114
List< face > faceList
Definition: faceListFwd.H:43
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
domainDecomposition(const IOobject &io)
Construct from IOobject.
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
Definition: ListOps.C:104
~domainDecomposition()
Destructor.
label nCells() const
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
const cellShapeList & cells
const pointField & points
Xfer< T > xferMove(T &)
Construct by transferring the contents of the arg.
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:201
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
static const char nl
Definition: Ostream.H:262
static autoPtr< dictionary > New(Istream &)
Construct top-level dictionary on freestore from Istream.
Definition: dictionaryIO.C:73
bool mkDir(const fileName &, mode_t=0777)
Make a directory and return an error if it could not be created.
Definition: POSIX.C:295
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:246
void setSize(const label)
Reset size of List.
Definition: List.C:295
label patchi
const fileName & pointsInstance() const
Return the current instance directory for points.
Definition: polyMesh.C:772
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.
const word & name() const
Return name.
Definition: IOobject.H:260
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:451
polyMesh(const IOobject &io)
Construct from IOobject.
Definition: polyMesh.C:140