triangle.H 11.4 KB
Newer Older
1 2 3 4
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
OpenFOAM bot's avatar
OpenFOAM bot committed
5
    \\  /    A nd           | www.openfoam.com
OpenFOAM bot's avatar
OpenFOAM bot committed
6 7
     \\/     M anipulation  |
-------------------------------------------------------------------------------
OpenFOAM bot's avatar
OpenFOAM bot committed
8
    Copyright (C) 2011-2017 OpenFOAM Foundation
9
    Copyright (C) 2018-2020 OpenCFD Ltd.
10 11 12 13
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

14 15 16 17
    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.
18 19 20 21 22 23 24

    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
25
    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42

Class
    Foam::triangle

Description
    A triangle primitive used to calculate face normals and swept volumes.

SourceFiles
    triangleI.H

\*---------------------------------------------------------------------------*/

#ifndef triangle_H
#define triangle_H

#include "intersection.H"
#include "vector.H"
43
#include "tensor.H"
44
#include "pointHit.H"
45
#include "Random.H"
46 47
#include "FixedList.H"
#include "UList.H"
mattijs's avatar
mattijs committed
48
#include "linePointRef.H"
49
#include "barycentric2D.H"
mattijs's avatar
mattijs committed
50

51 52 53 54 55
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

namespace Foam
{

56 57
// Forward declarations

58 59
class Istream;
class Ostream;
60 61
class triPoints;
class plane;
62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78

template<class Point, class PointRef> class triangle;

template<class Point, class PointRef>
inline Istream& operator>>
(
    Istream&,
    triangle<Point, PointRef>&
);

template<class Point, class PointRef>
inline Ostream& operator<<
(
    Ostream&,
    const triangle<Point, PointRef>&
);

79 80
typedef triangle<point, const point&> triPointRef;

81 82

/*---------------------------------------------------------------------------*\
83
                          Class triangle Declaration
84 85 86 87 88
\*---------------------------------------------------------------------------*/

template<class Point, class PointRef>
class triangle
{
89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143
public:

    // Public typedefs

        //- Storage type for triangles originating from intersecting triangle
        //  with another triangle
        typedef FixedList<triPoints, 27> triIntersectionList;

        //- Return types for classify
        enum proxType
        {
            NONE,
            POINT,  // Close to point
            EDGE    // Close to edge
        };


    // Public classes

        // Classes for use in sliceWithPlane. What to do with decomposition
        // of triangle.

            //- Dummy
            class dummyOp
            {
            public:
                inline void operator()(const triPoints&);
            };

            //- Sum resulting areas
            class sumAreaOp
            {
            public:
                scalar area_;

                inline sumAreaOp();

                inline void operator()(const triPoints&);
            };

            //- Store resulting tris
            class storeOp
            {
            public:
                triIntersectionList& tris_;
                label& nTris_;

                inline storeOp(triIntersectionList&, label&);

                inline void operator()(const triPoints&);
            };


private:

144 145 146 147 148
    // Private data

        PointRef a_, b_, c_;


149
   // Private Member Functions
150

151 152 153 154 155 156 157 158 159 160 161 162 163
        //- Helper: calculate intersection point
        inline static point planeIntersection
        (
            const FixedList<scalar, 3>& d,
            const triPoints& t,
            const label negI,
            const label posI
        );

        //- Helper: slice triangle with plane
        template<class AboveOp, class BelowOp>
        inline static void triSliceWithPlane
        (
164
            const plane& pln,
165 166 167 168
            const triPoints& tri,
            AboveOp& aboveOp,
            BelowOp& belowOp
        );
169 170


171 172
public:

173 174 175 176 177
    // Constructors

        //- Construct from three points
        inline triangle(const Point& a, const Point& b, const Point& c);

178
        //- Construct from three points
179
        inline triangle(const FixedList<Point, 3>& tri);
180

181 182 183 184
        //- Construct from three points in the list of points
        //  The indices could be from triFace etc.
        inline triangle
        (
185
            const UList<Point>& points,
186 187 188
            const FixedList<label, 3>& indices
        );

189
        //- Construct from Istream
190
        inline triangle(Istream& is);
191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211


    // Member Functions

        // Access

            //- Return first vertex
            inline const Point& a() const;

            //- Return second vertex
            inline const Point& b() const;

            //- Return third vertex
            inline const Point& c() const;


        // Properties

            //- Return centre (centroid)
            inline Point centre() const;

212 213 214 215 216 217
            //- The area normal - with magnitude equal to area of triangle
            inline vector areaNormal() const;

            //- Return unit normal
            inline vector unitNormal() const;

218 219
            //- Legacy name for areaNormal().
            //  \deprecated(2018-06) Deprecated for new use
220 221
            FOAM_DEPRECATED_FOR(2018-12, "areaNormal() or unitNormal()")
            vector normal() const
222 223 224 225
            {
                return areaNormal();
            }

226 227 228 229
            //- Return scalar magnitude
            inline scalar mag() const;

            //- Return circum-centre
230
            inline Point circumCentre() const;
231 232 233 234

