36 const Foam::scalar Foam::slidingInterface::pointMergeTolDefault_ = 0.05;
37 const Foam::scalar Foam::slidingInterface::edgeMergeTolDefault_ = 0.01;
38 const Foam::label Foam::slidingInterface::nFacesPerSlaveEdgeDefault_ = 5;
39 const Foam::label Foam::slidingInterface::edgeFaceEscapeLimitDefault_ = 10;
41 const Foam::scalar Foam::slidingInterface::integralAdjTolDefault_ = 0.05;
43 Foam::slidingInterface::edgeMasterCatchFractionDefault_ = 0.4;
44 const Foam::scalar Foam::slidingInterface::edgeEndCutoffTolDefault_ = 0.0001;
60 bool Foam::slidingInterface::projectPoints()
const 65 <<
": for object " <<
name() <<
" : " 66 <<
"Projecting slave points onto master surface." <<
endl;
96 mesh.faceZones()[masterFaceZoneID_.
index()]();
99 mesh.faceZones()[slaveFaceZoneID_.
index()]();
104 const pointField& masterLocalPoints = masterPatch.localPoints();
105 const faceList& masterLocalFaces = masterPatch.localFaces();
106 const edgeList& masterEdges = masterPatch.edges();
107 const labelListList& masterEdgeFaces = masterPatch.edgeFaces();
108 const labelListList& masterFaceEdges = masterPatch.faceEdges();
109 const labelListList& masterFaceFaces = masterPatch.faceFaces();
116 const pointField& slaveLocalPoints = slavePatch.localPoints();
117 const edgeList& slaveEdges = slavePatch.edges();
118 const labelListList& slaveEdgeFaces = slavePatch.edgeFaces();
119 const vectorField& slavePointNormals = slavePatch.pointNormals();
129 scalarField minMasterPointLength(masterLocalPoints.size(), GREAT);
130 scalarField minMasterFaceLength(masterPatch.size(), GREAT);
132 forAll(masterEdges, edgei)
134 const edge& curEdge = masterEdges[edgei];
136 const scalar curLength = curEdge.mag(masterLocalPoints);
139 minMasterPointLength[curEdge.start()] =
142 minMasterPointLength[curEdge.start()],
146 minMasterPointLength[curEdge.end()] =
149 minMasterPointLength[curEdge.end()],
154 const labelList& curFaces = masterEdgeFaces[edgei];
156 for (
const label facei : curFaces)
158 minMasterFaceLength[facei] =
161 minMasterFaceLength[facei],
171 scalarField minSlavePointLength(slaveLocalPoints.size(), GREAT);
172 scalarField minSlaveFaceLength(slavePatch.size(), GREAT);
176 const edge& curEdge = slaveEdges[edgei];
178 const scalar curLength = curEdge.mag(slaveLocalPoints);
181 minSlavePointLength[curEdge.start()] =
184 minSlavePointLength[curEdge.start()],
188 minSlavePointLength[curEdge.end()] =
191 minSlavePointLength[curEdge.end()],
196 const labelList& curFaces = slaveEdgeFaces[edgei];
198 for (
const label facei : curFaces)
200 minSlaveFaceLength[facei] =
203 minSlaveFaceLength[facei],
215 List<objectHit> slavePointFaceHits =
216 slavePatch.projectPoints
227 for (
const auto& hit : slavePointFaceHits)
235 Pout<<
"Number of hits in point projection: " << nHits
236 <<
" out of " << slavePointFaceHits.size() <<
" points." 241 projectedSlavePointsPtr_.reset
245 auto& projectedSlavePoints = *projectedSlavePointsPtr_;
249 label nAdjustedPoints = 0;
257 Pout<<
"bool slidingInterface::projectPoints() for object " 259 <<
"Adjusting point projection for integral match: ";
262 forAll(slavePointFaceHits, pointi)
264 if (slavePointFaceHits[pointi].hit())
267 projectedSlavePoints[pointi] =
269 [slavePointFaceHits[pointi].hitObject()].ray
271 slaveLocalPoints[pointi],
272 slavePointNormals[pointi],
281 masterLocalFaces[slavePointFaceHits[pointi].hitObject()].ray
283 slaveLocalPoints[pointi],
284 slavePointNormals[pointi],
289 const point nearPoint = missAdjust.missPoint();
290 const point missPoint =
291 slaveLocalPoints[pointi]
292 + slavePointNormals[pointi]*missAdjust.distance();
295 const scalar mergeTol =
296 integralAdjTol_*minSlavePointLength[pointi];
299 if (nearPoint.dist(missPoint) < mergeTol)
312 projectedSlavePoints[pointi] = nearPoint;
314 slavePointFaceHits[pointi] =
315 objectHit(
true, slavePointFaceHits[pointi].hitObject());
321 projectedSlavePoints[pointi] = slaveLocalPoints[pointi];
336 else if (matchType_ ==
PARTIAL)
338 forAll(slavePointFaceHits, pointi)
340 if (slavePointFaceHits[pointi].hit())
343 projectedSlavePoints[pointi] =
345 [slavePointFaceHits[pointi].hitObject()].ray
347 slaveLocalPoints[pointi],
348 slavePointNormals[pointi],
356 projectedSlavePoints[pointi] = slaveLocalPoints[pointi];
363 <<
" for object " <<
name()
369 Pout<<
"Number of adjusted points in projection: " 370 << nAdjustedPoints <<
endl;
373 scalar minEdgeLength = GREAT;
374 label nShortEdges = 0;
376 for (
const edge&
e : slaveEdges)
378 const scalar len =
e.mag(projectedSlavePoints);
379 minEdgeLength =
min(minEdgeLength, len);
383 Pout<<
"Point projection problems for edge: " 384 <<
e <<
". Length = " << len
394 << nShortEdges <<
" short projected edges " 395 <<
"after adjustment for object " <<
name()
400 Pout<<
" ... projection OK." <<
endl;
411 labelList slavePointPointHits(slaveLocalPoints.size(), -1);
412 labelList masterPointPointHits(masterLocalPoints.size(), -1);
425 label nMergedPoints = 0;
427 forAll(projectedSlavePoints, pointi)
429 if (slavePointFaceHits[pointi].hit())
432 point& curPoint = projectedSlavePoints[pointi];
435 const face& hitFace =
436 masterLocalFaces[slavePointFaceHits[pointi].hitObject()];
438 label mergePoint = -1;
439 scalar mergeDist = GREAT;
442 for (
const label hitPointi : hitFace)
445 mag(masterLocalPoints[hitPointi] - curPoint);
448 const scalar mergeTol =
452 minSlavePointLength[pointi],
453 minMasterPointLength[hitPointi]
456 if (dist < mergeTol && dist < mergeDist)
458 mergePoint = hitPointi;
479 slavePointPointHits[pointi] = mergePoint;
480 masterPointPointHits[mergePoint] = pointi;
483 curPoint = masterLocalPoints[mergePoint];
494 scalar minEdgeLength = GREAT;
496 for (
const edge&
e : slaveEdges)
498 const scalar len =
e.mag(projectedSlavePoints);
499 minEdgeLength =
min(minEdgeLength, len);
503 Pout<<
"Point projection problems for edge: " 504 <<
e <<
". Length = " << len
509 if (minEdgeLength < SMALL)
512 <<
" after point merge for object " <<
name()
517 Pout<<
" ... point merge OK." <<
endl;
523 labelList slavePointEdgeHits(slaveLocalPoints.size(), -1);
525 label nMovedPoints = 0;
527 forAll(projectedSlavePoints, pointi)
530 if (slavePointPointHits[pointi] < 0)
533 point& curPoint = projectedSlavePoints[pointi];
537 masterFaceEdges[slavePointFaceHits[pointi].hitObject()];
540 const scalar mergeLength =
543 minSlavePointLength[pointi],
544 minMasterFaceLength[slavePointFaceHits[pointi].hitObject()]
547 const scalar mergeTol = pointMergeTol_*mergeLength;
549 scalar minDistance = GREAT;
551 for (
const label edgei : hitFaceEdges)
553 const edge& curEdge = masterEdges[edgei];
556 curEdge.line(masterLocalPoints).nearestDist(curPoint);
561 edgeHit.point().dist(projectedSlavePoints[pointi]);
563 if (dist < mergeTol && dist < minDistance)
568 slavePointEdgeHits[pointi] = edgei;
569 projectedSlavePoints[pointi] = edgeHit.point();
588 if (slavePointEdgeHits[pointi] > -1)
592 point& curPoint = projectedSlavePoints[pointi];
594 const edge& hitMasterEdge =
595 masterEdges[slavePointEdgeHits[pointi]];
597 label mergePoint = -1;
598 scalar mergeDist = GREAT;
600 forAll(hitMasterEdge, hmeI)
602 const scalar hmeDist =
603 mag(masterLocalPoints[hitMasterEdge[hmeI]] - curPoint);
606 const scalar mergeTol =
610 minSlavePointLength[pointi],
611 minMasterPointLength[hitMasterEdge[hmeI]]
614 if (hmeDist < mergeTol && hmeDist < mergeDist)
616 mergePoint = hitMasterEdge[hmeI];
633 slavePointPointHits[pointi] = mergePoint;
634 curPoint = masterLocalPoints[mergePoint];
635 masterPointPointHits[mergePoint] = pointi;
637 slavePointFaceHits[pointi] =
638 objectHit(
true, slavePointFaceHits[pointi].hitObject());
642 slavePointEdgeHits[pointi] = -1;
650 Pout<<
"Number of merged master points: " << nMergedPoints <<
nl 651 <<
"Number of adjusted slave points: " << nMovedPoints <<
endl;
654 scalar minEdgeLength = GREAT;
656 for (
const edge&
e : slaveEdges)
658 const scalar len =
e.mag(projectedSlavePoints);
659 minEdgeLength =
min(minEdgeLength, len);
663 Pout<<
"Point projection problems for edge: " 664 <<
e <<
". Length = " << len
669 if (minEdgeLength < SMALL)
672 <<
" after correction for object " <<
name()
718 labelList masterPointEdgeHits(masterLocalPoints.size(), -1);
719 scalarField masterPointEdgeDist(masterLocalPoints.size(), GREAT);
729 Pout<<
"Processing slave edges " <<
endl;
742 const edge& curEdge = slaveEdges[edgei];
746 const label startFace =
747 slavePointFaceHits[curEdge.start()].hitObject();
749 const label endFace =
750 slavePointFaceHits[curEdge.end()].hitObject();
761 bool completed =
false;
772 curFaceMap.insert(startFace);
773 addedFaces.insert(startFace);
778 nSweeps < edgeFaceEscapeLimit_;
782 completed = addedFaces.found(endFace);
788 for (
const label cfi : cf)
790 const labelList& curNbrs = masterFaceFaces[cfi];
792 for (
const label nbri : curNbrs)
794 if (curFaceMap.insert(nbri))
796 addedFaces.insert(nbri);
801 if (completed)
break;
823 curFaceMap.insert(endFace);
824 addedFaces.insert(endFace);
829 nSweeps < edgeFaceEscapeLimit_;
833 completed = addedFaces.found(startFace);
839 for (
const label cfi : cf)
841 const labelList& curNbrs = masterFaceFaces[cfi];
843 for (
const label nbri : curNbrs)
845 if (curFaceMap.insert(nbri))
847 addedFaces.insert(nbri);
852 if (completed)
break;
884 for (
const label facei : curFaceMap)
886 const face&
f = masterLocalFaces[facei];
887 curPointMap.insert(
f);
890 const labelList curMasterPoints = curPointMap.toc();
894 linePointRef edgeLine = curEdge.line(projectedSlavePoints);
896 const vector edgeVec = edgeLine.vec();
897 const scalar edgeMag = edgeLine.
mag();
903 const scalar slaveCatchDist =
904 edgeMasterCatchFraction_*edgeMag
909 projectedSlavePoints[curEdge.start()]
910 - slaveLocalPoints[curEdge.start()]
914 projectedSlavePoints[curEdge.end()]
915 - slaveLocalPoints[curEdge.end()]
927 const vector edgeNormalInPlane =
932 slavePointNormals[curEdge.start()]
933 + slavePointNormals[curEdge.end()]
937 for (
const label cmp : curMasterPoints)
943 slavePointPointHits[curEdge.start()] == cmp
944 || slavePointPointHits[curEdge.end()] == cmp
945 || masterPointPointHits[cmp] > -1
954 edgeLine.nearestDist(masterLocalPoints[cmp]);
956 if (edgeLineHit.hit())
964 const scalar cutOnSlave =
965 ((edgeLineHit.point() - edgeLine.start()) & edgeVec)
968 const scalar distInEdgePlane =
971 edgeLineHit.distance(),
975 masterLocalPoints[cmp]
976 - edgeLineHit.point()
996 cutOnSlave > edgeEndCutoffTol_
997 && cutOnSlave < 1.0 - edgeEndCutoffTol_
998 && distInEdgePlane < edgeMergeTol_*edgeMag
999 && edgeLineHit.distance()
1003 masterPointEdgeDist[cmp]
1009 if (masterPointEdgeHits[cmp] == -1)
1022 masterPointEdgeHits[cmp] = edgei;
1023 masterPointEdgeDist[cmp] = edgeLineHit.distance();
1057 Pout<<
"bool slidingInterface::projectPoints() for object " 1059 <<
"Finished projecting points. Topology = ";
1085 slavePointPointHitsPtr_.reset
1087 new labelList(std::move(slavePointPointHits))
1089 slavePointEdgeHitsPtr_.reset
1091 new labelList(std::move(slavePointEdgeHits))
1093 slavePointFaceHitsPtr_.reset
1095 new List<objectHit>(std::move(slavePointFaceHits))
1097 masterPointEdgeHitsPtr_.reset
1099 new labelList(std::move(masterPointEdgeHits))
1104 Pout<<
"(Detached interface) changing." <<
endl;
1115 !slavePointPointHitsPtr_
1116 || !slavePointEdgeHitsPtr_
1117 || !slavePointFaceHitsPtr_
1118 || !masterPointEdgeHitsPtr_
1122 slavePointPointHitsPtr_.reset
1126 slavePointEdgeHitsPtr_.reset
1128 new labelList(std::move(slavePointEdgeHits))
1130 slavePointFaceHitsPtr_.reset
1132 new List<objectHit>(std::move(slavePointFaceHits))
1134 masterPointEdgeHitsPtr_.reset
1136 new labelList(std::move(masterPointEdgeHits))
1141 Pout<<
"(Attached interface restart) changing." <<
endl;
1148 if (slavePointPointHits != *slavePointPointHitsPtr_)
1152 Pout<<
"(Point projection) ";
1158 if (slavePointEdgeHits != *slavePointEdgeHitsPtr_)
1162 Pout<<
"(Edge projection) ";
1170 bool faceHitsDifferent =
false;
1172 const List<objectHit>& oldPointFaceHits = *slavePointFaceHitsPtr_;
1174 forAll(slavePointFaceHits, pointi)
1178 slavePointPointHits[pointi] < 0
1179 && slavePointEdgeHits[pointi] < 0
1183 if (slavePointFaceHits[pointi] != oldPointFaceHits[pointi])
1186 faceHitsDifferent =
true;
1192 if (faceHitsDifferent)
1196 Pout<<
"(Face projection) ";
1203 if (masterPointEdgeHits != *masterPointEdgeHitsPtr_)
1207 Pout<<
"(Master point projection) ";
1217 slavePointPointHitsPtr_.reset
1219 new labelList(std::move(slavePointPointHits))
1221 slavePointEdgeHitsPtr_.reset
1223 new labelList(std::move(slavePointEdgeHits))
1225 slavePointFaceHitsPtr_.reset
1227 new List<objectHit>(std::move(slavePointFaceHits))
1229 masterPointEdgeHitsPtr_.reset
1231 new labelList(std::move(masterPointEdgeHits))
1252 void Foam::slidingInterface::clearPointProjection()
const 1254 slavePointPointHitsPtr_.reset(
nullptr);
1255 slavePointEdgeHitsPtr_.reset(
nullptr);
1256 slavePointFaceHitsPtr_.reset(
nullptr);
1257 masterPointEdgeHitsPtr_.reset(
nullptr);
1259 projectedSlavePointsPtr_.reset(
nullptr);
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
List< edge > edgeList
List of edge.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
constexpr char nl
The newline '\n' character (0x0a)
Ostream & endl(Ostream &os)
Add newline and flush stream.
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
List< labelList > labelListList
List of labelList.
static const unsigned edgesPerFace_
Estimated number of edges per cell.
#define forAll(list, i)
Loop across all elements in list.
static const unsigned pointsPerFace_
Estimated number of points per face.
HashSet< label, Hash< label > > labelHashSet
A HashSet of labels, uses label hasher.
List< face > faceList
List of faces.
vectorField pointField
pointField is a vectorField.
const dimensionedScalar e
Elementary charge.
line< point, const point & > linePointRef
A line using referred points.
void clear()
Clear the list, i.e. set size to zero.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
errorManip< error > abort(error &err)
label index() const
The index of the first matching items, -1 if no matches.
int debug
Static debugging option.
scalar mag() const
The length (L2-norm) of the vector.
PrimitivePatch< List< face >, const pointField & > primitiveFacePatch
A PrimitivePatch with List storage for the faces, const reference for the point field.
const polyTopoChanger & topoChanger() const
Return reference to morph engine.
vector point
Point is a vector.
const word & name() const
Return name of this modifier.
Field< vector > vectorField
Specialisation of Field<T> for vector.
List< label > labelList
A List of labels.
const polyMesh & mesh() const
Return the mesh reference.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
PointHit< point > pointHit
A PointHit with a 3D point.
static constexpr const zero Zero
Global zero (0)