/*---------------------------------------------------------------------------*\ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | \\ / A nd | www.openfoam.com \\/ M anipulation | ------------------------------------------------------------------------------- Copyright (C) 2011-2018 OpenFOAM Foundation Copyright (C) 2019 OpenCFD Ltd. ------------------------------------------------------------------------------- License This file is part of OpenFOAM. OpenFOAM is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version. OpenFOAM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>. Class Foam::triSurfaceTools Description A collection of tools for triSurface. Note The curvature calculation is an implementation of the algorithm from: \verbatim "Estimating Curvatures and their Derivatives on Triangle Meshes" by S. Rusinkiewicz 3DPVT'04 Proceedings of the 3D Data Processing, Visualization, and Transmission, 2nd International Symposium Pages 486-493 http://gfx.cs.princeton.edu/pubs/_2004_ECA/curvpaper.pdf \endverbatim SourceFiles triSurfaceTools.C triSurfaceCloseness.C triSurfaceCurvature.C \*---------------------------------------------------------------------------*/ #ifndef triSurfaceTools_H #define triSurfaceTools_H #include "boolList.H" #include "pointField.H" #include "vectorField.H" #include "triadFieldFwd.H" #include "DynamicList.H" #include "HashSet.H" #include "Map.H" #include "FixedList.H" #include "Pair.H" #include "vector2D.H" #include "triPointRef.H" #include "surfaceLocation.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // namespace Foam { // Forward declaration of classes class boundBox; class edge; class labelledTri; class polyBoundaryMesh; class plane; class triSurface; class face; class Time; template<class Face> class MeshedSurface; /*---------------------------------------------------------------------------*\ Class triSurfaceTools Declaration \*---------------------------------------------------------------------------*/ class triSurfaceTools { // Private Member Functions // Refinement enum refineType { NONE, RED, GREEN }; static void calcRefineStatus ( const triSurface& surf, const label facei, List<refineType>& refine ); static void greenRefine ( const triSurface& surf, const label facei, const label edgeI, const label newPointi, DynamicList<labelledTri>& newFaces ); static triSurface doRefine ( const triSurface& surf, const List<refineType>& refineStatus ); // Coarsening static scalar faceCosAngle ( const point& pStart, const point& pEnd, const point& pLeft, const point& pRight ); static void protectNeighbours ( const triSurface& surf, const label vertI, labelList& faceStatus ); //- Faces to collapse because of edge collapse static labelHashSet getCollapsedFaces ( const triSurface& surf, label edgeI ); // Return value of faceUsed for faces using vertI (local numbering). // Used internally. static label vertexUsesFace ( const triSurface& surf, const labelHashSet& faceUsed, const label vertI ); // Get new connections between faces (because of edge collapse) // in form of tables: // - given edge get other edge // - given edge get other face // A face using point v1 on edge will get connected to a face using // point v2 if they share a common vertex // (but not a common edge since then the triangles collapse to // nothing) static void getMergedEdges ( const triSurface& surf, const label edgeI, const labelHashSet& collapsedFaces, Map<label>& edgeToEdge, Map<label>& edgeToFace ); //- Calculates (cos of) angle across edgeI of facei, // taking into account updated addressing (resulting from edge // collapse) static scalar edgeCosAngle ( const triSurface& surf, const label v1, const point& pt, const labelHashSet& collapsedFaces, const Map<label>& edgeToEdge, const Map<label>& edgeToFace, const label facei, const label edgeI ); //- Calculate minimum (cos of) edge angle using addressing from // collapsing // edge to v1 at pt. Returns 1 if v1 is on edge without neighbours // (and hence no edge angle can be defined) static scalar collapseMinCosAngle ( const triSurface& surf, const label v1, const point& pt, const labelHashSet& collapsedFaces, const Map<label>& edgeToEdge, const Map<label>& edgeToFace ); //- Like collapseMinCosAngle but return true for value < minCos bool collapseCreatesFold ( const triSurface& surf, const label v1, const point& pt, const labelHashSet& collapsedFaces, const Map<label>& edgeToEdge, const Map<label>& edgeToFace, const scalar minCos ); ////- Checks if edge collapse creates triangles on top of each //// other //static bool collapseCreatesDuplicates //( // const triSurface& surf, // const label edgeI, // const Map<bool>& collapsedFaces //); // Tracking //- Finds the triangle edge/point cut by the plane between // a point inside/on edge of a triangle and a point outside. // Returns // - location on edge/point and hit() // - or miss() if no intersection found static surfaceLocation cutEdge ( const triSurface& s, const label triI, const label excludeEdgeI, const label excludePointi, const point& triPoint, const plane& cutPlane, const point& toPoint ); //- Checks if current is on the same triangle as the endpoint // and shifts it there. If so updates current and sets a hit. static void snapToEnd ( const triSurface& s, const surfaceLocation& endInfo, surfaceLocation& current ); //- Visits faces eFaces around start. Does not visit triangle // start.triangle() nor edge excludeEdgeI. // Returns edge, triangle (if more than one choice) which gets // us nearer endpoint. // Returns // - hit() if triangle contains endpoint // - triangle()=-1 if no triangle found // - nearest triangle/edge otherwise static surfaceLocation visitFaces ( const triSurface& s, const labelList& eFaces, const surfaceLocation& start, const label excludeEdgeI, const label excludePointi, const surfaceLocation& end, const plane& cutPlane ); public: // OBJ writing //- Write pointField to OBJ format file static void writeOBJ ( const fileName& fName, const pointField& pts ); //- Write vertex subset to OBJ format file static void writeOBJ ( const triSurface& surf, const fileName& fName, const boolList& markedVerts ); // Additional addressing //- Get all triangles using edge endpoint static void getVertexTriangles ( const triSurface& surf, const label edgeI, labelList& edgeTris ); //- Get all vertices (local numbering) connected to vertices of edge static labelList getVertexVertices ( const triSurface& surf, const edge& e ); //- Get face connected to edge not facei static label otherFace ( const triSurface& surf, const label facei, const label edgeI ); //- Get the two edges on facei counterclockwise after edgeI static void otherEdges ( const triSurface& surf, const label facei, const label edgeI, label& e1, label& e2 ); //- Get the two vertices (local numbering) on facei counterclockwise // vertI static void otherVertices ( const triSurface& surf, const label facei, const label vertI, label& vert1I, label& vert2I ); //- Get edge opposite vertex (local numbering) static label oppositeEdge ( const triSurface& surf, const label facei, const label vertI ); //- Get vertex (local numbering) opposite edge static label oppositeVertex ( const triSurface& surf, const label facei, const label edgeI ); //- Returns edge label connecting v1, v2 (local numbering) static label getEdge ( const triSurface& surf, const label vert1I, const label vert2I ); //- Return index of triangle (or -1) using all three edges static label getTriangle ( const triSurface& surf, const label e0I, const label e1I, const label e2I ); // Coarsening //- Create new triSurface by collapsing edges to edge mids. static triSurface collapseEdges ( const triSurface& surf, const labelList& collapsableEdges ); //- Face collapse status. // anyEdge: any edge can be collapsed // noEdge: no edge can be collapsed // collapsed: already collapsed // >0: edge label that can be collapsed static const label ANYEDGE; static const label NOEDGE; static const label COLLAPSED; //- Create new triSurface by collapsing edges to specified // positions. faceStatus allows // explicit control over which faces need to be protected (see above). // faceStatus gets updated to protect collapsing already collapsed // faces. static triSurface collapseEdges ( const triSurface& surf, const labelList& collapsableEdges, const pointField& edgeMids, labelList& faceStatus ); // Refinement //- Refine edges by splitting to opposite vertex static triSurface greenRefine ( const triSurface& surf, const labelList& refineEdges ); //- Refine face by splitting all edges. Neighbouring face is // greenRefine'd. static triSurface redGreenRefine ( const triSurface& surf, const labelList& refineFaces ); // Geometric //- Returns element in edgeIndices with minimum length static label minEdge ( const triSurface& surf, const labelList& edgeIndices ); //- Returns element in edgeIndices with minimum length static label maxEdge ( const triSurface& surf, const labelList& edgeIndices ); //- Merge points within distance static triSurface mergePoints ( const triSurface& surf, const scalar mergeTol ); //- Triangle (unit) normal. If nearest point to triangle on edge use // edge normal (calculated on the fly); if on vertex use vertex normal. // Uses planarTol. static vector surfaceNormal ( const triSurface& surf, const label nearestFacei, const point& nearestPt ); //- On which side of surface enum sideType { UNKNOWN, // cannot be determined (e.g. non-manifold) INSIDE, // inside OUTSIDE // outside }; //- If nearest point is on edgeI, determine on which side of surface // sample is. static sideType edgeSide ( const triSurface& surf, const point& sample, const point& nearestPoint, const label edgeI ); //- Given nearest point (to sample) on surface determines which side // sample is. Uses either face normal, edge normal or point normal // (non-trivial). Uses triangle::classify. static sideType surfaceSide ( const triSurface& surf, const point& sample, const label nearestFacei ); // Triangulation of faces //- Simple triangulation of (selected patches of) boundaryMesh. Needs // polyMesh (or polyBoundaryMesh) since only at this level are the // triangles on neighbouring patches connected. // Return faceMap from triI to faceI static triSurface triangulate ( const polyBoundaryMesh& mBesh, const labelHashSet& includePatches, labelList& faceMap, const bool verbose = false ); static triSurface triangulate ( const polyBoundaryMesh& bMesh, const labelHashSet& includePatches, const boundBox& bBox, const bool verbose = false ); //- Face-centre triangulation of (selected patches of) boundaryMesh. // Needs // polyMesh (or polyBoundaryMesh) since only at this level are the // triangles on neighbouring patches connected. triSurface triangulateFaceCentre ( const polyBoundaryMesh& mBesh, const labelHashSet& includePatches, const bool verbose = false ); // Triangulation and interpolation //- Do unconstrained Delaunay of points. Returns triSurface with 3D // points with z=0. All triangles in region 0. static triSurface delaunay2D(const List<vector2D>&); //- Calculate linear interpolation weights for point (guaranteed to be // inside triangle) static void calcInterpolationWeights ( const triPointRef& tri, const point& p, FixedList<scalar, 3>& weights ); // Calculate weighting factors from samplePts to triangle it is in. // Uses linear search to find triangle. // Vertices are: // (a b c) : vertices of the triangle abc the point is in // or if the point is outside all triangles: // (a b -1) : the edge ab the point is nearest to. // (a -1 -1) : the vertex a the point is nearest to static void calcInterpolationWeights ( const triSurface& s, const pointField& samplePts, List<FixedList<label, 3>>& allVerts, List<FixedList<scalar, 3>>& allWeights ); // Curvature //- Weighting for normals of faces attached to vertex static scalar vertexNormalWeight ( const triFace& f, const label pI, const vector& fN, const UList<point>& points ); //- Weighted average of normals of attached faces static tmp<vectorField> vertexNormals(const triSurface& surf); //- Local coordinate-system for each point normal static tmp<triadField> vertexTriads ( const triSurface& surf, const vectorField& pointNormals ); //- Surface curvatures at the vertex points static tmp<scalarField> curvatures ( const triSurface& surf, const vectorField& pointNormals, const triadField& pointTriads ); //- Surface curvatures at the vertex points static tmp<scalarField> curvatures ( const triSurface& surf ); //- Write surface curvature at the vertex points and return the field static tmp<scalarField> writeCurvature ( const Time& runTime, const word& basename, const triSurface& surf ); // Closeness //- Check and write internal/external closeness fields static Pair<tmp<scalarField>> writeCloseness ( const Time& runTime, const word& basename, const triSurface& surf, const scalar internalAngleTolerance = 45, const scalar externalAngleTolerance = 10 ); // Feature Proximity //- Calculate feature proximity scalarField featureProximity ( const triSurface& surf, const scalar searchDistance ); //- Check and write internal/external closeness fields static void writeFeatureProximity ( const Time& runTime, const word& basename, const triSurface& surf, const bool writeVTK, const scalar searchDistance ); // Surface checking functionality //- Check single triangle for (topological) validity static bool validTri ( const triSurface&, const label facei, const bool verbose = true ); //- Check single triangle for (topological) validity static bool validTri ( const MeshedSurface<face>&, const label facei, const bool verbose = true ); // Tracking //- Test point on plane of triangle to see if on edge or point or inside static surfaceLocation classify ( const triSurface&, const label triI, const point& trianglePoint ); //- Track on surface to get closer to point. // Possible situations: // - 1. reached endpoint // - 2. reached edge (normal situation) // - 3. reached end of surface (edge on single face) // Input: // - starting position+triangle/edge/point (so has to be on surface!) // - (optional) previous position+triangle/edge/point to prevent // going back. Set index (of triangle/edge/point) to -1 if not // used. // - end position+triangle/edge/point (so has to be on surface!) // - plane to follow. Has to go through end point! // Returns: // - true if end point reached (situation 1) // - new position+triangle/edge/point // Caller has to check for situation 3 by checking that triangle() // is not set. static surfaceLocation trackToEdge ( const triSurface&, const surfaceLocation& start, const surfaceLocation& end, const plane& cutPlane ); //- Track from edge to edge across surface. Uses trackToEdge. // Not really useful by itself, more example of how to use trackToEdge. // endInfo should be location on surface. // hitInfo should be initialised to starting location (on surface as // well). Upon return is set to end location. static void track ( const triSurface&, const surfaceLocation& endInfo, const plane& cutPlane, surfaceLocation& hitInfo ); }; // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // } // End namespace Foam // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // #endif // ************************************************************************* //