51 if (
findIndex(elems2, elems1[elemI]) != -1)
61 bool Foam::meshCutter::isIn
67 label index =
findIndex(cuts, twoCuts[0]);
76 cuts[cuts.fcIndex(index)] == twoCuts[1]
77 || cuts[cuts.rcIndex(index)] == twoCuts[1]
85 Foam::label Foam::meshCutter::findCutCell
93 label cellI = cellLabels[
labelI];
95 if (cuts.cellLoops()[cellI].size())
108 Foam::label Foam::meshCutter::findInternalFacePoint
115 label pointI = pointLabels[
labelI];
121 label faceI = pFaces[pFaceI];
123 if (
mesh().isInternalFace(faceI))
130 if (pointLabels.empty())
132 FatalErrorIn(
"meshCutter::findInternalFacePoint(const labelList&)")
142 void Foam::meshCutter::faceCells
144 const cellCuts& cuts,
157 if (cellLoops[own].size() && uses(f, anchorPts[own]))
159 own = addedCells_[own];
164 if (
mesh().isInternalFace(faceI))
168 if (cellLoops[nei].size() && uses(f, anchorPts[nei]))
170 nei = addedCells_[nei];
176 void Foam::meshCutter::getFaceInfo
186 if (!
mesh().isInternalFace(faceI))
199 zoneFlip = fZone.flipMap()[fZone.whichFace(faceI)];
205 void Foam::meshCutter::addFace
207 polyTopoChange& meshMod,
214 label patchID, zoneID, zoneFlip;
216 getFaceInfo(faceI, patchID, zoneID, zoneFlip);
218 if ((nei == -1) || (own < nei))
223 Pout<<
"Adding face " << newFace
224 <<
" with new owner:" << own
225 <<
" with new neighbour:" << nei
226 <<
" patchID:" << patchID
227 <<
" zoneID:" << zoneID
228 <<
" zoneFlip:" << zoneFlip
254 Pout<<
"Adding (reversed) face " << newFace.reverseFace()
255 <<
" with new owner:" << nei
256 <<
" with new neighbour:" << own
257 <<
" patchID:" << patchID
258 <<
" zoneID:" << zoneID
259 <<
" zoneFlip:" << zoneFlip
267 newFace.reverseFace(),
284 void Foam::meshCutter::modFace
286 polyTopoChange& meshMod,
293 label patchID, zoneID, zoneFlip;
295 getFaceInfo(faceI, patchID, zoneID, zoneFlip);
299 (own !=
mesh().faceOwner()[faceI])
301 mesh().isInternalFace(faceI)
302 && (nei !=
mesh().faceNeighbour()[faceI])
304 || (newFace !=
mesh().faces()[faceI])
309 Pout<<
"Modifying face " << faceI
310 <<
" old vertices:" <<
mesh().
faces()[faceI]
311 <<
" new vertices:" << newFace
312 <<
" new owner:" << own
313 <<
" new neighbour:" << nei
314 <<
" new zoneID:" << zoneID
315 <<
" new zoneFlip:" << zoneFlip
319 if ((nei == -1) || (own < nei))
343 newFace.reverseFace(),
360 void Foam::meshCutter::copyFace
374 newFace[newFp++] = f[fp];
376 fp = (fp + 1) % f.size();
378 newFace[newFp] = f[fp];
386 void Foam::meshCutter::splitFace
403 "meshCutter::splitFace"
404 ", const face&, const label, const label, face&, face&)"
405 ) <<
"Cannot find vertex (new numbering) " << v0
416 "meshCutter::splitFace("
417 ", const face&, const label, const label, face&, face&)"
418 ) <<
"Cannot find vertex (new numbering) " << v1
424 f0.setSize((endFp + 1 + f.size() - startFp) % f.size());
425 f1.setSize(f.size() - f0.size() + 2);
427 copyFace(f, startFp, endFp, f0);
428 copyFace(f, endFp, startFp, f1);
434 Foam::face Foam::meshCutter::addEdgeCutsToFace(
const label faceI)
const
438 face newFace(2 * f.size());
445 newFace[newFp++] = f[fp];
448 label fp1 = f.fcIndex(fp);
450 HashTable<label, edge, Hash<edge> >::const_iterator fnd =
451 addedPoints_.find(edge(f[fp], f[fp1]));
453 if (fnd != addedPoints_.end())
456 newFace[newFp++] = fnd();
460 newFace.setSize(newFp);
475 face newFace(2*loop.size());
481 label cut = loop[fp];
485 label edgeI = getEdge(cut);
489 label vertI = addedPoints_[
e];
491 newFace[newFaceI++] = vertI;
496 label vertI = getVertex(cut);
498 newFace[newFaceI++] = vertI;
500 label nextCut = loop[loop.fcIndex(fp)];
502 if (!isEdge(nextCut))
505 label nextVertI = getVertex(nextCut);
512 HashTable<label, edge, Hash<edge> >::const_iterator fnd =
513 addedPoints_.find(
mesh().edges()[edgeI]);
515 if (fnd != addedPoints_.end())
517 newFace[newFaceI++] = fnd();
523 newFace.setSize(newFaceI);
558 addedCells_.resize(cuts.
nLoops());
561 addedFaces_.resize(cuts.
nLoops());
563 addedPoints_.clear();
564 addedPoints_.resize(cuts.
nLoops());
585 if (debug && findCutCell(cuts,
mesh().edgeCells()[edgeI]) == -1)
589 "meshCutter::setRefinement(const cellCuts&"
591 ) <<
"Problem: cut edge but none of the cells using it is\n"
592 <<
"edge:" << edgeI <<
" verts:" << e
597 label masterPointI = e.
start();
604 point newPt = weight*v1 + (1.0-weight)*v0;
619 addedPoints_.insert(e, addedPointI);
623 Pout<<
"Added point " << addedPointI
625 << masterPointI <<
" of edge " << edgeI
626 <<
" vertices " << e <<
endl;
637 if (cellLoops[cellI].size())
653 addedCells_.insert(cellI, addedCellI);
657 Pout<<
"Added cell " << addedCells_[cellI] <<
" to cell "
670 const labelList& loop = cellLoops[cellI];
678 face newFace(loopToFace(cellI, loop));
681 label masterPointI = findInternalFacePoint(anchorPts[cellI]);
701 addedFaces_.insert(cellI, addedFaceI);
719 Pout<<
"Added splitting face " << newFace <<
" index:"
721 <<
" to owner " << cellI
722 <<
" neighbour " << addedCells_[cellI]
724 writeCuts(
Pout, loop, weights);
745 iter != faceSplitCuts.
end();
749 label faceI = iter.key();
752 face newFace(addEdgeCutsToFace(faceI));
755 const edge& splitEdge = iter();
757 label cut0 = splitEdge[0];
762 label edgeI = getEdge(cut0);
763 v0 = addedPoints_[
mesh().
edges()[edgeI]];
767 v0 = getVertex(cut0);
770 label cut1 = splitEdge[1];
774 label edgeI = getEdge(cut1);
775 v1 = addedPoints_[
mesh().
edges()[edgeI]];
779 v1 = getVertex(cut1);
784 splitFace(newFace, v0, v1, f0, f1);
790 if (
mesh().isInternalFace(faceI))
798 <<
" own:" << own <<
" nei:" << nei
800 <<
" and f1:" << f1 <<
endl;
816 if (cellLoops[own].empty())
821 else if (isIn(splitEdge, cellLoops[own]))
825 if (uses(f0, anchorPts[own]))
827 f0Owner = addedCells_[own];
833 f1Owner = addedCells_[own];
840 if (uses(f, anchorPts[own]))
842 label newCellI = addedCells_[own];
854 label f0Neighbour = -1;
855 label f1Neighbour = -1;
859 if (cellLoops[nei].empty())
864 else if (isIn(splitEdge, cellLoops[nei]))
868 if (uses(f0, anchorPts[nei]))
870 f0Neighbour = addedCells_[nei];
876 f1Neighbour = addedCells_[nei];
883 if (uses(f, anchorPts[nei]))
885 label newCellI = addedCells_[nei];
886 f0Neighbour = newCellI;
887 f1Neighbour = newCellI;
898 addFace(meshMod, faceI, f0, f0Owner, f0Neighbour);
900 modFace(meshMod, faceI, f1, f1Owner, f1Neighbour);
902 faceUptodate[faceI] =
true;
915 if (edgeIsCut[edgeI])
921 label faceI = eFaces[i];
923 if (!faceUptodate[faceI])
926 face newFace(addEdgeCutsToFace(faceI));
930 Pout<<
"Added edge cuts to face " << faceI
932 <<
" newFace:" << newFace <<
endl;
937 faceCells(cuts, faceI, own, nei);
939 modFace(meshMod, faceI, newFace, own, nei);
941 faceUptodate[faceI] =
true;
954 if (cellLoops[cellI].size())
958 forAll(cllFaces, cllFaceI)
960 label faceI = cllFaces[cllFaceI];
962 if (!faceUptodate[faceI])
967 if (debug && (f != addEdgeCutsToFace(faceI)))
971 "meshCutter::setRefinement(const cellCuts&"
973 ) <<
"Problem: edges added to face which does "
974 <<
" not use a marked cut" <<
endl
975 <<
"faceI:" << faceI <<
endl
976 <<
"face:" << f <<
endl
977 <<
"newFace:" << addEdgeCutsToFace(faceI)
983 faceCells(cuts, faceI, own, nei);
994 faceUptodate[faceI] =
true;
1002 Pout<<
"meshCutter:" <<
nl
1003 <<
" cells split:" << addedCells_.size() <<
nl
1004 <<
" faces added:" << addedFaces_.size() <<
nl
1005 <<
" points added on edges:" << addedPoints_.size() <<
nl
1018 Map<label> newAddedCells(addedCells_.size());
1023 iter != addedCells_.end();
1027 label cellI = iter.key();
1031 label addedCellI = iter();
1035 if (newCellI >= 0 && newAddedCellI >= 0)
1040 && (newCellI != cellI || newAddedCellI != addedCellI)
1043 Pout<<
"meshCutter::updateMesh :"
1044 <<
" updating addedCell for cell " << cellI
1045 <<
" from " << addedCellI
1046 <<
" to " << newAddedCellI <<
endl;
1048 newAddedCells.insert(newCellI, newAddedCellI);
1053 addedCells_.transfer(newAddedCells);
1057 Map<label> newAddedFaces(addedFaces_.size());
1062 iter != addedFaces_.end();
1066 label cellI = iter.key();
1070 label addedFaceI = iter();
1074 if ((newCellI >= 0) && (newAddedFaceI >= 0))
1079 && (newCellI != cellI || newAddedFaceI != addedFaceI)
1082 Pout<<
"meshCutter::updateMesh :"
1083 <<
" updating addedFace for cell " << cellI
1084 <<
" from " << addedFaceI
1085 <<
" to " << newAddedFaceI
1088 newAddedFaces.insert(newCellI, newAddedFaceI);
1093 addedFaces_.transfer(newAddedFaces);
1102 addedPoints_.begin();
1103 iter != addedPoints_.end();
1107 const edge& e = iter.key();
1113 label addedPointI = iter();
1117 if ((newStart >= 0) && (newEnd >= 0) && (newAddedPointI >= 0))
1119 edge newE =
edge(newStart, newEnd);
1124 && (e != newE || newAddedPointI != addedPointI)
1127 Pout<<
"meshCutter::updateMesh :"
1128 <<
" updating addedPoints for edge " << e
1129 <<
" from " << addedPointI
1130 <<
" to " << newAddedPointI
1134 newAddedPoints.insert(newE, newAddedPointI);
1139 addedPoints_.transfer(newAddedPoints);