diff --git a/src/OpenFOAM/meshes/meshShapes/face/face.H b/src/OpenFOAM/meshes/meshShapes/face/face.H index 81ce90d8a9d537d2021ab3393bfb50235094f93c..dd5fc4a832b4265e0cd9aa0dd5842c224c13ba22 100644 --- a/src/OpenFOAM/meshes/meshShapes/face/face.H +++ b/src/OpenFOAM/meshes/meshShapes/face/face.H @@ -126,6 +126,14 @@ class face public: + //- Return types for classify + enum proxType + { + NONE, + POINT, // Close to point + EDGE // Close to edge + }; + // Static data members static const char* const typeName; @@ -249,6 +257,20 @@ public: const pointField& meshPoints ) const; + //- Return nearest point to face and classify it: + // + near point (nearType=POINT, nearLabel=0, 1, 2) + // + near edge (nearType=EDGE, nearLabel=0, 1, 2) + // Note: edges are counted from starting vertex so + // e.g. edge n is from f[n] to f[0], where the face has n + 1 + // points + pointHit nearestPointClassify + ( + const point& p, + const pointField& meshPoints, + label& nearType, + label& nearLabel + ) const; + //- Return contact sphere diameter scalar contactSphereDiameter ( diff --git a/src/OpenFOAM/meshes/meshShapes/face/faceIntersection.C b/src/OpenFOAM/meshes/meshShapes/face/faceIntersection.C index d4184e766aa6a3dd225942da99c4d1e964330e83..f9f02eb61293307a29975ff878f75ffaef138ec0 100644 --- a/src/OpenFOAM/meshes/meshShapes/face/faceIntersection.C +++ b/src/OpenFOAM/meshes/meshShapes/face/faceIntersection.C @@ -181,6 +181,22 @@ Foam::pointHit Foam::face::nearestPoint const point& p, const pointField& meshPoints ) const +{ + // Dummy labels + label nearType = -1; + label nearLabel = -1; + + return nearestPointClassify(p, meshPoints, nearType, nearLabel); +} + + +Foam::pointHit Foam::face::nearestPointClassify +( + const point& p, + const pointField& meshPoints, + label& nearType, + label& nearLabel +) const { const face& f = *this; point ctr = centre(meshPoints); @@ -188,6 +204,9 @@ Foam::pointHit Foam::face::nearestPoint // Initialize to miss, distance=GREAT pointHit nearest(p); + nearType = -1; + nearLabel = -1; + label nPoints = f.size(); point nextPoint = ctr; @@ -196,8 +215,10 @@ Foam::pointHit Foam::face::nearestPoint { nextPoint = meshPoints[f[fcIndex(pI)]]; + label tmpNearType = -1; + label tmpNearLabel = -1; + // Note: for best accuracy, centre point always comes last - // triPointRef tri ( meshPoints[f[pI]], @@ -205,12 +226,42 @@ Foam::pointHit Foam::face::nearestPoint ctr ); - pointHit curHit = tri.nearestPoint(p); + pointHit curHit = tri.nearestPointClassify + ( + p, + tmpNearType, + tmpNearLabel + ); if (Foam::mag(curHit.distance()) < Foam::mag(nearest.distance())) { nearest.setDistance(curHit.distance()); + // Assume at first that the near type is NONE on the + // triangle (i.e. on the face of the triangle) then it is + // therefore also for the face. + + nearType = NONE; + + if (tmpNearType == triPointRef::EDGE && tmpNearLabel == 0) + { + // If the triangle edge label is 0, then this is also + // an edge of the face, if not, it is on the face + + nearType = EDGE; + + nearLabel = pI; + } + else if (tmpNearType == triPointRef::POINT && tmpNearLabel < 2) + { + // If the triangle point label is 0 or 1, then this is + // also a point of the face, if not, it is on the face + + nearType = POINT; + + nearLabel = pI + tmpNearLabel; + } + if (curHit.hit()) { nearest.setHit(); diff --git a/src/meshTools/indexedOctree/treeDataTriSurface.C b/src/meshTools/indexedOctree/treeDataTriSurface.C index b019a85304642b522af2c4df75b2c55f0cad2036..f380167174481939ae61a056bae233485239d815 100644 --- a/src/meshTools/indexedOctree/treeDataTriSurface.C +++ b/src/meshTools/indexedOctree/treeDataTriSurface.C @@ -360,8 +360,8 @@ void Foam::treeDataTriSurface::findNearest // scalar b(E0& E1); // scalar c(E1& E1); - // // Get nearest point in s,t coordinates (s is along E0, t is along - // // E1) + // // Get nearest point in s,t coordinates (s is along E0, t + // // is along E1) // scalar s; // scalar t; diff --git a/src/meshTools/triSurface/triSurfaceTools/triSurfaceTools.C b/src/meshTools/triSurface/triSurfaceTools/triSurfaceTools.C index c2a88bf7e25598c7fc99c9e756f9a426f4d5e47b..8ade80cbf8127aeeb7f4ccf79330a8bc95ccbbec 100644 --- a/src/meshTools/triSurface/triSurfaceTools/triSurfaceTools.C +++ b/src/meshTools/triSurface/triSurfaceTools/triSurfaceTools.C @@ -2220,40 +2220,40 @@ Foam::triSurfaceTools::sideType Foam::triSurfaceTools::surfaceSide if (nearType == triPointRef::NONE) { + vector sampleNearestVec = (sample - nearestPoint); + // Nearest to face interior. Use faceNormal to determine side -<<<<<<< HEAD - scalar c = (sample - nearestPoint) & surf.faceNormals()[nearestFaceI]; -======= scalar c = sampleNearestVec & surf.faceNormals()[nearestFaceI]; - // If the sample is essentially on the face, do not check for - // it being perpendicular. - - if (magSampleNearestVec > SMALL) - { - c /= magSampleNearestVec*mag(surf.faceNormals()[nearestFaceI]); - - if (mag(c) < 0.99) - { - FatalErrorIn("triSurfaceTools::surfaceSide") - << "nearestPoint identified as being on triangle face " - << "but vector from nearestPoint to sample is not " - << "perpendicular to the normal." << nl - << "sample: " << sample << nl - << "nearestPoint: " << nearestPoint << nl - << "sample - nearestPoint: " << sample - nearestPoint << nl - << "normal: " << surf.faceNormals()[nearestFaceI] << nl - << "mag(sample - nearestPoint): " - << mag(sample - nearestPoint) << nl - << "normalised dot product: " << c << nl - << "triangle vertices: " << nl - << " " << points[f[0]] << nl - << " " << points[f[1]] << nl - << " " << points[f[2]] << nl - << abort(FatalError); - } - } ->>>>>>> 0bb6ebd... ENH: Making nearestPointClassify query for triangle. + // // If the sample is essentially on the face, do not check for + // // it being perpendicular. + + // scalar magSampleNearestVec = mag(sampleNearestVec); + + // if (magSampleNearestVec > SMALL) + // { + // c /= magSampleNearestVec*mag(surf.faceNormals()[nearestFaceI]); + + // if (mag(c) < 0.99) + // { + // FatalErrorIn("triSurfaceTools::surfaceSide") + // << "nearestPoint identified as being on triangle face " + // << "but vector from nearestPoint to sample is not " + // << "perpendicular to the normal." << nl + // << "sample: " << sample << nl + // << "nearestPoint: " << nearestPoint << nl + // << "sample - nearestPoint: " << sample - nearestPoint << nl + // << "normal: " << surf.faceNormals()[nearestFaceI] << nl + // << "mag(sample - nearestPoint): " + // << mag(sample - nearestPoint) << nl + // << "normalised dot product: " << c << nl + // << "triangle vertices: " << nl + // << " " << points[f[0]] << nl + // << " " << points[f[1]] << nl + // << " " << points[f[2]] << nl + // << abort(FatalError); + // } + // } if (c > 0) { @@ -2273,26 +2273,26 @@ Foam::triSurfaceTools::sideType Foam::triSurfaceTools::surfaceSide label edgeI = surf.faceEdges()[nearestFaceI][nearLabel]; // if (debug) - { - // Check order of faceEdges same as face vertices. - const edge& e = surf.edges()[edgeI]; - const labelList& meshPoints = surf.meshPoints(); - const edge meshEdge(meshPoints[e[0]], meshPoints[e[1]]); - - if - ( - meshEdge - != edge(f[nearLabel], f[f.fcIndex(nearLabel)]) - ) - { - FatalErrorIn("triSurfaceTools::surfaceSide") - << "Edge:" << edgeI << " local vertices:" << e - << " mesh vertices:" << meshEdge - << " not at position " << nearLabel - << " in face " << f - << abort(FatalError); - } - } + // { + // // Check order of faceEdges same as face vertices. + // const edge& e = surf.edges()[edgeI]; + // const labelList& meshPoints = surf.meshPoints(); + // const edge meshEdge(meshPoints[e[0]], meshPoints[e[1]]); + + // if + // ( + // meshEdge + // != edge(f[nearLabel], f[f.fcIndex(nearLabel)]) + // ) + // { + // FatalErrorIn("triSurfaceTools::surfaceSide") + // << "Edge:" << edgeI << " local vertices:" << e + // << " mesh vertices:" << meshEdge + // << " not at position " << nearLabel + // << " in face " << f + // << abort(FatalError); + // } + // } return edgeSide(surf, sample, nearestPoint, edgeI); }