            //- Return circum-radius
            inline scalar circumRadius() const;

235 236 237
            //- Return quality: Ratio of triangle and circum-circle
            //  area, scaled so that an equilateral triangle has a
            //  quality of 1
238 239 240 241 242
            inline scalar quality() const;

            //- Return swept-volume
            inline scalar sweptVol(const triangle& t) const;

243 244 245 246
            //- Return the inertia tensor, with optional reference
            //  point and density specification
            inline tensor inertia
            (
Andrew Heather's avatar
Andrew Heather committed
247
                PointRef refPt = Zero,
248 249 250
                scalar density = 1.0
            ) const;

251 252 253 254
            //- Return a random point on the triangle from a uniform
            //  distribution
            inline Point randomPoint(Random& rndGen) const;

255 256 257 258 259 260 261 262 263
            //- Calculate the point from the given barycentric coordinates.
            inline Point barycentricToPoint(const barycentric2D& bary) const;

            //- Calculate the barycentric coordinates from the given point
            inline barycentric2D pointToBarycentric(const point& pt) const;

            //- Calculate the barycentric coordinates from the given point.
            //  Returns the determinant.
            inline scalar pointToBarycentric
264 265
            (
                const point& pt,
266
                barycentric2D& bary
267 268
            ) const;

269
            //- Return point intersection with a ray.
270
            //  For a hit, the distance is signed. Positive number
271 272 273 274 275 276 277 278 279 280 281 282
            //  represents the point in front of triangle.
            //  In case of miss pointHit is set to nearest point
            //  on triangle and its distance to the distance between
            //  the original point and the plane intersection point
            inline pointHit ray
            (
                const point& p,
                const vector& q,
                const intersection::algorithm = intersection::FULL_RAY,
                const intersection::direction dir = intersection::VECTOR
            ) const;

mattijs's avatar
mattijs committed
283 284 285
            //- Fast intersection with a ray.
            //  For a hit, the pointHit.distance() is the line parameter t :
            //  intersection=p+t*q. Only defined for VISIBLE, FULL_RAY or
mattijs's avatar
mattijs committed
286
            //  HALF_RAY. tol increases the virtual size of the triangle
mattijs's avatar
mattijs committed
287
            // by a relative factor.
mattijs's avatar
mattijs committed
288 289 290 291
            inline pointHit intersection
            (
                const point& p,
                const vector& q,
mattijs's avatar
mattijs committed
292 293
                const intersection::algorithm alg,
                const scalar tol = 0.0
mattijs's avatar
mattijs committed
294 295
            ) const;

296 297 298 299 300 301
            //- Find the nearest point to p on the triangle 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 2 is from f[2] to f[0]
            pointHit nearestPointClassify
302
            (
303 304 305
                const point& p,
                label& nearType,
                label& nearLabel
306 307
            ) const;

308 309 310 311 312 313
            //- Return nearest point to p on triangle
            inline pointHit nearestPoint(const point& p) const;

            //- Classify nearest point to p in triangle plane
            //  w.r.t. triangle edges and points.  Returns inside
            //  (true)/outside (false).
314 315 316 317 318 319 320
            bool classify
            (
                const point& p,
                label& nearType,
                label& nearLabel
            ) const;

mattijs's avatar
mattijs committed
321 322 323 324 325 326 327 328 329
            //- Return nearest point to line on triangle. Returns hit if
            //  point is inside triangle. Sets edgePoint to point on edge
            //  (hit if nearest is inside line)
            inline pointHit nearestPoint
            (
                const linePointRef& edge,
                pointHit& edgePoint
            ) const;

330 331 332 333 334 335 336 337
            //- The sign for which side of the face plane the point is on.
            //  Uses the supplied tolerance for rounding around zero.
            //  \return
            //  -  0: on plane
            //  - +1: front-side
            //  - -1: back-side
            inline int sign(const point& p, const scalar tol = SMALL) const;

338 339 340 341
            //- Decompose triangle into triangles above and below plane
            template<class AboveOp, class BelowOp>
            inline void sliceWithPlane
            (
342
                const plane& pln,
343 344 345 346 347 348 349 350 351 352 353 354 355 356 357 358
                AboveOp& aboveOp,
                BelowOp& belowOp
            ) const;

            //- Decompose triangle into triangles inside and outside
            //  (with respect to user provided normal) other
            //  triangle.
            template<class InsideOp, class OutsideOp>
            inline void triangleOverlap
            (
                const vector& n,
                const triangle<Point, PointRef>& tri,
                InsideOp& insideOp,
                OutsideOp& outsideOp
            ) const;

359 360 361 362 363 364 365 366 367 368 369 370 371 372 373 374 375 376 377 378 379 380 381 382 383 384 385

    // IOstream operators

        friend Istream& operator>> <Point, PointRef>
        (
            Istream&,
            triangle&
        );

        friend Ostream& operator<< <Point, PointRef>
        (
            Ostream&,
            const triangle&
        );
};


// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

} // End namespace Foam

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#include "triangleI.H"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

386 387 388 389 390 391
#ifdef NoRepository
#   include "triangle.C"
#endif

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

392 393 394
#endif

// ************************************************************************* //