45 void Foam::domainDecomposition::mark
54 label pointi = zoneElems[i];
56 if (elementToZone[pointi] == -1)
59 elementToZone[pointi] = zoneI;
61 else if (elementToZone[pointi] >= 0)
64 elementToZone[pointi] = -2;
75 const fileName& dictFile
79 facesInstancePointsPtr_
99 decompositionModel::
New 103 ).
lookup<int>(
"numberOfSubdomains")
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_)
136 Info<<
"\nConstructing processor meshes" <<
endl;
149 forAll(pointZones(), zoneI)
151 mark(pointZones()[zoneI], zoneI, pointToZone);
155 labelList faceToZone(faces().size(), -1);
157 forAll(faceZones(), zoneI)
159 mark(faceZones()[zoneI], zoneI, faceToZone);
165 forAll(cellZones(), zoneI)
167 mark(cellZones()[zoneI], zoneI, cellToZone);
171 PtrList<const cellSet> cellSets;
172 PtrList<const faceSet> faceSets;
173 PtrList<const pointSet> pointSets;
177 IOobjectList
objects(*
this, facesInstance(),
"polyMesh/sets");
179 IOobjectList cSets(
objects.lookupClass(cellSet::typeName));
182 cellSets.append(
new cellSet(*iter()));
186 IOobjectList fSets(
objects.lookupClass(faceSet::typeName));
189 faceSets.append(
new faceSet(*iter()));
193 IOobjectList pSets(
objects.lookupClass(pointSet::typeName));
196 pointSets.append(
new pointSet(*iter()));
203 hexRef8Data baseMeshData
219 label maxProcCells = 0;
220 label totProcFaces = 0;
221 label maxProcPatches = 0;
222 label totProcPatches = 0;
223 label maxProcFaces = 0;
227 for (
label proci = 0; proci < nProcs_; proci++)
230 const labelList& curPointLabels = procPointAddressing_[proci];
238 forAll(curPointLabels, pointi)
240 procPoints[pointi] = meshPoints[curPointLabels[pointi]];
242 pointLookup[curPointLabels[pointi]] = pointi;
246 const labelList& curFaceLabels = procFaceAddressing_[proci];
248 const faceList& meshFaces = faces();
252 faceList procFaces(curFaceLabels.size());
254 forAll(curFaceLabels, facei)
259 label curF =
mag(curFaceLabels[facei]) - 1;
261 faceLookup[curF] = facei;
266 if (curFaceLabels[facei] >= 0)
269 origFaceLabels = meshFaces[curF];
273 origFaceLabels = meshFaces[curF].reverseFace();
277 face& procFaceLabels = procFaces[facei];
279 procFaceLabels.
setSize(origFaceLabels.size());
281 forAll(origFaceLabels, pointi)
283 procFaceLabels[pointi] = pointLookup[origFaceLabels[pointi]];
288 const labelList& curCellLabels = procCellAddressing_[proci];
292 cellList procCells(curCellLabels.size());
294 forAll(curCellLabels, celli)
296 const labelList& origCellLabels = meshCells[curCellLabels[celli]];
298 cell& curCell = procCells[celli];
300 curCell.
setSize(origCellLabels.size());
302 forAll(origCellLabels, cellFacei)
304 curCell[cellFacei] = faceLookup[origCellLabels[cellFacei]];
310 fileName processorCasePath
312 time().caseName()/fileName(word(
"processor") +
Foam::name(proci))
324 processorDb.setTime(time());
336 autoPtr<polyMesh> procMeshPtr;
338 if (facesInstancePointsPtr_.valid())
343 facesInstancePointsPtr_(),
357 move(facesInstancePoints),
381 polyMesh& procMesh = procMeshPtr();
385 const labelList& curPatchSizes = procPatchSize_[proci];
387 const labelList& curPatchStarts = procPatchStartIndex_[proci];
389 const labelList& curNeighbourProcessors =
390 procNeighbourProcessors_[proci];
392 const labelList& curProcessorPatchSizes =
393 procProcessorPatchSize_[proci];
395 const labelList& curProcessorPatchStarts =
396 procProcessorPatchStartIndex_[proci];
399 procProcessorPatchSubPatchIDs_[proci];
402 procProcessorPatchSubPatchStarts_[proci];
408 label nInterProcPatches = 0;
409 forAll(curSubPatchIDs, procPatchi)
411 nInterProcPatches += curSubPatchIDs[procPatchi].size();
414 List<polyPatch*> procPatches
416 curPatchSizes.size() + nInterProcPatches,
417 reinterpret_cast<polyPatch*
>(0)
422 forAll(curPatchSizes, patchi)
426 const polyPatch& meshPatch = meshPatches[
patchi];
428 fvFieldDecomposer::patchFieldDecomposer patchMapper
433 curPatchSizes[patchi],
434 curPatchStarts[patchi]
440 procPatches[
nPatches] = meshPatch.clone
442 procMesh.boundaryMesh(),
444 patchMapper.addressing(),
451 forAll(curProcessorPatchSizes, procPatchi)
453 const labelList& subPatchID = curSubPatchIDs[procPatchi];
454 const labelList& subStarts = curSubStarts[procPatchi];
456 label curStart = curProcessorPatchStarts[procPatchi];
462 i < subPatchID.size()-1
463 ? subStarts[i+1] - subStarts[i]
464 : curProcessorPatchSizes[procPatchi] - subStarts[i]
467 if (subPatchID[i] == -1)
471 new processorPolyPatch
476 procMesh.boundaryMesh(),
478 curNeighbourProcessors[procPatchi]
483 const coupledPolyPatch& pcPatch
484 = refCast<const coupledPolyPatch>
486 boundaryMesh()[subPatchID[i]]
490 new processorCyclicPolyPatch
495 procMesh.boundaryMesh(),
497 curNeighbourProcessors[procPatchi],
509 procMesh.addPatches(procPatches);
520 List<DynamicList<label>> zonePoints(pz.size());
525 zonePoints[zoneI].setCapacity(pz[zoneI].size() / nProcs_);
530 forAll(curPointLabels, pointi)
532 label curPoint = curPointLabels[pointi];
534 label zoneI = pointToZone[curPoint];
539 zonePoints[zoneI].append(pointi);
541 else if (zoneI == -2)
546 label index = pz[zoneI].whichPoint(curPoint);
550 zonePoints[zoneI].append(pointi);
556 procMesh.pointZones().clearAddressing();
557 procMesh.pointZones().setSize(zonePoints.size());
560 procMesh.pointZones().set
565 procMesh.pointZones(),
567 zonePoints[zoneI].shrink()
586 List<DynamicList<label>> zoneFaces(fz.size());
587 List<DynamicList<bool>> zoneFaceFlips(fz.size());
592 label procSize = fz[zoneI].size() / nProcs_;
594 zoneFaces[zoneI].setCapacity(procSize);
595 zoneFaceFlips[zoneI].setCapacity(procSize);
601 forAll(curFaceLabels, facei)
605 label curF =
mag(curFaceLabels[facei]) - 1;
607 label zoneI = faceToZone[curF];
612 zoneFaces[zoneI].append(facei);
614 label index = fz[zoneI].whichFace(curF);
616 bool flip = fz[zoneI].flipMap()[index];
618 if (curFaceLabels[facei] < 0)
623 zoneFaceFlips[zoneI].append(flip);
625 else if (zoneI == -2)
630 label index = fz[zoneI].whichFace(curF);
634 zoneFaces[zoneI].append(facei);
636 bool flip = fz[zoneI].flipMap()[index];
638 if (curFaceLabels[facei] < 0)
643 zoneFaceFlips[zoneI].append(flip);
649 procMesh.faceZones().clearAddressing();
650 procMesh.faceZones().setSize(zoneFaces.size());
653 procMesh.faceZones().set
658 zoneFaces[zoneI].shrink(),
659 zoneFaceFlips[zoneI].shrink(),
680 List<DynamicList<label>> zoneCells(cz.size());
685 zoneCells[zoneI].setCapacity(cz[zoneI].size() / nProcs_);
688 forAll(curCellLabels, celli)
690 label curCelli = curCellLabels[celli];
692 label zoneI = cellToZone[curCelli];
697 zoneCells[zoneI].append(celli);
699 else if (zoneI == -2)
704 label index = cz[zoneI].whichCell(curCelli);
708 zoneCells[zoneI].append(celli);
714 procMesh.cellZones().clearAddressing();
715 procMesh.cellZones().setSize(zoneCells.size());
718 procMesh.cellZones().set
723 zoneCells[zoneI].shrink(),
743 if (facesInstancePointsPtr_.valid())
759 pointsInstancePoints.write();
768 const cellSet& cs = cellSets[i];
769 cellSet
set(procMesh, cs.name(), cs.size()/nProcs_);
772 if (cs.found(curCellLabels[i]))
781 const faceSet& cs = faceSets[i];
782 faceSet
set(procMesh, cs.name(), cs.size()/nProcs_);
785 if (cs.found(
mag(curFaceLabels[i])-1))
794 const pointSet& cs = pointSets[i];
795 pointSet
set(procMesh, cs.name(), cs.size()/nProcs_);
798 if (cs.found(curPointLabels[i]))
822 procCellAddressing_[proci],
823 procPointAddressing_[proci]
830 <<
"Processor " << proci <<
nl 831 <<
" Number of cells = " << procMesh.nCells()
834 maxProcCells =
max(maxProcCells, procMesh.nCells());
836 label nBoundaryFaces = 0;
837 label nProcPatches = 0;
838 label nProcFaces = 0;
842 if (isA<processorPolyPatch>(procMesh.boundaryMesh()[
patchi]))
844 const processorPolyPatch& ppp =
845 refCast<const processorPolyPatch>
847 procMesh.boundaryMesh()[
patchi]
850 Info<<
" Number of faces shared with processor " 851 << ppp.neighbProcNo() <<
" = " << ppp.size() <<
endl;
854 nProcFaces += ppp.size();
858 nBoundaryFaces += procMesh.boundaryMesh()[
patchi].size();
862 Info<<
" Number of processor patches = " << nProcPatches <<
nl 863 <<
" Number of processor faces = " << nProcFaces <<
nl 864 <<
" Number of boundary faces = " << nBoundaryFaces <<
endl;
866 totProcFaces += nProcFaces;
867 totProcPatches += nProcPatches;
868 maxProcPatches =
max(maxProcPatches, nProcPatches);
869 maxProcFaces =
max(maxProcFaces, nProcFaces);
876 "pointProcAddressing",
877 procMesh.facesInstance(),
883 procPointAddressing_[proci]
885 pointProcAddressing.write();
891 "faceProcAddressing",
892 procMesh.facesInstance(),
898 procFaceAddressing_[proci]
906 "cellProcAddressing",
907 procMesh.facesInstance(),
913 procCellAddressing_[proci]
915 cellProcAddressing.write();
919 label nMeshPatches = curPatchSizes.size();
921 procBoundaryAddressing.setSize(nMeshPatches+nProcPatches, -1);
927 "boundaryProcAddressing",
928 procMesh.facesInstance(),
934 procBoundaryAddressing
936 boundaryProcAddressing.write();
939 scalar avgProcCells = scalar(nCells())/nProcs_;
940 scalar avgProcPatches = scalar(totProcPatches)/nProcs_;
941 scalar avgProcFaces = scalar(totProcFaces)/nProcs_;
944 if (totProcPatches == 0)
948 if (totProcFaces == 0)
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
PtrList< polyPatch > polyPatchList
container classes for polyPatch
List< labelList > labelListList
A List of labelList.
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.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
const word & name() const
Return name.
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.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
vectorIOField pointIOField
pointIOField is a vectorIOField.
PtrList< labelIOList > & faceProcAddressing
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
static unsigned int defaultPrecision()
Return the default precision.
IOobject(const word &name, const fileName &instance, const objectRegistry ®istry, readOption r=NO_READ, writeOption w=NO_WRITE, bool registerObject=true)
Construct from name, instance, registry, io options.
Ostream & endl(Ostream &os)
Add newline and flush stream.
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
~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.
const fileName & pointsInstance() const
Return the current instance directory for points.
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.
ZoneMesh< pointZone, polyMesh > pointZoneMesh
A ZoneMesh with the type pointZone.
static word controlDictName
The default control dictionary name (normally "controlDict")
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
static autoPtr< dictionary > New(Istream &)
Construct top-level dictionary on freestore from Istream.
word name(const complex &)
Return a string representation of a complex.
fvMesh(const IOobject &io)
Construct from IOobject.
void setSize(const label)
Reset size of List.
dimensioned< scalar > mag(const dimensioned< Type > &)
List< cell > cellList
list of cells
ZoneMesh< cellZone, polyMesh > cellZoneMesh
A ZoneMesh with the type cellZone.
IOList< label > labelIOList
Label container classes.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
polyMesh(const IOobject &io)
Construct from IOobject.