43 const scalar faceCoupleInfo::angleTol_ = 1
e-3;
49void Foam::faceCoupleInfo::writeOBJ
51 const fileName& fName,
52 const edgeList& edges,
69 const edge&
e = edges[edgeI];
73 const label pointi =
e[eI];
75 if (pointMap[pointi] == -1)
96 const edge&
e = edges[edgeI];
98 str<<
"l " << pointMap[
e[0]]+1 <<
' ' << pointMap[
e[1]]+1 <<
nl;
103void Foam::faceCoupleInfo::writeOBJ
105 const fileName& fName,
110 Pout<<
"Writing connections as edges to " << fName <<
endl;
122 str <<
"l " << vertI-1 <<
' ' << vertI <<
nl;
127void Foam::faceCoupleInfo::writePointsFaces()
const
135 OFstream str(
"masterPatch.obj");
140 OFstream str(
"slavePatch.obj");
145 OFstream str(
"cutFaces.obj");
152 Pout<<
"Writing cutToMasterPoints to cutToMasterPoints.obj" <<
endl;
156 "cutToMasterPoints.obj",
161 Pout<<
"Writing cutToSlavePoints to cutToSlavePoints.obj" <<
endl;
165 "cutToSlavePoints.obj",
173 Pout<<
"Writing cutToMasterFaces to cutToMasterFaces.obj" <<
endl;
177 forAll(cutToMasterFaces(), cutFacei)
179 label masterFacei = cutToMasterFaces()[cutFacei];
181 if (masterFacei != -1)
183 equivMasterFaces[cutFacei] = m[masterFacei].centre(m.points());
188 <<
"No master face for cut face " << cutFacei
189 <<
" at position " <<
c[cutFacei].centre(
c.points())
192 equivMasterFaces[cutFacei] =
Zero;
198 "cutToMasterFaces.obj",
199 calcFaceCentres<List>(c, cutPoints(), 0,
c.size()),
205 Pout<<
"Writing cutToSlaveFaces to cutToSlaveFaces.obj" <<
endl;
209 forAll(cutToSlaveFaces(), cutFacei)
211 label slaveFacei = cutToSlaveFaces()[cutFacei];
213 equivSlaveFaces[cutFacei] =
s[slaveFacei].centre(
s.points());
218 "cutToSlaveFaces.obj",
219 calcFaceCentres<List>(c, cutPoints(), 0,
c.size()),
228void Foam::faceCoupleInfo::writeEdges
240 OFstream str(
"cutToMasterEdges.obj");
241 Pout<<
"Writing cutToMasterEdges to " << str.
name() <<
endl;
245 forAll(cutToMasterEdges, cutEdgeI)
247 if (cutToMasterEdges[cutEdgeI] != -1)
249 const edge& masterEdge = m.edges()[cutToMasterEdges[cutEdgeI]];
250 const edge& cutEdge =
c.edges()[cutEdgeI];
260 str <<
"l " << vertI-3 <<
' ' << vertI-2 <<
nl;
261 str <<
"l " << vertI-3 <<
' ' << vertI-1 <<
nl;
262 str <<
"l " << vertI-3 <<
' ' << vertI <<
nl;
263 str <<
"l " << vertI-2 <<
' ' << vertI-1 <<
nl;
264 str <<
"l " << vertI-2 <<
' ' << vertI <<
nl;
265 str <<
"l " << vertI-1 <<
' ' << vertI <<
nl;
270 OFstream str(
"cutToSlaveEdges.obj");
271 Pout<<
"Writing cutToSlaveEdges to " << str.
name() <<
endl;
279 if (slaveToCut[edgeI] != -1)
281 const edge& slaveEdge =
s.edges()[edgeI];
282 const edge& cutEdge =
c.edges()[slaveToCut[edgeI]];
292 str <<
"l " << vertI-3 <<
' ' << vertI-2 <<
nl;
293 str <<
"l " << vertI-3 <<
' ' << vertI-1 <<
nl;
294 str <<
"l " << vertI-3 <<
' ' << vertI <<
nl;
295 str <<
"l " << vertI-2 <<
' ' << vertI-1 <<
nl;
296 str <<
"l " << vertI-2 <<
' ' << vertI <<
nl;
297 str <<
"l " << vertI-1 <<
' ' << vertI <<
nl;
315 forAll(toPatchEdges, edgeI)
317 const edge&
e = edges[edgeI];
319 label v0 = pointMap[
e[0]];
320 label v1 = pointMap[
e[1]];
322 toPatchEdges[edgeI] =
326 patch.pointEdges()[v0],
335bool Foam::faceCoupleInfo::regionEdge
337 const polyMesh& slaveMesh,
338 const label slaveEdgeI
341 const labelList& eFaces = slavePatch().edgeFaces()[slaveEdgeI];
343 if (eFaces.size() == 1)
355 label facei = eFaces[i];
357 label meshFacei = slavePatch().addressing()[facei];
359 label patchi = slaveMesh.boundaryMesh().whichPatch(meshFacei);
365 else if (patchi != patch0)
377Foam::label Foam::faceCoupleInfo::mostAlignedCutEdge
380 const polyMesh& slaveMesh,
381 const bool patchDivision,
385 const label edgeStart,
393 const pointField& localPoints = cutFaces().localPoints();
395 const labelList& pEdges = cutFaces().pointEdges()[pointi];
399 Pout<<
"mostAlignedEdge : finding nearest edge among "
400 << UIndirectList<edge>(cutFaces().edges(), pEdges)
401 <<
" connected to point " << pointi
402 <<
" coord:" << localPoints[pointi]
403 <<
" running between " << edgeStart <<
" coord:"
404 << localPoints[edgeStart]
405 <<
" and " << edgeEnd <<
" coord:"
406 << localPoints[edgeEnd]
413 scalar maxCos = -GREAT;
417 label edgeI = pEdges[i];
423 && cutToMasterEdges[edgeI] == -1
427 && regionEdge(slaveMesh, cutToSlaveEdges[edgeI])
431 const edge&
e = cutFaces().edges()[edgeI];
433 label otherPointi =
e.otherVertex(pointi);
435 if (otherPointi == edgeEnd)
440 Pout<<
" mostAlignedEdge : found end point " << edgeEnd
448 vector eVec(localPoints[otherPointi] - localPoints[pointi]);
450 scalar magEVec =
mag(eVec);
452 if (magEVec < VSMALL)
455 <<
"Crossing zero sized edge " << edgeI
456 <<
" coords:" << localPoints[otherPointi]
457 << localPoints[pointi]
458 <<
" when walking from " << localPoints[edgeStart]
459 <<
" to " << localPoints[edgeEnd]
466 const vector eToEndPoint =
469 localPoints[edgeEnd] - localPoints[otherPointi]
472 scalar cosAngle = eVec & eToEndPoint;
476 Pout<<
" edge:" <<
e <<
" points:" << localPoints[pointi]
477 << localPoints[otherPointi]
479 <<
" vecToEnd:" << eToEndPoint
480 <<
" cosAngle:" << cosAngle
484 if (cosAngle > maxCos)
492 if (maxCos > 1 - angleTol_)
503void Foam::faceCoupleInfo::setCutEdgeToPoints(
const labelList& cutToMasterEdges)
514 const edgeList& cutEdges = cutFaces().edges();
524 forAll(masterToCutEdges, masterEdgeI)
526 const edge& masterE = masterPatch().edges()[masterEdgeI];
531 const labelList& stringedEdges = masterToCutEdges[masterEdgeI];
533 if (stringedEdges.empty())
536 <<
"Did not match all of master edges to cutFace edges"
538 <<
"First unmatched edge:" << masterEdgeI <<
" endPoints:"
539 << masterPatch().localPoints()[masterE[0]]
540 << masterPatch().localPoints()[masterE[1]]
542 <<
"This usually means that the slave patch is not a"
543 <<
" subdivision of the master patch"
546 else if (stringedEdges.size() > 1)
551 DynamicList<label> splitPoints(stringedEdges.size()-1);
554 const edge unsplitEdge
556 masterToCutPoints_[masterE[0]],
557 masterToCutPoints_[masterE[1]]
560 label startVertI = unsplitEdge[0];
561 label startEdgeI = -1;
563 while (startVertI != unsplitEdge[1])
571 label oldStart = startVertI;
575 label edgeI = stringedEdges[i];
577 if (edgeI != startEdgeI)
579 const edge&
e = cutEdges[edgeI];
585 if (
e[0] == startVertI)
589 if (
e[1] != unsplitEdge[1])
595 else if (
e[1] == startVertI)
599 if (
e[0] != unsplitEdge[1])
609 if (oldStart == startVertI)
612 <<
" unsplitEdge:" << unsplitEdge
613 <<
" does not correspond to split edges "
614 << UIndirectList<edge>(cutEdges, stringedEdges)
629 cutEdgeToPoints_.insert(unsplitEdge, splitPoints.shrink());
635Foam::label Foam::faceCoupleInfo::matchFaces
642 const bool sameOrientation
645 if (f0.size() != f1.size())
648 <<
"Different sizes for supposedly matching faces." <<
nl
649 <<
"f0:" << f0 <<
" coords:" << UIndirectList<point>(
points0, f0)
651 <<
"f1:" << f1 <<
" coords:" << UIndirectList<point>(points1, f1)
655 const scalar absTolSqr =
sqr(absTol);
663 bool fullMatch =
true;
672 if (distSqr > absTolSqr)
678 fp0 = f0.fcIndex(fp0);
682 fp1 = f1.fcIndex(fp1);
686 fp1 = f1.rcIndex(fp1);
700 <<
"No unique match between two faces" <<
nl
701 <<
"Face " << f0 <<
" coords "
702 << UIndirectList<point>(
points0, f0) <<
nl
703 <<
"Face " << f1 <<
" coords "
704 << UIndirectList<point>(points1, f1)
705 <<
"when using tolerance " << absTol
706 <<
" and forwardMatching:" << sameOrientation
714bool Foam::faceCoupleInfo::matchPointsThroughFaces
721 const bool sameOrientation,
734 patchToCutPoints.setSize(patchPoints.size());
735 patchToCutPoints = -1;
739 labelList cutPointRegion(cutPoints.size(), -1);
740 DynamicList<label> cutPointRegionMaster;
742 forAll(patchFaces, patchFacei)
744 const face& patchF = patchFaces[patchFacei];
747 const face& cutF = cutFaces[patchFacei];
750 label patchFp = matchFaces
762 label cutPointi = cutF[cutFp];
763 label patchPointi = patchF[patchFp];
775 if (patchToCutPoints[patchPointi] == -1)
777 patchToCutPoints[patchPointi] = cutPointi;
779 else if (patchToCutPoints[patchPointi] != cutPointi)
783 label otherCutPointi = patchToCutPoints[patchPointi];
792 if (cutPointRegion[otherCutPointi] != -1)
795 label region = cutPointRegion[otherCutPointi];
796 cutPointRegion[cutPointi] = region;
799 cutPointRegionMaster[region] =
min
801 cutPointRegionMaster[region],
808 label region = cutPointRegionMaster.size();
809 cutPointRegionMaster.append
811 min(cutPointi, otherCutPointi)
813 cutPointRegion[cutPointi] = region;
814 cutPointRegion[otherCutPointi] = region;
820 patchFp = patchF.fcIndex(patchFp);
824 patchFp = patchF.rcIndex(patchFp);
830 compactToCut.setSize(cutPointRegion.size());
831 cutToCompact.setSize(cutPointRegion.size());
833 label compactPointi = 0;
837 if (cutPointRegion[i] == -1)
840 cutToCompact[i] = compactPointi;
841 compactToCut[compactPointi] = i;
848 label masterPointi = cutPointRegionMaster[cutPointRegion[i]];
850 if (cutToCompact[masterPointi] == -1)
852 cutToCompact[masterPointi] = compactPointi;
853 compactToCut[compactPointi] = masterPointi;
856 cutToCompact[i] = cutToCompact[masterPointi];
859 compactToCut.setSize(compactPointi);
861 return compactToCut.size() != cutToCompact.size();
873 scalar maxDist = -GREAT;
877 const point& cutPt = cutPoints[cutF[fp]];
879 pointHit pHit = masterF.nearestPoint(cutPt, masterPoints);
881 maxDist =
max(maxDist, pHit.distance());
887void Foam::faceCoupleInfo::findPerfectMatchingFaces
889 const primitiveMesh& mesh0,
890 const primitiveMesh& mesh1,
898 if (!mesh0.nFaces() || !mesh1.nFaces())
911 calcFaceCentres<List>
915 mesh0.nInternalFaces(),
916 mesh0.nBoundaryFaces()
922 calcFaceCentres<List>
926 mesh1.nInternalFaces(),
927 mesh1.nBoundaryFaces()
934 Pout<<
"Face matching tolerance : " << absTol <<
endl;
952 <<
"Matched ALL " << fc1.size()
953 <<
" boundary faces of mesh0 to boundary faces of mesh1." <<
endl
954 <<
"This is only valid if the mesh to add is fully"
955 <<
" enclosed by the mesh it is added to." <<
endl;
962 mesh0Faces.setSize(fc0.size());
963 mesh1Faces.setSize(fc1.size());
967 if (from1To0[i] != -1)
969 mesh1Faces[nMatched] = i + mesh1.nInternalFaces();
970 mesh0Faces[nMatched] = from1To0[i] + mesh0.nInternalFaces();
976 mesh0Faces.setSize(nMatched);
977 mesh1Faces.setSize(nMatched);
981void Foam::faceCoupleInfo::findSlavesCoveringMaster
983 const primitiveMesh& mesh0,
984 const primitiveMesh& mesh1,
994 identity(mesh0.nBoundaryFaces(), mesh0.nInternalFaces())
997 treeBoundBox overallBb(mesh0.points());
1001 indexedOctree<treeDataFace> tree
1009 overallBb.extend(
rndGen, 1
e-4),
1017 Pout<<
"findSlavesCoveringMaster :"
1018 <<
" constructed octree for mesh0 boundary faces" <<
endl;
1028 label mesh1Facei = mesh1.nInternalFaces();
1029 mesh1Facei < mesh1.nFaces();
1033 const face& f1 = mesh1.faces()[mesh1Facei];
1036 point fc(f1.centre(mesh1.points()));
1042 label mesh0Facei = tree.shapes().faceLabels()[nearInfo.index()];
1055 mesh0.faces()[mesh0Facei],
1061 mesh0Set.insert(mesh0Facei);
1062 mesh1Set.insert(mesh1Facei);
1069 Pout<<
"findSlavesCoveringMaster :"
1070 <<
" matched " << mesh1Set.size() <<
" mesh1 faces to "
1071 << mesh0Set.size() <<
" mesh0 faces" <<
endl;
1074 mesh0Faces = mesh0Set.toc();
1075 mesh1Faces = mesh1Set.toc();
1079Foam::label Foam::faceCoupleInfo::growCutFaces
1082 Map<labelList>& candidates
1087 Pout<<
"growCutFaces :"
1088 <<
" growing cut faces to masterPatch" <<
endl;
1091 label nTotChanged = 0;
1100 forAll(cutToMasterFaces_, cutFacei)
1102 const label masterFacei = cutToMasterFaces_[cutFacei];
1104 if (masterFacei != -1)
1110 const labelList& fEdges = cutFaceEdges[cutFacei];
1114 const label cutEdgeI = fEdges[i];
1116 if (cutToMasterEdges[cutEdgeI] == -1)
1124 const labelList& eFaces = cutEdgeFaces[cutEdgeI];
1128 const label facei = eFaces[j];
1130 if (cutToMasterFaces_[facei] == -1)
1132 cutToMasterFaces_[facei] = masterFacei;
1133 candidates.erase(facei);
1136 else if (cutToMasterFaces_[facei] != masterFacei)
1139 cutFaces().points();
1141 masterPatch().points();
1143 const edge&
e = cutFaces().edges()[cutEdgeI];
1145 label myMaster = cutToMasterFaces_[facei];
1146 const face& myF = masterPatch()[myMaster];
1148 const face& nbrF = masterPatch()[masterFacei];
1152 << cutFaces()[facei].points(cutPoints)
1155 <<
" but also connects to nbr face "
1156 << cutFaces()[cutFacei].points(cutPoints)
1157 <<
" with master " << masterFacei
1160 << myF.points(masterPoints)
1161 <<
" nbrMasterFace:"
1162 << nbrF.points(masterPoints) <<
nl
1163 <<
"Across cut edge " << cutEdgeI
1165 << cutFaces().localPoints()[
e[0]]
1166 << cutFaces().localPoints()[
e[1]]
1177 Pout<<
"growCutFaces : Grown an additional " << nChanged
1178 <<
" cut to master face" <<
" correspondence" <<
endl;
1181 nTotChanged += nChanged;
1193void Foam::faceCoupleInfo::checkMatch(
const labelList& cutToMasterEdges)
const
1195 const pointField& cutLocalPoints = cutFaces().localPoints();
1197 const pointField& masterLocalPoints = masterPatch().localPoints();
1198 const faceList& masterLocalFaces = masterPatch().localFaces();
1200 forAll(cutToMasterEdges, cutEdgeI)
1202 const edge&
e = cutFaces().edges()[cutEdgeI];
1204 if (cutToMasterEdges[cutEdgeI] == -1)
1207 const labelList& cutEFaces = cutFaces().edgeFaces()[cutEdgeI];
1209 label masterFacei = -1;
1213 label cutFacei = cutEFaces[i];
1215 if (cutToMasterFaces_[cutFacei] != -1)
1217 if (masterFacei == -1)
1219 masterFacei = cutToMasterFaces_[cutFacei];
1221 else if (masterFacei != cutToMasterFaces_[cutFacei])
1223 label myMaster = cutToMasterFaces_[cutFacei];
1224 const face& myF = masterLocalFaces[myMaster];
1226 const face& nbrF = masterLocalFaces[masterFacei];
1229 <<
"Internal CutEdge " <<
e
1231 << cutLocalPoints[
e[0]]
1232 << cutLocalPoints[
e[1]]
1233 <<
" connects to master " << myMaster
1234 <<
" and to master " << masterFacei <<
nl
1236 << myF.points(masterLocalPoints)
1237 <<
" nbrMasterFace:"
1238 << nbrF.points(masterLocalPoints)
1248Foam::label Foam::faceCoupleInfo::matchEdgeFaces
1251 Map<labelList>& candidates
1260 candidates.resize(cutFaces().size());
1264 forAll(cutToMasterEdges, cutEdgeI)
1266 label masterEdgeI = cutToMasterEdges[cutEdgeI];
1268 if (masterEdgeI != -1)
1273 const labelList& cutEFaces = cutFaces().edgeFaces()[cutEdgeI];
1275 masterPatch().edgeFaces()[masterEdgeI];
1279 label cutFacei = cutEFaces[i];
1281 if (cutToMasterFaces_[cutFacei] == -1)
1290 auto fnd = candidates.find(cutFacei);
1296 candidates.insert(cutFacei, masterEFaces);
1303 const labelList& masterFaces = fnd.val();
1305 DynamicList<label> newCandidates(masterFaces.size());
1309 if (masterFaces.found(masterEFaces[j]))
1311 newCandidates.append(masterEFaces[j]);
1315 if (newCandidates.size() == 1)
1318 cutToMasterFaces_[cutFacei] = newCandidates[0];
1319 candidates.erase(cutFacei);
1332 fnd() = newCandidates.shrink();
1343 Pout<<
"matchEdgeFaces : Found " << nChanged
1344 <<
" faces where there was"
1345 <<
" only one remaining choice for cut-master correspondence"
1353Foam::label Foam::faceCoupleInfo::geometricMatchEdgeFaces
1355 Map<labelList>& candidates
1358 const pointField& cutPoints = cutFaces().points();
1368 masterPatch().size(),
1375 label cutFacei = iter.key();
1376 const labelList& masterFaces = iter.val();
1378 const face& cutF = cutFaces()[cutFacei];
1380 if (cutToMasterFaces_[cutFacei] == -1)
1383 scalar minDist = GREAT;
1384 label minMasterFacei = -1;
1388 label masterFacei = masterFaces[i];
1390 if (masterToCutFaces[masterFaces[i]].empty())
1392 scalar dist = maxDistance
1396 masterPatch()[masterFacei],
1403 minMasterFacei = masterFacei;
1408 if (minMasterFacei != -1)
1410 cutToMasterFaces_[cutFacei] = minMasterFacei;
1411 masterToCutFaces[minMasterFacei] = cutFacei;
1418 forAll(cutToMasterFaces_, cutFacei)
1420 if (cutToMasterFaces_[cutFacei] != -1)
1422 candidates.erase(cutFacei);
1429 Pout<<
"geometricMatchEdgeFaces : Found " << nChanged
1430 <<
" faces where there was"
1431 <<
" only one remaining choice for cut-master correspondence"
1439void Foam::faceCoupleInfo::perfectPointMatch
1441 const scalar absTol,
1442 const bool slaveFacesOrdered
1450 Pout<<
"perfectPointMatch :"
1451 <<
" Matching master and slave to cut."
1452 <<
" Master and slave faces are identical" <<
nl;
1454 if (slaveFacesOrdered)
1456 Pout<<
"and master and slave faces are ordered"
1457 <<
" (on coupled patches)" <<
endl;
1461 Pout<<
"and master and slave faces are not ordered"
1462 <<
" (on coupled patches)" <<
endl;
1466 cutToMasterFaces_ =
identity(masterPatch().size());
1467 cutPoints_ = masterPatch().localPoints();
1472 masterPatch().localFaces(),
1476 masterToCutPoints_ =
identity(cutPoints_.size());
1480 bool matchedAllFaces =
false;
1482 if (slaveFacesOrdered)
1484 cutToSlaveFaces_ =
identity(cutFaces().size());
1485 matchedAllFaces = (cutFaces().size() == slavePatch().size());
1495 calcFaceCentres<List>
1502 calcFaceCentres<IndirectList>
1518 if (!matchedAllFaces)
1520 labelList cutToSlaveFacesTemp(cutToSlaveFaces_.size(), -1);
1524 calcFacePointAverages<List>
1531 calcFacePointAverages<IndirectList>
1543 cutToSlaveFaces_ =
max(cutToSlaveFaces_, cutToSlaveFacesTemp);
1545 matchedAllFaces =
min(cutToSlaveFaces_) != -1;
1550 if (!matchedAllFaces)
1553 <<
"Did not match all of the master faces to the slave faces"
1555 <<
"This usually means that the slave patch and master patch"
1556 <<
" do not align to within " << absTol <<
" metre."
1566 matchPointsThroughFaces
1569 cutFaces().localPoints(),
1570 reorder(cutToSlaveFaces_, cutFaces().localFaces()),
1571 slavePatch().localPoints(),
1572 slavePatch().localFaces(),
1582 cutPoints_ = UIndirectList<point>(cutPoints_, compactToCut)();
1584 const faceList& cutLocalFaces = cutFaces().localFaces();
1586 faceList compactFaces(cutLocalFaces.size());
1589 compactFaces[i] =
renumber(cutToCompact, cutLocalFaces[i]);
1605void Foam::faceCoupleInfo::subDivisionMatch
1607 const polyMesh& slaveMesh,
1608 const bool patchDivision,
1614 Pout<<
"subDivisionMatch :"
1615 <<
" Matching master and slave to cut."
1616 <<
" Slave can be subdivision of master but all master edges"
1617 <<
" have to be covered by slave edges." <<
endl;
1622 cutPoints_ = slavePatch().localPoints();
1624 faceList cutFaces(slavePatch().size());
1628 cutFaces[i] = slavePatch().localFaces()[i].reverseFace();
1634 cutToSlaveFaces_ =
identity(cutFaces().size());
1651 OFstream str(
"regionEdges.obj");
1653 Pout<<
"subDivisionMatch :"
1654 <<
" addressing for slave patch fully done."
1655 <<
" Dumping region edges to " << str.
name() <<
endl;
1659 forAll(slavePatch().edges(), slaveEdgeI)
1661 if (regionEdge(slaveMesh, slaveEdgeI))
1663 const edge&
e = slavePatch().edges()[slaveEdgeI];
1669 str<<
"l " << vertI-1 <<
' ' << vertI <<
nl;
1682 Pout<<
"subDivisionMatch :"
1683 <<
" matching master points to cut points" <<
endl;
1688 masterPatch().localPoints(),
1695 if (!matchedAllPoints)
1698 <<
"Did not match all of the master points to the slave points"
1700 <<
"This usually means that the slave patch is not a"
1701 <<
" subdivision of the master patch"
1713 Pout<<
"subDivisionMatch :"
1714 <<
" matching cut edges to master edges" <<
endl;
1717 const edgeList& masterEdges = masterPatch().edges();
1718 const pointField& masterPoints = masterPatch().localPoints();
1720 const edgeList& cutEdges = cutFaces().edges();
1724 forAll(masterEdges, masterEdgeI)
1726 const edge& masterEdge = masterEdges[masterEdgeI];
1728 label cutPoint0 = masterToCutPoints_[masterEdge[0]];
1729 label cutPoint1 = masterToCutPoints_[masterEdge[1]];
1733 label cutPointi = cutPoint0;
1766 Pout<<
"Dumping unmatched pointEdges to errorEdges.obj"
1777 cutFaces().pointEdges()[cutPointi]
1780 cutFaces().localPoints(),
1785 <<
"Problem in finding cut edges corresponding to"
1786 <<
" master edge " << masterEdge
1787 <<
" points:" << masterPoints[masterEdge[0]]
1788 <<
' ' << masterPoints[masterEdge[1]]
1789 <<
" corresponding cut edge: (" << cutPoint0
1790 <<
") " << cutPoint1
1794 cutToMasterEdges[cutEdgeI] = masterEdgeI;
1796 cutPointi = cutEdges[cutEdgeI].otherVertex(cutPointi);
1798 }
while (cutPointi != cutPoint1);
1804 writeEdges(cutToMasterEdges, cutToSlaveEdges);
1809 setCutEdgeToPoints(cutToMasterEdges);
1823 Pout<<
"subDivisionMatch :"
1824 <<
" matching (topological) cut faces to masterPatch" <<
endl;
1828 Map<labelList> candidates(cutFaces().size());
1830 cutToMasterFaces_.setSize(cutFaces().size());
1831 cutToMasterFaces_ = -1;
1837 label nChanged = matchEdgeFaces(cutToMasterEdges, candidates);
1839 checkMatch(cutToMasterEdges);
1846 nChanged += growCutFaces(cutToMasterEdges, candidates);
1848 checkMatch(cutToMasterEdges);
1858 Pout<<
"subDivisionMatch :"
1859 <<
" matching (geometric) cut faces to masterPatch" <<
endl;
1868 label nChanged = geometricMatchEdgeFaces(candidates);
1870 checkMatch(cutToMasterEdges);
1872 nChanged += growCutFaces(cutToMasterEdges, candidates);
1874 checkMatch(cutToMasterEdges);
1884 forAll(cutToMasterFaces_, cutFacei)
1886 if (cutToMasterFaces_[cutFacei] == -1)
1888 const face& cutF = cutFaces()[cutFacei];
1891 <<
"Did not match all of cutFaces to a master face" <<
nl
1892 <<
"First unmatched cut face:" << cutFacei <<
" with points:"
1893 << UIndirectList<point>(cutFaces().
points(), cutF)
1895 <<
"This usually means that the slave patch is not a"
1896 <<
" subdivision of the master patch"
1903 Pout<<
"subDivisionMatch :"
1904 <<
" finished matching master and slave to cut" <<
endl;
1920 const scalar absTol,
1921 const bool perfectMatch
1924 masterPatchPtr_(nullptr),
1925 slavePatchPtr_(nullptr),
1927 cutFacesPtr_(nullptr),
1928 cutToMasterFaces_(0),
1929 masterToCutPoints_(0),
1930 cutToSlaveFaces_(0),
1931 slaveToCutPoints_(0),
1943 findPerfectMatchingFaces
1956 findSlavesCoveringMaster
1967 masterPatchPtr_.reset
1976 slavePatchPtr_.reset
1989 perfectPointMatch(absTol,
false);
1994 subDivisionMatch(slaveMesh,
false, absTol);
2010 const scalar absTol,
2011 const bool perfectMatch,
2012 const bool orderedFaces,
2013 const bool patchDivision
2033 cutFacesPtr_(nullptr),
2034 cutToMasterFaces_(0),
2035 masterToCutPoints_(0),
2036 cutToSlaveFaces_(0),
2037 slaveToCutPoints_(0),
2040 if (perfectMatch && (masterAddressing.
size() != slaveAddressing.
size()))
2043 <<
"Perfect match specified but number of master and slave faces"
2044 <<
" differ." <<
endl
2045 <<
"master:" << masterAddressing.
size()
2046 <<
" slave:" << slaveAddressing.
size()
2052 masterAddressing.
size()
2057 <<
"Supplied internal face on master mesh to couple." <<
nl
2058 <<
"Faces to be coupled have to be boundary faces."
2063 slaveAddressing.
size()
2068 <<
"Supplied internal face on slave mesh to couple." <<
nl
2069 <<
"Faces to be coupled have to be boundary faces."
2076 perfectPointMatch(absTol, orderedFaces);
2081 subDivisionMatch(slaveMesh, patchDivision, absTol);
2097 label facei = pp.
start();
bool insert(const Key &key, const T &obj)
Copy insert a new entry, not overwriting existing entries.
A List with indirect addressing.
void append(const T &val)
Append an element at the end of the list.
void resize(const label len)
Adjust allocated size of list.
A HashTable to objects of type <T> with a label key.
virtual const fileName & name() const
Get the name of the stream.
A list of faces which address into the list of points.
void size(const label n)
Older name for setAddressableSize.
Container for information needed to couple to meshes. When constructed from two meshes and a geometri...
static Map< label > makeMap(const labelList &)
Create Map from List.
A face is a list of labels corresponding to mesh vertices.
Mesh consisting of general polyhedral cells.
virtual const faceList & faces() const
Return raw faces.
virtual const pointField & points() const
Return raw points.
A patch is a list of labels that address the faces in the global face list.
label start() const
Return start label of this patch in the polyMesh face list.
label nInternalFaces() const noexcept
Number of internal faces.
scalar maxDistance() const
Highest distance of all features.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
const labelList & faceLabels() const noexcept
Values for "faces" (XML only)
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const labelList nEdges(UPstream::listGatherValues< label >(aMesh.nEdges()))
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Determine correspondence between points. See below.
#define WarningInFunction
Report a warning using Foam::Warning.
const dimensionedScalar c
Speed of light in a vacuum.
const std::string patch
OpenFOAM patch number as a std::string.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i)
void inplaceRenumber(const labelUList &oldToNew, IntListType &input)
Inplace renumber the values (not the indices) of a list.
List< label > labelList
A List of labels.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
PointIndexHit< point > pointIndexHit
A PointIndexHit for 3D points.
IntListType renumber(const labelUList &oldToNew, const IntListType &input)
Renumber the values (not the indices) of a list.
vectorField pointField
pointField is a vectorField.
vector point
Point is a vector.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
List< labelList > labelListList
A List of labelList.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
errorManip< error > abort(error &err)
PointHit< point > pointHit
A PointIndexHit for 3D points.
static constexpr const zero Zero
Global zero (0)
bool matchPoints(const UList< point > &pts0, const UList< point > &pts1, const UList< scalar > &matchDistance, const bool verbose, labelList &from0To1, const point &origin=point::zero)
Determine correspondence between pointFields. Gets passed.
HashSet< label, Hash< label > > labelHashSet
A HashSet of labels, uses label hasher.
PrimitivePatch< List< face >, const pointField & > primitiveFacePatch
A PrimitivePatch with List storage for the faces, const reference for the point field.
List< edge > edgeList
A List of edges.
labelList invert(const label len, const labelUList &map)
Create an inverse one-to-one mapping.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
labelListList invertOneToMany(const label len, const labelUList &map)
Invert one-to-many map. Unmapped elements will be size 0.
List< face > faceList
A List of faces.
PrimitivePatch< IndirectList< face >, const pointField & > indirectPrimitivePatch
A PrimitivePatch with an IndirectList for the faces, const reference for the point field.
ListType reorder(const labelUList &oldToNew, const ListType &input, const bool prune=false)
Reorder the elements of a list.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
constexpr char nl
The newline '\n' character (0x0a)
#define forAll(list, i)
Loop across all elements in list.
#define forAllConstIters(container, iter)
Iterate across all elements of the container object with const access.
pointField points0(pointIOField(IOobject("points", mesh.time().constant(), polyMesh::meshSubDir, mesh, IOobject::MUST_READ, IOobject::NO_WRITE, false)))