32 void Foam::pairPatchAgglomeration::compactLevels(
const label nCreatedLevels)
40 bool Foam::pairPatchAgglomeration::continueAgglomerating
42 const label nCoarseFaces
46 label localnCoarseFaces = nCoarseFaces;
47 bool contAgg = localnCoarseFaces >= nFacesInCoarsestLevel_;
52 void Foam::pairPatchAgglomeration::setBasedEdgeWeights()
54 const bPatch& coarsePatch = patchLevels_[0];
55 forAll(coarsePatch.edges(), i)
57 if (coarsePatch.isInternalEdge(i))
60 coarsePatch.edges()[i].mag(coarsePatch.localPoints());
62 const labelList& eFaces = coarsePatch.edgeFaces()[i];
64 if (eFaces.size() == 2)
67 coarsePatch.faceNormals()[eFaces[0]]
68 & coarsePatch.faceNormals()[eFaces[1]];
70 const edge edgeCommon = edge(eFaces[0], eFaces[1]);
72 if (facePairWeight_.found(edgeCommon))
74 facePairWeight_[edgeCommon] += edgeLength;
78 facePairWeight_.insert(edgeCommon, edgeLength);
83 facePairWeight_[edgeCommon] = -1.0;
90 for (
label k = j+1;
k<eFaces.size();
k++)
92 facePairWeight_.insert
94 edge(eFaces[j], eFaces[
k]),
105 void Foam::pairPatchAgglomeration::setEdgeWeights
107 const label fineLevelIndex
111 const bPatch& coarsePatch = patchLevels_[fineLevelIndex];
113 const labelList& fineToCoarse = restrictAddressing_[fineLevelIndex];
114 const label nCoarseI =
max(fineToCoarse) + 1;
117 HashSet<edge, Hash<edge>> fineFeaturedFaces(coarsePatch.nEdges()/10);
124 const edge
e = iter.key();
125 const edge edgeFeatured
130 fineFeaturedFaces.insert(edgeFeatured);
135 facePairWeight_.clear();
136 facePairWeight_.resize(coarsePatch.nEdges());
138 forAll(coarsePatch.edges(), i)
140 if (coarsePatch.isInternalEdge(i))
143 coarsePatch.edges()[i].mag(coarsePatch.localPoints());
145 const labelList& eFaces = coarsePatch.edgeFaces()[i];
147 if (eFaces.size() == 2)
149 const edge edgeCommon = edge(eFaces[0], eFaces[1]);
151 if (facePairWeight_.found(edgeCommon))
153 facePairWeight_[edgeCommon] += edgeLength;
157 facePairWeight_.insert(edgeCommon, edgeLength);
162 if (fineFeaturedFaces.found(edgeCommon))
164 facePairWeight_[edgeCommon] = -1.0;
172 for (
label k = j+1;
k<eFaces.size();
k++)
174 facePairWeight_.insert
176 edge(eFaces[j], eFaces[
k]),
193 const bool additionalWeights
201 nFacesInCoarsestLevel_
210 restrictAddressing_(maxLevels_),
211 restrictTopBottomAddressing_(
identityMap(patch.size())),
212 patchLevels_(maxLevels_),
213 facePairWeight_(patch.size())
230 setBasedEdgeWeights();
247 return patchLevels_[i];
251 void Foam::pairPatchAgglomeration::mapBaseToTopAgglom
253 const label fineLevelIndex
256 const labelList& fineToCoarse = restrictAddressing_[fineLevelIndex];
257 forAll(restrictTopBottomAddressing_, i)
259 restrictTopBottomAddressing_[i] =
260 fineToCoarse[restrictTopBottomAddressing_[i]];
265 bool Foam::pairPatchAgglomeration::agglomeratePatch
269 const label fineLevelIndex
272 if (
min(fineToCoarse) == -1)
278 if (fineToCoarse.size() == 0)
283 if (fineToCoarse.size() != patch.size())
286 <<
"restrict map does not correspond to fine level. " <<
endl
287 <<
" Sizes: restrictMap: " << fineToCoarse.size()
288 <<
" nEqns: " << patch.size()
292 const label nCoarseI =
max(fineToCoarse)+1;
293 List<face> patchFaces(nCoarseI);
299 for (
label coarseI = 0; coarseI < nCoarseI; coarseI++)
301 const labelList& fineFaces = coarseToFine[coarseI];
306 IndirectList<face>(patch, fineFaces),
310 if (upp.edgeLoops().size() != 1)
312 if (fineFaces.size() == 2)
314 const edge
e(fineFaces[0], fineFaces[1]);
315 facePairWeight_[
e] = -1.0;
317 else if (fineFaces.size() == 3)
319 const edge
e(fineFaces[0], fineFaces[1]);
320 const edge e1(fineFaces[0], fineFaces[2]);
321 const edge e2(fineFaces[2], fineFaces[1]);
322 facePairWeight_[
e] = -1.0;
323 facePairWeight_[e1] = -1.0;
324 facePairWeight_[e2] = -1.0;
330 patchFaces[coarseI] = face
345 SubList<face>(patchFaces, nCoarseI, 0),
356 label nPairLevels = 0;
357 label nCreatedLevels = 1;
358 label nCoarseFaces = 0;
359 label nCoarseFacesOld = 0;
361 while (nCreatedLevels < maxLevels_)
363 const bPatch& patch = patchLevels_[nCreatedLevels - 1];
365 bool agglomOK =
false;
369 label nCoarseFacesPrev = nCoarseFaces;
371 finalAgglomPtr = agglomerateOneLevel
377 if (nCoarseFaces > 0 && nCoarseFaces != nCoarseFacesPrev)
382 agglomOK = agglomeratePatch
391 restrictAddressing_.set(nCreatedLevels, finalAgglomPtr);
392 mapBaseToTopAgglom(nCreatedLevels);
393 setEdgeWeights(nCreatedLevels);
395 if (nPairLevels % mergeLevels_)
397 combineLevels(nCreatedLevels);
417 nFaces_[nCreatedLevels] = nCoarseFaces;
421 !continueAgglomerating(nCoarseFaces)
422 || (nCoarseFacesOld == nCoarseFaces)
428 nCoarseFacesOld = nCoarseFaces;
431 compactLevels(nCreatedLevels);
441 const label nFineFaces = patch.size();
444 labelField& coarseCellMap = tcoarseCellMap.ref();
452 const labelList& fFaces = faceFaces[facei];
454 if (coarseCellMap[facei] < 0)
456 label matchFaceNo = -1;
457 label matchFaceNeibNo = -1;
458 scalar maxFaceWeight = -great;
463 label faceNeig = fFaces[i];
464 const edge edgeCommon =
edge(facei, faceNeig);
467 facePairWeight_[edgeCommon] > maxFaceWeight
468 && coarseCellMap[faceNeig] < 0
469 && facePairWeight_[edgeCommon] != -1.0
474 matchFaceNeibNo = faceNeig;
475 maxFaceWeight = facePairWeight_[edgeCommon];
479 if (matchFaceNo >= 0)
482 coarseCellMap[matchFaceNo] = nCoarseFaces;
483 coarseCellMap[matchFaceNeibNo] = nCoarseFaces;
490 label clusterMatchFaceNo = -1;
491 scalar clusterMaxFaceCoeff = -great;
495 label faceNeig = fFaces[i];
496 const edge edgeCommon = edge(facei, faceNeig);
499 facePairWeight_[edgeCommon] > clusterMaxFaceCoeff
500 && facePairWeight_[edgeCommon] != -1.0
501 && coarseCellMap[faceNeig] > 0
504 clusterMatchFaceNo = faceNeig;
505 clusterMaxFaceCoeff = facePairWeight_[edgeCommon];
509 if (clusterMatchFaceNo >= 0)
512 coarseCellMap[facei] = coarseCellMap[clusterMatchFaceNo];
517 coarseCellMap[facei] = nCoarseFaces;
526 for (
label facei=0; facei<nFineFaces; facei++)
528 if (coarseCellMap[facei] < 0)
532 <<
" is not part of a cluster"
537 return tcoarseCellMap;
541 void Foam::pairPatchAgglomeration::combineLevels(
const label curLevel)
543 label prevLevel = curLevel - 1;
546 nFaces_[prevLevel] = nFaces_[curLevel];
551 const labelList& curResAddr = restrictAddressing_[curLevel];
552 labelList& prevResAddr = restrictAddressing_[prevLevel];
556 prevResAddr[i] = curResAddr[prevResAddr[i]];
560 restrictAddressing_.set(curLevel,
nullptr);
562 patchLevels_.set(prevLevel, patchLevels_.set(curLevel,
nullptr));
#define forAll(list, i)
Loop across all elements in list.
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
void setSize(const label)
Reset size of List.
A list of faces which address into the list of points.
const List< FaceType > & localFaces() const
Return patch faces addressing into local point list.
const labelListList & faceFaces() const
Return face-face addressing.
const Field< PointType > & localPoints() const
Return pointField of points in patch.
label size() const
Return the number of elements in the UList.
A list of keyword definitions, which are a keyword followed by any number of values (e....
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
const bPatch & patchLevel(const label leveli) const
Return primitivePatch of given level.
~pairPatchAgglomeration()
PtrList< bPatch > patchLevels_
Hierarchy of patch addressing.
labelList nFaces_
The number of faces in each level.
void agglomerate()
Agglomerate patch.
pairPatchAgglomeration(const polyPatch &patch, const dictionary &controlDict, const bool additionalWeights=false)
Construct given mesh and controls.
PtrList< labelField > restrictAddressing_
Cell restriction addressing array.
A patch is a list of labels that address the faces in the global face list.
A class for managing temporary objects.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
errorManipArg< error, int > exit(error &err, const int errNo=1)
List< label > labelList
A List of labels.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Ostream & endl(Ostream &os)
Add newline and flush stream.
errorManip< error > abort(error &err)
ListType renumber(const labelUList &oldToNew, const ListType &)
Renumber the values (not the indices) of a list.
labelListList invertOneToMany(const label len, const labelUList &)
Invert one-to-many map. Unmapped elements will be size 0.
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
List< labelList > labelListList
A List of labelList.
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
Field< label > labelField
Specialisation of Field<T> for label.
labelList identityMap(const label len)
Create identity map (map[i] == i) of given length.
PrimitivePatch< IndirectList< face >, const pointField & > indirectPrimitivePatch
Foam::indirectPrimitivePatch.
PrimitivePatch< faceList, const pointField > bPatch
dimensionedScalar cos(const dimensionedScalar &ds)
const unitConversion unitDegrees
runTime controlDict().lookup("adjustTimeStep") >> adjustTimeStep