38bool Foam::featurePointConformer::createSpecialisedFeaturePoint
64 if (debug)
Info<<
"nExternal == 2 && nInternal == 1" <<
endl;
77 label nVert = foamyHexMesh_.number_of_vertices();
79 const label initialNumOfPoints = pts.
size();
87 label concaveEdgeI = -1;
97 concaveEdgeI = pEds[i];
101 convexEdgesI[nConvex++] = pEds[i];
106 <<
"Edge " << eS <<
" is flat"
112 <<
"Edge " << eS <<
" not concave/convex"
117 const vector& concaveEdgePlaneANormal =
118 normals[edgeNormals[concaveEdgeI][0]];
120 const vector& concaveEdgePlaneBNormal =
121 normals[edgeNormals[concaveEdgeI][1]];
127 featPt + ppDist*concaveEdgePlaneANormal,
128 concaveEdgePlaneANormal
133 featPt + ppDist*concaveEdgePlaneBNormal,
134 concaveEdgePlaneBNormal
145 featPt + ppDist*concaveEdgeDir;
150 plane planeF(concaveEdgeLocalFeatPt, concaveEdgeDir);
152 const Foam::point concaveEdgeExternalPt = planeF.planePlaneIntersect
160 planeA =
plane(featPt, concaveEdgePlaneANormal);
162 planeB =
plane(featPt, concaveEdgePlaneBNormal);
165 concaveEdgeExternalPt
166 - 2.0*planeA.distance(concaveEdgeExternalPt)
167 *concaveEdgePlaneANormal;
180 const label internalPtAIndex(pts.
last().index());
183 concaveEdgeExternalPt
184 - 2.0*planeB.distance(concaveEdgeExternalPt)
185 *concaveEdgePlaneBNormal;
198 const label internalPtBIndex(pts.
last().index());
208 const labelList& concaveEdgeNormals = edgeNormals[concaveEdgeI];
209 const labelList& convexEdgeANormals = edgeNormals[convexEdgesI[0]];
210 const labelList& convexEdgeBNormals = edgeNormals[convexEdgesI[1]];
212 forAll(concaveEdgeNormals, edgeNormalI)
214 bool convexEdgeA =
false;
215 bool convexEdgeB =
false;
217 forAll(convexEdgeANormals, edgeAnormalI)
219 const vector& concaveNormal
220 = normals[concaveEdgeNormals[edgeNormalI]];
221 const vector& convexNormal
222 = normals[convexEdgeANormals[edgeAnormalI]];
226 Info<<
"Angle between vectors = "
232 if (
areParallel(concaveNormal, convexNormal, tolParallel))
238 forAll(convexEdgeBNormals, edgeBnormalI)
240 const vector& concaveNormal
241 = normals[concaveEdgeNormals[edgeNormalI]];
242 const vector& convexNormal
243 = normals[convexEdgeBNormals[edgeBnormalI]];
247 Info<<
"Angle between vectors = "
253 if (
areParallel(concaveNormal, convexNormal, tolParallel))
259 if ((convexEdgeA && convexEdgeB) || (!convexEdgeA && !convexEdgeB))
262 <<
"Both or neither of the convex edges share the concave "
264 <<
" convexEdgeA = " << convexEdgeA
265 <<
" convexEdgeB = " << convexEdgeB
269 for (label i = 0; i < 2; ++i)
280 forAll(convexEdgeANormals, edgeAnormalI)
282 const vector& concaveNormal
283 = normals[concaveEdgeNormals[edgeNormalI]];
284 const vector& convexNormal
285 = normals[convexEdgeANormals[edgeAnormalI]];
289 !
areParallel(concaveNormal, convexNormal, tolParallel)
292 convexEdgePlaneCNormal = convexNormal;
294 plane planeC(featPt, convexEdgePlaneCNormal);
298 + 2.0*planeC.distance(internalPtA)
299 *convexEdgePlaneCNormal;
323 forAll(convexEdgeBNormals, edgeBnormalI)
325 const vector& concaveNormal
326 = normals[concaveEdgeNormals[edgeNormalI]];
327 const vector& convexNormal
328 = normals[convexEdgeBNormals[edgeBnormalI]];
332 !
areParallel(concaveNormal, convexNormal, tolParallel)
335 convexEdgePlaneDNormal = convexNormal;
337 plane planeD(featPt, convexEdgePlaneDNormal);
341 + 2.0*planeD.distance(internalPtB)
342 *convexEdgePlaneDNormal;
369 concaveEdgeExternalPt,
388 const label concaveEdgeExternalPtIndex(pts.
last().index());
402 vector convexEdgesPlaneNormal =
403 0.5*(convexEdgePlaneCNormal + convexEdgePlaneDNormal);
405 plane planeM(featPt, convexEdgesPlaneNormal);
430 concaveEdgeExternalPt
432 + 2.0*(concaveEdgeLocalFeatPt - concaveEdgeExternalPt);
445 const label internalPtFIndex(pts.
last().index());
449 concaveEdgeExternalPtIndex,
455 + 2.0*planeM.distance(internalPtF)*convexEdgesPlaneNormal;
477 for (label ptI = initialNumOfPoints; ptI < pts.
size(); ++ptI)
479 Info<<
"Point " << ptI <<
" : ";
495 Info<<
"nExternal == 1 && nInternal == 2" <<
endl;
509 label nVert = foamyHexMesh_.number_of_vertices();
511 const label initialNumOfPoints = pts.
size();
519 label convexEdgeI = -1;
529 convexEdgeI = pEds[i];
533 concaveEdgesI[nConcave++] = pEds[i];
538 <<
"Edge " << eS <<
" is flat"
544 <<
"Edge " << eS <<
" not concave/convex"
549 const vector& convexEdgePlaneANormal =
550 normals[edgeNormals[convexEdgeI][0]];
552 const vector& convexEdgePlaneBNormal =
553 normals[edgeNormals[convexEdgeI][1]];
559 featPt - ppDist*convexEdgePlaneANormal,
560 convexEdgePlaneANormal
565 featPt - ppDist*convexEdgePlaneBNormal,
566 convexEdgePlaneBNormal
577 featPt + ppDist*convexEdgeDir;
582 plane planeF(convexEdgeLocalFeatPt, convexEdgeDir);
584 const Foam::point convexEdgeExternalPt = planeF.planePlaneIntersect
592 planeA =
plane(featPt, convexEdgePlaneANormal);
594 planeB =
plane(featPt, convexEdgePlaneBNormal);
598 + 2.0*planeA.distance(convexEdgeExternalPt)
599 *convexEdgePlaneANormal;
612 const label internalPtAIndex(pts.
last().index());
616 + 2.0*planeB.distance(convexEdgeExternalPt)
617 *convexEdgePlaneBNormal;
630 const label internalPtBIndex(pts.
last().index());
640 const labelList& convexEdgeNormals = edgeNormals[convexEdgeI];
641 const labelList& concaveEdgeANormals = edgeNormals[concaveEdgesI[0]];
642 const labelList& concaveEdgeBNormals = edgeNormals[concaveEdgesI[1]];
644 forAll(convexEdgeNormals, edgeNormalI)
646 bool concaveEdgeA =
false;
647 bool concaveEdgeB =
false;
649 forAll(concaveEdgeANormals, edgeAnormalI)
651 const vector& convexNormal
652 = normals[convexEdgeNormals[edgeNormalI]];
653 const vector& concaveNormal
654 = normals[concaveEdgeANormals[edgeAnormalI]];
658 Info<<
"Angle between vectors = "
664 if (
areParallel(convexNormal, concaveNormal, tolParallel))
670 forAll(concaveEdgeBNormals, edgeBnormalI)
672 const vector& convexNormal
673 = normals[convexEdgeNormals[edgeNormalI]];
674 const vector& concaveNormal
675 = normals[concaveEdgeBNormals[edgeBnormalI]];
679 Info<<
"Angle between vectors = "
685 if (
areParallel(convexNormal, concaveNormal, tolParallel))
693 (concaveEdgeA && concaveEdgeB)
694 || (!concaveEdgeA && !concaveEdgeB)
698 <<
"Both or neither of the concave edges share the convex "
700 <<
" concaveEdgeA = " << concaveEdgeA
701 <<
" concaveEdgeB = " << concaveEdgeB
705 for (label i = 0; i < 2; ++i)
716 forAll(concaveEdgeANormals, edgeAnormalI)
718 const vector& convexNormal
719 = normals[convexEdgeNormals[edgeNormalI]];
720 const vector& concaveNormal
721 = normals[concaveEdgeANormals[edgeAnormalI]];
725 !
areParallel(convexNormal, concaveNormal, tolParallel)
728 concaveEdgePlaneCNormal = concaveNormal;
730 plane planeC(featPt, concaveEdgePlaneCNormal);
734 - 2.0*planeC.distance(internalPtA)
735 *concaveEdgePlaneCNormal;
759 forAll(concaveEdgeBNormals, edgeBnormalI)
761 const vector& convexNormal
762 = normals[convexEdgeNormals[edgeNormalI]];
763 const vector& concaveNormal
764 = normals[concaveEdgeBNormals[edgeBnormalI]];
768 !
areParallel(convexNormal, concaveNormal, tolParallel)
771 concaveEdgePlaneDNormal = concaveNormal;
773 plane planeD(featPt, concaveEdgePlaneDNormal);
777 - 2.0*planeD.distance(internalPtB)
778 *concaveEdgePlaneDNormal;
805 convexEdgeExternalPt,
836 vector convexEdgesPlaneNormal =
837 0.5*(concaveEdgePlaneCNormal + concaveEdgePlaneDNormal);
839 plane planeM(featPt, convexEdgesPlaneNormal);
866 + 2.0*(convexEdgeLocalFeatPt - convexEdgeExternalPt);
881 pts[pts.
size() - 2].index(),
887 - 2.0*planeM.distance(internalPtF)*convexEdgesPlaneNormal;
902 pts[pts.
size() - 2].index(),
909 for (label ptI = initialNumOfPoints; ptI < pts.
size(); ++ptI)
911 Info<<
"Point " << ptI <<
" "
An indexed form of CGAL::Triangulation_vertex_base_3<K> used to keep track of the Delaunay vertices i...
label vertexCount() const
Return the vertex count (the next unique vertex index)
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
T remove()
Remove and return the last element. Fatal on an empty list.
void append(const T &val)
Copy append an element to the end of this list.
bool found(const Key &key) const
Return true if hashed entry is found in table.
void size(const label n)
Older name for setAddressableSize.
T & last()
Return the last element of the list.
static bool & parRun() noexcept
Test if this a parallel run.
bool positionOnThisProcessor(const point &pt) const
Is the given position inside the domain of this decomposition.
scalar maxQuadAngle() const
Return the maxQuadAngle.
const pointField & points() const noexcept
Return points.
const labelListList & edgeNormals() const
Return the indices of the normals that are adjacent to the.
@ FLAT
Neither concave or convex, on a flat surface.
const vectorField & normals() const
Return the normals of the surfaces adjacent to the feature edges.
vector edgeDirection(label edgeI, label ptI) const
Return the direction of edgeI, pointing away from ptI.
static const Enum< vertexType > vertexTypeNames_
Geometric class that creates a 3D plane and can return the intersection point between a line and the ...
Hold the types of feature edges attached to the point.
bool addPointPair(const labelPair &vA, const labelPair &vB)
int myProcNo() const noexcept
Return processor number.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
#define WarningInFunction
Report a warning using Foam::Warning.
constexpr scalar pi(M_PI)
pointFromPoint topoint(const Point &P)
messageStream Info
Information stream (stdout output on master, null elsewhere)
vector point
Point is a vector.
Ostream & endl(Ostream &os)
Add newline and flush stream.
constexpr scalar radToDeg() noexcept
Multiplication factor for radians to degrees conversion.
static constexpr const zero Zero
Global zero (0)
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define forAll(list, i)
Loop across all elements in list.