triangle.H
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2017 OpenFOAM Foundation
9 Copyright (C) 2018-2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 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
27Class
28 Foam::triangle
29
30Description
31 A triangle primitive used to calculate face normals and swept volumes.
32
33SourceFiles
34 triangleI.H
35
36\*---------------------------------------------------------------------------*/
37
38#ifndef triangle_H
39#define triangle_H
40
41#include "intersection.H"
42#include "vector.H"
43#include "tensor.H"
44#include "pointHit.H"
45#include "Random.H"
46#include "FixedList.H"
47#include "UList.H"
48#include "linePointRef.H"
49#include "barycentric2D.H"
50
51// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
52
53namespace Foam
54{
55
56// Forward Declarations
57class plane;
58class triPoints;
59
60template<class Point, class PointRef> class triangle;
61
62template<class Point, class PointRef>
64
65template<class Point, class PointRef>
67
68
69// Common Typedefs
70
71//- A triangle using referred points
73
74
75/*---------------------------------------------------------------------------*\
76 Class triangle Declaration
77\*---------------------------------------------------------------------------*/
78
79template<class Point, class PointRef>
80class triangle
81{
82public:
83
84 // Public Typedefs
85
86 //- The point type
87 typedef Point point_type;
88
89 //- Storage type for triangles originating from intersecting triangle
90 //- with another triangle
92
93 //- The proximity classifications
94 enum proxType
95 {
96 NONE = 0,
97 POINT,
99 };
100
101
102 // Public classes
103
104 // Classes for use in sliceWithPlane. What to do with decomposition
105 // of triangle.
106
107 //- Dummy
108 class dummyOp
109 {
110 public:
111 inline void operator()(const triPoints&);
112 };
113
114 //- Sum resulting areas
115 class sumAreaOp
116 {
117 public:
118 scalar area_;
119
120 inline sumAreaOp();
121
122 inline void operator()(const triPoints&);
123 };
124
125 //- Store resulting tris
126 class storeOp
127 {
128 public:
130 label& nTris_;
131
132 inline storeOp(triIntersectionList&, label&);
133
134 inline void operator()(const triPoints&);
135 };
136
137
138private:
139
140 // Private Data
141
142 PointRef a_, b_, c_;
143
144
145 // Private Member Functions
146
147 //- Helper: calculate intersection point
148 inline static point planeIntersection
149 (
150 const FixedList<scalar, 3>& d,
151 const triPoints& t,
152 const label negI,
153 const label posI
154 );
155
156 //- Helper: slice triangle with plane
157 template<class AboveOp, class BelowOp>
158 inline static void triSliceWithPlane
159 (
160 const plane& pln,
161 const triPoints& tri,
162 AboveOp& aboveOp,
163 BelowOp& belowOp
164 );
165
166
167public:
168
169 // Constructors
170
171 //- Construct from three points
172 inline triangle(const Point& a, const Point& b, const Point& c);
173
174 //- Construct from three points
175 inline triangle(const FixedList<Point, 3>& tri);
176
177 //- Construct from three points in the list of points
178 // The indices could be from triFace etc.
179 inline triangle
180 (
181 const UList<Point>& points,
182 const FixedList<label, 3>& indices
183 );
184
185 //- Construct from Istream
186 inline explicit triangle(Istream& is);
187
188
189 // Member Functions
190
191 // Access
192
193 //- Return first vertex
194 inline const Point& a() const;
195
196 //- Return second vertex
197 inline const Point& b() const;
198
199 //- Return third vertex
200 inline const Point& c() const;
201
202
203 // Properties
204
205 //- Return centre (centroid)
206 inline Point centre() const;
207
208 //- The area normal - with magnitude equal to area of triangle
209 inline vector areaNormal() const;
210
211 //- Return unit normal
212 inline vector unitNormal() const;
213
214 //- Legacy name for areaNormal().
215 // \deprecated(2018-06) Deprecated for new use
216 FOAM_DEPRECATED_FOR(2018-12, "areaNormal() or unitNormal()")
217 vector normal() const
218 {
219 return areaNormal();
220 }
221
222 //- Return scalar magnitude
223 inline scalar mag() const;
224
225 //- Return circum-centre
226 inline Point circumCentre() const;
227
228 //- Return circum-radius
229 inline scalar circumRadius() const;
230
231 //- Return quality: Ratio of triangle and circum-circle
232 // area, scaled so that an equilateral triangle has a
233 // quality of 1
234 inline scalar quality() const;
235
236 //- Return swept-volume
237 inline scalar sweptVol(const triangle& t) const;
238
239 //- Return the inertia tensor, with optional reference
240 // point and density specification
241 inline tensor inertia
242 (
243 PointRef refPt = Zero,
244 scalar density = 1.0
245 ) const;
246
247 //- Return a random point on the triangle from a uniform
248 // distribution
249 inline Point randomPoint(Random& rndGen) const;
250
251 //- Calculate the point from the given barycentric coordinates.
252 inline Point barycentricToPoint(const barycentric2D& bary) const;
253
254 //- Calculate the barycentric coordinates from the given point
255 inline barycentric2D pointToBarycentric(const point& pt) const;
256
257 //- Calculate the barycentric coordinates from the given point.
258 // Returns the determinant.
259 inline scalar pointToBarycentric
260 (
261 const point& pt,
262 barycentric2D& bary
263 ) const;
264
265 //- Return point intersection with a ray.
266 // For a hit, the distance is signed. Positive number
267 // represents the point in front of triangle.
268 // In case of miss pointHit is set to nearest point
269 // on triangle and its distance to the distance between
270 // the original point and the plane intersection point
271 inline pointHit ray
272 (
273 const point& p,
274 const vector& q,
277 ) const;
278
279 //- Fast intersection with a ray.
280 // For a hit, the pointHit.distance() is the line parameter t :
281 // intersection=p+t*q. Only defined for VISIBLE, FULL_RAY or
282 // HALF_RAY. tol increases the virtual size of the triangle
283 // by a relative factor.
285 (
286 const point& p,
287 const vector& q,
288 const intersection::algorithm alg,
289 const scalar tol = 0.0
290 ) const;
291
292 //- Find the nearest point to p on the triangle and classify it:
293 // + near point (nearType=POINT, nearLabel=0, 1, 2)
294 // + near edge (nearType=EDGE, nearLabel=0, 1, 2)
295 // Note: edges are counted from starting
296 // vertex so e.g. edge 2 is from f[2] to f[0]
298 (
299 const point& p,
300 label& nearType,
301 label& nearLabel
302 ) const;
303
304 //- Return nearest point to p on triangle
305 inline pointHit nearestPoint(const point& p) const;
306
307 //- Classify nearest point to p in triangle plane
308 // w.r.t. triangle edges and points. Returns inside
309 // (true)/outside (false).
310 bool classify
311 (
312 const point& p,
313 label& nearType,
314 label& nearLabel
315 ) const;
316
317 //- Return nearest point to line on triangle. Returns hit if
318 // point is inside triangle. Sets edgePoint to point on edge
319 // (hit if nearest is inside line)
321 (
322 const linePointRef& edge,
323 pointHit& edgePoint
324 ) const;
325
326 //- The sign for which side of the face plane the point is on.
327 // Uses the supplied tolerance for rounding around zero.
328 // \return
329 // - 0: on plane
330 // - +1: front-side
331 // - -1: back-side
332 inline int sign(const point& p, const scalar tol = SMALL) const;
333
334 //- Decompose triangle into triangles above and below plane
335 template<class AboveOp, class BelowOp>
336 inline void sliceWithPlane
337 (
338 const plane& pln,
339 AboveOp& aboveOp,
340 BelowOp& belowOp
341 ) const;
342
343 //- Decompose triangle into triangles inside and outside
344 // (with respect to user provided normal) other
345 // triangle.
346 template<class InsideOp, class OutsideOp>
347 inline void triangleOverlap
348 (
349 const vector& n,
350 const triangle<Point, PointRef>& tri,
351 InsideOp& insideOp,
352 OutsideOp& outsideOp
353 ) const;
354
355
356 // IOstream Operators
358 friend Istream& operator>> <Point, PointRef>
359 (
360 Istream&,
361 triangle&
362 );
364 friend Ostream& operator<< <Point, PointRef>
365 (
366 Ostream&,
367 const triangle&
368 );
369};
370
371
372// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
373
374} // End namespace Foam
375
376// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
377
378#include "triangleI.H"
379
380// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
381
382#ifdef NoRepository
383# include "triangle.C"
384#endif
385
386// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
387
388#endif
389
390// ************************************************************************* //
CGAL::Point_3< K > Point
label n
Templated 2D Barycentric derived from VectorSpace. Has 3 components, one of which is redundant.
Definition: Barycentric2D.H:57
A 1D vector of objects of type <T> with a fixed length <N>.
Definition: FixedList.H:81
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:64
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
Describes the interaction of a face and a point. It carries the info of a successful hit and (if succ...
Definition: PointHit.H:54
Random number generator.
Definition: Random.H:60
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: UList.H:94
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:66
Foam::intersection.
Definition: intersection.H:53
A line primitive.
Definition: line.H:68
Geometric class that creates a 3D plane and can return the intersection point between a line and the ...
Definition: plane.H:95
Triangle storage. Null constructable (unfortunately triangle<point, point> is not)
Definition: triPoints.H:55
void operator()(const triPoints &)
Definition: triangleI.H:836
Store resulting tris.
Definition: triangle.H:126
triIntersectionList & tris_
Definition: triangle.H:128
void operator()(const triPoints &)
Definition: triangleI.H:873
Sum resulting areas.
Definition: triangle.H:115
void operator()(const triPoints &)
Definition: triangleI.H:851
A triangle primitive used to calculate face normals and swept volumes.
Definition: triangle.H:80
Point circumCentre() const
Return circum-centre.
Definition: triangleI.H:135
barycentric2D pointToBarycentric(const point &pt) const
Calculate the barycentric coordinates from the given point.
Definition: triangleI.H:272
void sliceWithPlane(const plane &pln, AboveOp &aboveOp, BelowOp &belowOp) const
Decompose triangle into triangles above and below plane.
Definition: triangle.C:115
scalar circumRadius() const
Return circum-radius.
Definition: triangleI.H:162
Point barycentricToPoint(const barycentric2D &bary) const
Calculate the point from the given barycentric coordinates.
Definition: triangleI.H:262
scalar sweptVol(const triangle &t) const
Return swept-volume.
Definition: triangleI.H:199
vector areaNormal() const
The area normal - with magnitude equal to area of triangle.
Definition: triangleI.H:112
const Point & a() const
Return first vertex.
Definition: triangleI.H:86
Point centre() const
Return centre (centroid)
Definition: triangleI.H:105
FixedList< triPoints, 27 > triIntersectionList
Definition: triangle.H:90
void triangleOverlap(const vector &n, const triangle< Point, PointRef > &tri, InsideOp &insideOp, OutsideOp &outsideOp) const
Decompose triangle into triangles inside and outside.
Definition: triangle.C:128
pointHit nearestPointClassify(const point &p, label &nearType, label &nearLabel) const
Find the nearest point to p on the triangle and classify it:
Definition: triangleI.H:522
proxType
The proximity classifications.
Definition: triangle.H:94
@ POINT
Close to point.
Definition: triangle.H:96
@ EDGE
Close to edge.
Definition: triangle.H:97
@ NONE
Unknown proximity.
Definition: triangle.H:95
pointHit ray(const point &p, const vector &q, const intersection::algorithm=intersection::FULL_RAY, const intersection::direction dir=intersection::VECTOR) const
Return point intersection with a ray.
Definition: triangleI.H:323
vector unitNormal() const
Return unit normal.
Definition: triangleI.H:119
const Point & c() const
Return third vertex.
Definition: triangleI.H:98
tensor inertia(PointRef refPt=Zero, scalar density=1.0) const
Return the inertia tensor, with optional reference.
Definition: triangleI.H:218
scalar mag() const
Return scalar magnitude.
Definition: triangleI.H:128
Point point_type
The point type.
Definition: triangle.H:86
pointHit nearestPoint(const point &p) const
Return nearest point to p on triangle.
Definition: triangleI.H:671
Point randomPoint(Random &rndGen) const
Return a random point on the triangle from a uniform.
Definition: triangleI.H:254
bool classify(const point &p, label &nearType, label &nearLabel) const
Classify nearest point to p in triangle plane.
Definition: triangleI.H:685
int sign(const point &p, const scalar tol=SMALL) const
The sign for which side of the face plane the point is on.
Definition: triangleI.H:823
const Point & b() const
Return second vertex.
Definition: triangleI.H:92
scalar quality() const
Return quality: Ratio of triangle and circum-circle.
Definition: triangleI.H:183
volScalarField & p
const pointField & points
Namespace for OpenFOAM.
Ostream & operator<<(Ostream &, const boundaryPatch &p)
Write boundaryPatch as dictionary entries (without surrounding braces)
Definition: boundaryPatch.C:83
Istream & operator>>(Istream &, directionInfo &)
triangle< point, const point & > triPointRef
A triangle using referred points.
Definition: triangle.H:71
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
#define FOAM_DEPRECATED_FOR(since, replacement)
Definition: stdFoam.H:52
Random rndGen
Definition: createFields.H:23