54 void Foam::edgeIntersections::checkEdges(
const triSurface& surf)
56 const pointField& localPoints = surf.localPoints();
57 const edgeList& edges = surf.edges();
59 treeBoundBox bb(localPoints);
61 scalar minSize = SMALL * bb.minDim();
65 const edge&
e = edges[edgeI];
67 scalar eMag =
e.mag(localPoints);
72 <<
"Edge " << edgeI <<
" vertices " <<
e 73 <<
" coords:" << localPoints[
e[0]] <<
' ' 74 << localPoints[
e[1]] <<
" is very small compared to bounding" 75 <<
" box dimensions " << bb <<
endl 76 <<
"This might lead to problems in intersection" 84 void Foam::edgeIntersections::intersectEdges
86 const triSurface& surf1,
88 const triSurfaceSearch& querySurf2,
93 const triSurface& surf2 = querySurf2.surface();
96 const labelList& meshPoints = surf1.meshPoints();
100 Pout<<
"Calculating intersection of " << edgeLabels.size() <<
" edges" 101 <<
" out of " << surf1.nEdges() <<
" with " << surf2.size()
102 <<
" triangles ..." <<
endl;
112 label edgeI = edgeLabels[i];
116 Pout<<
"Intersecting edge " << edgeI <<
" with surface" <<
endl;
119 const edge&
e = surf1.edges()[edgeI];
121 const point& pStart = points1[meshPoints[
e.start()]];
122 const point& pEnd = points1[meshPoints[
e.end()]];
131 start[i] = pStart - 0.5*surf1PointTol[
e[0]]*
n;
132 end[i] = pEnd + 0.5*surf1PointTol[
e[1]]*
n;
137 List<List<pointIndexHit>> edgeIntersections;
138 querySurf2.findLineAll
150 const label edgeI = edgeLabels[i];
152 labelList& intersectionTypes = classification_[edgeI];
153 intersectionTypes.
setSize(edgeIntersections[i].size(), -1);
155 this->operator[](edgeI).transfer(edgeIntersections[i]);
157 forAll(intersectionTypes, hitI)
160 label& hitType = intersectionTypes[hitI];
167 const edge&
e = surf1.edges()[edgeI];
170 if (pHit.point().dist(start[i]) < surf1PointTol[
e[0]])
175 else if (pHit.point().dist(
end[i]) < surf1PointTol[
e[1]])
180 else if (
mag(edgeDirs[i] & normals2[pHit.index()]) < alignedCos_)
184 Pout<<
"Flat angle edge:" << edgeI
185 <<
" face:" << pHit.index()
186 <<
" cos:" <<
mag(edgeDirs[i] & normals2[pHit.index()])
194 Info<<
" hit " << pHit <<
" classify = " << hitType <<
endl;
203 Pout<<
"Found " << nHits <<
" intersections of edges with surface ..." 215 bool Foam::edgeIntersections::inlinePerturb
217 const triSurface& surf1,
225 bool hasPerturbed =
false;
231 const labelList& edgeEnds = classification_[edgeI];
235 bool perturbStart =
false;
236 bool perturbEnd =
false;
239 if (edgeEnds.first() == 0)
244 if (edgeEnds.last() == 1)
250 if (perturbStart || perturbEnd)
252 const edge&
e = surf1.edges()[edgeI];
254 label v0 = surf1.meshPoints()[
e[0]];
255 label v1 = surf1.meshPoints()[
e[1]];
263 points1[v0] += t*surf1PointTol[
e[0]]*
n;
265 const labelList& pEdges = surf1.pointEdges()[
e[0]];
269 affectedEdges[pEdges[i]] =
true;
276 points1[v1] += t*surf1PointTol[
e[1]]*
n;
278 const labelList& pEdges = surf1.pointEdges()[
e[1]];
282 affectedEdges[pEdges[i]] =
true;
294 bool Foam::edgeIntersections::rotatePerturb
296 const triSurface& surf1,
305 const labelList& meshPoints = surf1.meshPoints();
307 const labelList& edgeEnds = classification_[edgeI];
309 bool hasPerturbed =
false;
313 if (edgeEnds[i] == 2)
315 const edge&
e = surf1.edges()[edgeI];
325 vector n(points1[meshPoints[
e[1]]] - points1[meshPoints[
e[0]]]);
326 scalar magN =
mag(
n) + VSMALL;
329 rndVec.removeCollinear(
n);
335 Pout<<
"rotating: shifting endpoint " << meshPoints[pointi]
336 <<
" of edge:" << edgeI <<
" verts:" 337 << points1[meshPoints[
e[0]]] <<
' ' 338 << points1[meshPoints[
e[1]]]
340 <<
" tol:" << surf1PointTol[pointi] <<
endl;
342 points1[meshPoints[pointi]] += rndVec;
345 const labelList& pEdges = surf1.pointEdges()[pointi];
349 affectedEdges[pEdges[i]] =
true;
365 bool Foam::edgeIntersections::offsetPerturb
367 const triSurface& surf1,
368 const triSurface& surf2,
376 const labelList& meshPoints = surf1.meshPoints();
378 const edge&
e = surf1.edges()[edgeI];
380 const List<pointIndexHit>& hits = operator[](edgeI);
383 bool hasPerturbed =
false;
389 const label surf2Facei = pHit.index();
392 const pointField& surf2Pts = surf2.localPoints();
394 const point ctr = f2.centre(surf2Pts);
396 label nearType, nearLabel;
398 f2.nearestPointClassify(pHit.hitPoint(), surf2Pts, nearType, nearLabel);
407 points1[meshPoints[
e[0]]] += offset;
410 const labelList& pEdges0 = surf1.pointEdges()[
e[0]];
414 affectedEdges[pEdges0[i]] =
true;
418 points1[meshPoints[
e[1]]] += offset;
421 const labelList& pEdges1 = surf1.pointEdges()[
e[1]];
425 affectedEdges[pEdges1[i]] =
true;
456 classification_(surf1.nEdges())
482 classification_(classification)
498 const labelList& pEdges = pointEdges[pointi];
500 scalar minDist = GREAT;
504 minDist =
min(minDist, edges[pEdges[i]].
mag(localPoints));
507 minLen[pointi] = minDist;
540 for (; iter < nIters; iter++)
552 label edgeI = edgesToTest[i];
555 if (!affectedEdges[edgeI])
559 bool shiftedEdgeEndPoints =
570 nShifted += (shiftedEdgeEndPoints ? 1 : 0);
572 if (!shiftedEdgeEndPoints)
585 nRotated += (rotatedEdge ? 1 : 0);
593 bool offsetEdgePoints =
604 nOffset += (offsetEdgePoints ? 1 : 0);
612 Pout<<
"Edges to test : " <<
nl 613 <<
" total:" << edgesToTest.size() <<
nl 614 <<
" resolved by:" <<
nl 615 <<
" shifting : " << nShifted <<
nl 616 <<
" rotating : " << nRotated <<
nl 617 <<
" offsetting : " << nOffset <<
nl 622 if (nShifted == 0 && nRotated == 0 && nOffset == 0)
633 forAll(affectedEdges, edgeI)
635 if (affectedEdges[edgeI])
637 newEdgesToTest[newEdgeI++] = edgeI;
640 newEdgesToTest.setSize(newEdgeI);
644 Pout<<
"Edges to test:" <<
nl 645 <<
" was : " << edgesToTest.size() <<
nl 646 <<
" is : " << newEdgesToTest.size() <<
nl 651 edgesToTest.transfer(newEdgesToTest);
653 if (edgesToTest.empty())
675 const edgeIntersections& subInfo,
683 const List<pointIndexHit>& subHits = subInfo[subI];
684 const labelList& subClass = subInfo.classification()[subI];
686 label edgeI = edgeMap[subI];
687 List<pointIndexHit>& intersections = operator[](edgeI);
688 labelList& intersectionTypes = classification_[edgeI];
694 sz = intersections.
size();
702 bool foundFace =
false;
703 for (label interI = 0; interI < sz; interI++)
705 if (intersections[interI].index() ==
faceMap[subHit.index()])
718 intersections.setSize(sz+nNew);
719 intersectionTypes.setSize(sz+nNew);
726 bool foundFace =
false;
727 for (label interI = 0; interI < sz; interI++)
729 if (intersections[interI].index() ==
faceMap[subHit.index()])
743 intersectionTypes[nNew] = subClass[i];
const labelListList & pointEdges() const
Return point-edge addressing.
void size(const label n)
Older name for setAddressableSize.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const Field< point_type > & localPoints() const
Return pointField of points in patch.
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.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Unit conversion functions.
constexpr char nl
The newline '\n' character (0x0a)
void merge(const edgeIntersections &, const labelList &edgeMap, const labelList &faceMap, const bool merge=true)
Merge (or override) edge intersection for a subset.
Ostream & endl(Ostream &os)
Add newline and flush stream.
PointIndexHit< point > pointIndexHit
A PointIndexHit with a 3D point.
This class describes the interaction of an object (often a face) and a point. It carries the info of ...
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
#define forAll(list, i)
Loop across all elements in list.
Helper class to search on triSurface.
static scalarField minEdgeLength(const triSurface &surf)
Calculate min edge length for every surface point.
vectorField pointField
pointField is a vectorField.
Type sample01()
Return a sample whose components lie in the range [0,1].
const dimensionedScalar e
Elementary charge.
void setSize(const label n)
Alias for resize()
dimensionedScalar cos(const dimensionedScalar &ds)
List< edge > edgeList
A List of edges.
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
labelledTri face_type
The face type (same as the underlying PrimitivePatch)
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
const Field< point_type > & points() const noexcept
Return reference to global points.
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 nEdges() const
Number of edges in patch.
int debug
Static debugging option.
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
defineTypeNameAndDebug(combustionModel, 0)
label removeDegenerates(const label nIters, const triSurface &surf1, const triSurfaceSearch &query2, const scalarField &surf1PointTol, pointField &points1)
Resolve ties. Shuffles points so all edge - face intersections.
const triSurface & surface() const
Return reference to the surface.
vector point
Point is a vector.
#define WarningInFunction
Report a warning using Foam::Warning.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Field< vector > vectorField
Specialisation of Field<T> for vector.
List< label > labelList
A List of labels.
static scalar alignedCos_
Cosine between edge and face normal when considered parallel.
Triangulated surface description with patch information.
constexpr scalar degToRad(const scalar deg) noexcept
Conversion from degrees to radians.
edgeIntersections()
Construct null.
List< bool > boolList
A List of bools.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
int bit()
Return a random bit.