47 scalar edgeIntersections::alignedCos_ =
cos(
degToRad(89.0));
53 void Foam::edgeIntersections::checkEdges(
const triSurface& surf)
55 const pointField& localPoints = surf.localPoints();
56 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" 80 if (edgeFaces[edgeI].size() == 1)
83 <<
"Edge " << edgeI <<
" vertices " << e
84 <<
" coords:" << localPoints[e[0]] <<
' ' 85 << localPoints[e[1]] <<
" has only one face connected to it:" 86 << edgeFaces[edgeI] <<
endl 87 <<
"This might lead to problems in intersection" 95 void Foam::edgeIntersections::intersectEdges
97 const triSurface& surf1,
99 const triSurfaceSearch& querySurf2,
104 const triSurface& surf2 = querySurf2.surface();
107 const labelList& meshPoints = surf1.meshPoints();
111 Pout<<
"Calculating intersection of " << edgeLabels.size() <<
" edges" 112 <<
" out of " << surf1.nEdges() <<
" with " << surf2.size()
113 <<
" triangles ..." <<
endl;
123 label edgeI = edgeLabels[i];
127 Pout<<
"Intersecting edge " << edgeI <<
" with surface" <<
endl;
130 const edge& e = surf1.edges()[edgeI];
132 const point& pStart = points1[meshPoints[e.start()]];
133 const point& pEnd = points1[meshPoints[e.end()]];
135 const vector eVec(pEnd - pStart);
143 start[i] = pStart - 0.5*surf1PointTol[e[0]]*
n;
144 end[i] = pEnd + 0.5*surf1PointTol[e[1]]*
n;
149 List<List<pointIndexHit>> edgeIntersections;
150 querySurf2.findLineAll
162 const label edgeI = edgeLabels[i];
164 labelList& intersectionTypes = classification_[edgeI];
165 intersectionTypes.
setSize(edgeIntersections[i].size(), -1);
167 this->operator[](edgeI).transfer(edgeIntersections[i]);
169 forAll(intersectionTypes, hitI)
172 label& hitType = intersectionTypes[hitI];
179 const edge& e = surf1.edges()[edgeI];
182 if (
mag(pHit.hitPoint() - start[i]) < surf1PointTol[e[0]])
187 else if (
mag(pHit.hitPoint() - end[i]) < surf1PointTol[e[1]])
192 else if (
mag(edgeDirs[i] & normals2[pHit.index()]) < alignedCos_)
196 Pout<<
"Flat angle edge:" << edgeI
197 <<
" face:" << pHit.index()
198 <<
" cos:" <<
mag(edgeDirs[i] & normals2[pHit.index()])
206 Info<<
" hit " << pHit <<
" classify = " << hitType <<
endl;
215 Pout<<
"Found " << nHits <<
" intersections of edges with surface ..." 227 bool Foam::edgeIntersections::inlinePerturb
229 const triSurface& surf1,
237 bool hasPerturbed =
false;
243 const labelList& edgeEnds = classification_[edgeI];
247 bool perturbStart =
false;
248 bool perturbEnd =
false;
251 if (edgeEnds.first() == 0)
256 if (edgeEnds.last() == 1)
262 if (perturbStart || perturbEnd)
264 const edge& e = surf1.edges()[edgeI];
266 label v0 = surf1.meshPoints()[e[0]];
267 label v1 = surf1.meshPoints()[e[1]];
269 vector eVec(points1[v1] - points1[v0]);
275 scalar t = 4.0*(rndGen.scalar01() - 0.5);
276 points1[v0] += t*surf1PointTol[e[0]]*
n;
278 const labelList& pEdges = surf1.pointEdges()[e[0]];
282 affectedEdges[pEdges[i]] =
true;
288 scalar t = 4.0*(rndGen.scalar01() - 0.5);
289 points1[v1] += t*surf1PointTol[e[1]]*
n;
291 const labelList& pEdges = surf1.pointEdges()[e[1]];
295 affectedEdges[pEdges[i]] =
true;
307 bool Foam::edgeIntersections::rotatePerturb
309 const triSurface& surf1,
318 const labelList& meshPoints = surf1.meshPoints();
320 const labelList& edgeEnds = classification_[edgeI];
322 bool hasPerturbed =
false;
326 if (edgeEnds[i] == 2)
328 const edge& e = surf1.edges()[edgeI];
331 label pointi = e[rndGen.sampleAB<
label>(0, 2)];
338 vector n(points1[meshPoints[e[1]]] - points1[meshPoints[e[0]]]);
339 scalar magN =
mag(n) + vSmall;
342 rndVec -= n*(n & rndVec);
345 rndVec /=
mag(rndVec) + vSmall;
350 Pout<<
"rotating: shifting endpoint " << meshPoints[pointi]
351 <<
" of edge:" << edgeI <<
" verts:" 352 << points1[meshPoints[e[0]]] <<
' ' 353 << points1[meshPoints[e[1]]]
355 <<
" tol:" << surf1PointTol[pointi] <<
endl;
357 points1[meshPoints[pointi]] += rndVec;
360 const labelList& pEdges = surf1.pointEdges()[pointi];
364 affectedEdges[pEdges[i]] =
true;
380 bool Foam::edgeIntersections::offsetPerturb
382 const triSurface& surf1,
383 const triSurface& surf2,
391 const labelList& meshPoints = surf1.meshPoints();
393 const edge& e = surf1.edges()[edgeI];
395 const List<pointIndexHit>& hits = operator[](edgeI);
398 bool hasPerturbed =
false;
406 label surf2Facei = pHit.index();
409 const pointField& surf2Pts = surf2.localPoints();
413 label nearType, nearLabel;
415 f2.nearestPointClassify(pHit.hitPoint(), surf2Pts, nearType, nearLabel);
420 vector offset = 0.01*rndGen.scalar01()*(ctr - pHit.hitPoint());
423 points1[meshPoints[e[0]]] += offset;
426 const labelList& pEdges0 = surf1.pointEdges()[e[0]];
430 affectedEdges[pEdges0[i]] =
true;
434 points1[meshPoints[e[1]]] += offset;
437 const labelList& pEdges1 = surf1.pointEdges()[e[1]];
441 affectedEdges[pEdges1[i]] =
true;
474 classification_(surf1.
nEdges())
508 classification_(classification)
524 const labelList& pEdges = pointEdges[pointi];
530 minDist =
min(minDist, edges[pEdges[i]].
mag(localPoints));
566 for (; iter < nIters; iter++)
578 label edgeI = edgesToTest[i];
581 if (!affectedEdges[edgeI])
585 bool shiftedEdgeEndPoints =
596 nShifted += (shiftedEdgeEndPoints ? 1 : 0);
598 if (!shiftedEdgeEndPoints)
611 nRotated += (rotatedEdge ? 1 : 0);
619 bool offsetEdgePoints =
630 nOffset += (offsetEdgePoints ? 1 : 0);
638 Pout<<
"Edges to test : " <<
nl 639 <<
" total:" << edgesToTest.size() <<
nl 640 <<
" resolved by:" <<
nl 641 <<
" shifting : " << nShifted <<
nl 642 <<
" rotating : " << nRotated <<
nl 643 <<
" offsetting : " << nOffset <<
nl 648 if (nShifted == 0 && nRotated == 0 && nOffset == 0)
659 forAll(affectedEdges, edgeI)
661 if (affectedEdges[edgeI])
663 newEdgesToTest[newEdgeI++] = edgeI;
666 newEdgesToTest.setSize(newEdgeI);
670 Pout<<
"Edges to test:" <<
nl 671 <<
" was : " << edgesToTest.size() <<
nl 672 <<
" is : " << newEdgesToTest.size() <<
nl 677 edgesToTest.transfer(newEdgesToTest);
679 if (edgesToTest.empty())
701 const edgeIntersections& subInfo,
712 label edgeI = edgeMap[subI];
714 labelList& intersectionTypes = classification_[edgeI];
720 sz = intersections.
size();
728 bool foundFace =
false;
729 for (
label interI = 0; interI < sz; interI++)
731 if (intersections[interI].index() == faceMap[subHit.
index()])
744 intersections.
setSize(sz+nNew);
745 intersectionTypes.
setSize(sz+nNew);
752 bool foundFace =
false;
753 for (
label interI = 0; interI < sz; interI++)
755 if (intersections[interI].index() == faceMap[subHit.
index()])
768 faceMap[subHit.
index()]
770 intersectionTypes[nNew] = subClass[i];
const labelListList & pointEdges() const
Return point-edge addressing.
List< labelList > labelListList
A List of labelList.
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
T & operator[](const label)
Return element of UList.
#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.
void size(const label)
Override size to be inconsistent with allocated storage.
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.
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
PointIndexHit< point > pointIndexHit
Vector< scalar > vector
A scalar version of the templated Vector.
const Field< PointType > & localPoints() const
Return pointField of points in patch.
scalar minDist(const List< pointIndexHit > &hitList)
scalar degToRad(const scalar deg)
Conversion from degrees to radians.
List< bool > boolList
Bool container classes.
const Vector< Cmpt > & centre(const Foam::List< Vector< Cmpt >> &) const
Return *this (used for point which is a typedef to Vector<scalar>.
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.
dimensionedScalar cos(const dimensionedScalar &ds)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
bool hit() const
Is there a hit.
const Field< PointType > & points() const
Return reference to global points.
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
List< label > labelList
A List of labels.
const labelListList & classification() const
For every intersection the classification status.
errorManip< error > abort(error &err)
const Point & rawPoint() const
Return point with no checking.
label nEdges() const
Return number of edges in patch.
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.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
const triSurface & surface() const
Return reference to the surface.
void setSize(const label)
Reset size of List.
std::remove_reference< ::Foam::List< labelledTri > >::type::value_type FaceType
vector point
Point is a vector.
#define WarningInFunction
Report a warning using Foam::Warning.
prefixOSstream Pout(cout, "Pout")
dimensioned< scalar > mag(const dimensioned< Type > &)
Field< vector > vectorField
Specialisation of Field<T> for vector.
Triangulated surface description with patch information.
edgeIntersections()
Construct null.
label index() const
Return index.