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