triFaceI.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-2013 OpenFOAM Foundation
9 Copyright (C) 2017-2022 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
27\*---------------------------------------------------------------------------*/
28
29#include "IOstreams.H"
30#include "face.H"
31#include "triPointRef.H"
32#include <algorithm> // For std::swap
33
34// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
35
36inline int Foam::triFace::compare(const triFace& a, const triFace& b)
37{
38 if
39 (
40 (a[0] == b[0] && a[1] == b[1] && a[2] == b[2])
41 || (a[0] == b[1] && a[1] == b[2] && a[2] == b[0])
42 || (a[0] == b[2] && a[1] == b[0] && a[2] == b[1])
43 )
44 {
45 // identical
46 return 1;
47 }
48 else if
49 (
50 (a[0] == b[2] && a[1] == b[1] && a[2] == b[0])
51 || (a[0] == b[1] && a[1] == b[0] && a[2] == b[2])
52 || (a[0] == b[0] && a[1] == b[2] && a[2] == b[1])
53 )
54 {
55 // same face, but reversed orientation
56 return -1;
57 }
58
59 return 0;
60}
61
62
63// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
64
66:
67 FixedList<label, 3>(-1)
68{}
69
70
72(
73 const label a,
74 const label b,
75 const label c
76)
77{
78 operator[](0) = a;
79 operator[](1) = b;
80 operator[](2) = c;
81}
82
83
84inline Foam::triFace::triFace(std::initializer_list<label> list)
85:
86 FixedList<label, 3>(list)
87{}
88
89
91:
92 FixedList<label, 3>(list)
93{}
94
95
97(
98 const labelUList& list,
99 const FixedList<label, 3>& indices
100)
101:
102 FixedList<label, 3>(list, indices)
103{}
104
105
107:
108 FixedList<label, 3>(is)
109{}
110
111
112// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
113
114inline bool Foam::triFace::valid() const
115{
116 return
117 (
118 operator[](0) >= 0 && operator[](0) != operator[](1)
119 && operator[](1) >= 0 && operator[](1) != operator[](2)
120 && operator[](2) >= 0 && operator[](0) != operator[](2)
121 );
122}
123
124
125inline Foam::label Foam::triFace::collapse()
126{
127 // Cannot resize FixedList, so mark duplicates with '-1'
128 // (the lower vertex is retained)
129 // catch any '-1' (eg, if called multiple times)
130
131 label n = 3;
132 if (operator[](0) == operator[](1) || operator[](1) == -1)
133 {
134 operator[](1) = -1;
135 n--;
136 }
137 else if (operator[](1) == operator[](2) || operator[](2) == -1)
138 {
139 operator[](2) = -1;
140 n--;
141 }
142 if (operator[](0) == operator[](2))
143 {
144 operator[](2) = -1;
145 n--;
146 }
147
148 return n;
149}
150
151
153{
154 std::swap(operator[](1), operator[](2));
155}
156
157
159{
160 pointField p(3);
161
162 p[0] = pts[operator[](0)];
163 p[1] = pts[operator[](1)];
164 p[2] = pts[operator[](2)];
165
166 return p;
167}
168
169
171{
172 return Foam::face(*this);
173}
174
175
177{
178 return triPointRef
179 (
180 points[operator[](0)],
181 points[operator[](1)],
182 points[operator[](2)]
183 );
184}
185
186
188{
189 return (1.0/3.0)*
190 (
191 points[operator[](0)]
192 + points[operator[](1)]
193 + points[operator[](2)]
194 );
195}
196
197
199{
200 // As per triPointRef(...).areaNormal()
201 return 0.5*
202 (
203 (points[operator[](1)] - points[operator[](0)])
204 ^(points[operator[](2)] - points[operator[](0)])
205 );
206}
207
208
210{
211 const vector n(areaNormal(points));
212 const scalar s(Foam::mag(n));
213 return s < ROOTVSMALL ? Zero : n/s;
214}
215
216
217inline Foam::scalar Foam::triFace::mag(const UList<point>& points) const
218{
219 return ::Foam::mag(areaNormal(points));
220}
221
222
223inline Foam::label Foam::triFace::nTriangles() const noexcept
224{
225 return 1;
226}
227
228
230{
231 // The starting points of the original and reverse face are identical.
232 return triFace(operator[](0), operator[](2), operator[](1));
233}
234
235
236inline Foam::label Foam::triFace::which(const label pointLabel) const
237{
238 return FixedList<label, 3>::find(pointLabel);
239}
240
241
242inline Foam::label Foam::triFace::thisLabel(const label i) const
243{
244 return operator[](i);
245}
246
247
248inline Foam::label Foam::triFace::nextLabel(const label i) const
249{
250 return operator[]((i == 2 ? 0 : i+1));
251}
252
253
254inline Foam::label Foam::triFace::prevLabel(const label i) const
255{
256 return operator[]((i ? i-1 : 2));
257}
258
259
260inline Foam::scalar Foam::triFace::sweptVol
261(
262 const UList<point>& opts,
263 const UList<point>& npts
264) const
265{
266 return (1.0/6.0)*
267 (
268 (
269 (npts[operator[](0)] - opts[operator[](0)])
270 & (
271 (opts[operator[](1)] - opts[operator[](0)])
272 ^ (opts[operator[](2)] - opts[operator[](0)])
273 )
274 )
275 + (
276 (npts[operator[](1)] - opts[operator[](1)])
277 & (
278 (opts[operator[](2)] - opts[operator[](1)])
279 ^ (npts[operator[](0)] - opts[operator[](1)])
280 )
281 )
282 + (
283 (opts[operator[](2)] - npts[operator[](2)])
284 & (
285 (npts[operator[](1)] - npts[operator[](2)])
286 ^ (npts[operator[](0)] - npts[operator[](2)])
287 )
288 )
289 );
290}
291
292
294(
295 const UList<point>& points,
296 const point& refPt,
297 scalar density
298) const
299{
300 // a triangle, do a direct calculation
301 return this->tri(points).inertia(refPt, density);
302}
303
304
306(
307 const point& p,
308 const vector& q,
309 const UList<point>& points,
310 const intersection::algorithm alg,
312) const
313{
314 return this->tri(points).ray(p, q, alg, dir);
315}
316
317
318
320(
321 const point& p,
322 const vector& q,
323 const UList<point>& points,
324 const intersection::algorithm alg,
325 const scalar tol
326) const
327{
328 return this->tri(points).intersection(p, q, alg, tol);
329}
330
331
333(
334 const point& p,
335 const vector& q,
336 const point& ctr,
337 const UList<point>& points,
338 const intersection::algorithm alg,
339 const scalar tol
340) const
341{
342 return intersection(p, q, points, alg, tol);
343}
344
345
347(
348 const point& p,
349 const UList<point>& points
350) const
351{
352 return this->tri(points).nearestPoint(p);
353}
354
355
357(
358 const point& p,
359 const UList<point>& points,
360 label& nearType,
361 label& nearLabel
362) const
363{
364 return this->tri(points).nearestPointClassify(p, nearType, nearLabel);
365}
366
367
369(
370 const point& p,
371 const UList<point>& points,
372 const scalar tol
373) const
374{
375 return this->tri(points).sign(p, tol);
376}
377
378
379inline Foam::label Foam::triFace::nEdges() const noexcept
380{
381 return 3;
382}
383
384
385inline Foam::edge Foam::triFace::edge(const label edgei) const
386{
387 return Foam::edge(thisLabel(edgei), nextLabel(edgei));
388}
389
390
392(
393 const label edgei,
394 const UList<point>& pts
395) const
396{
397 return vector(pts[nextLabel(edgei)] - pts[thisLabel(edgei)]);
398}
399
400
401inline Foam::edge Foam::triFace::rcEdge(const label edgei) const
402{
403 // Edge 0 (forward and reverse) always starts at [0]
404 // for consistency with face flipping
405 const label pointi = edgei ? (3 - edgei) : 0;
406 return Foam::edge(thisLabel(pointi), prevLabel(pointi));
407}
408
409
411(
412 const label edgei,
413 const UList<point>& pts
414) const
415{
416 // Edge 0 (forward and reverse) always starts at [0]
417 // for consistency with face flipping
418 const label pointi = edgei ? (3 - edgei) : 0;
419 return vector(pts[prevLabel(pointi)] - pts[thisLabel(pointi)]);
420}
421
422
424{
425 edgeList theEdges(3);
426
427 theEdges[0].first() = operator[](0);
428 theEdges[0].second() = operator[](1);
429
430 theEdges[1].first() = operator[](1);
431 theEdges[1].second() = operator[](2);
432
433 theEdges[2].first() = operator[](2);
434 theEdges[2].second() = operator[](0);
435
436 return theEdges;
437}
438
439
441{
442 edgeList theEdges(3);
443
444 theEdges[0].first() = operator[](0);
445 theEdges[0].second() = operator[](2);
446
447 theEdges[1].first() = operator[](2);
448 theEdges[1].second() = operator[](1);
449
450 theEdges[2].first() = operator[](1);
451 theEdges[2].second() = operator[](0);
452
453 return theEdges;
454}
455
456
458{
459 if (e.first() == operator[](0))
460 {
461 if (e.second() == operator[](1)) return 1; // Forward
462 if (e.second() == operator[](2)) return -1; // Reverse
463 }
464 if (e.first() == operator[](1))
465 {
466 if (e.second() == operator[](2)) return 1; // Forward
467 if (e.second() == operator[](0)) return -1; // Reverse
468 }
469 if (e.first() == operator[](2))
470 {
471 if (e.second() == operator[](0)) return 1; // Forward
472 if (e.second() == operator[](1)) return -1; // Reverse
473 }
474
475 return 0; // Not found
476}
477
478
479// * * * * * * * * * * * * * * * Global Operators * * * * * * * * * * * * * //
480
481inline bool Foam::operator==(const triFace& a, const triFace& b)
482{
483 return triFace::compare(a,b) != 0;
484}
485
486
487inline bool Foam::operator!=(const triFace& a, const triFace& b)
488{
489 return triFace::compare(a,b) == 0;
490}
491
492
493// ************************************************************************* //
Useful combination of include files which define Sin, Sout and Serr and the use of IO streams general...
label n
A 1D vector of objects of type <T> with a fixed length <N>.
Definition: FixedList.H:81
label find(const T &val, label pos=0) const
Find index of the first occurrence of the value.
Definition: FixedList.C:50
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:64
scalar centre() const
Mid-point location, zero for an empty list.
Definition: PDRblockI.H:67
Describes the interaction of a face and a point. It carries the info of a successful hit and (if succ...
Definition: PointHit.H:54
T & first()
Return the first element of the list.
Definition: UListI.H:202
friend complex sign(const complex &c)
sgn() https://en.wikipedia.org/wiki/Sign_function#Complex_signum
Definition: complexI.H:228
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:66
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
Computes the magnitude of an input field.
Definition: mag.H:153
Foam::intersection.
Definition: intersection.H:53
A reference point and direction.
Definition: plane.H:110
ray(const point &pt, const vector &dir)
Definition: plane.H:116
A triangular face using a FixedList of labels corresponding to mesh vertices.
Definition: triFace.H:72
scalar sweptVol(const UList< point > &opts, const UList< point > &npts) const
Return swept-volume from old-points to new-points.
Definition: triFaceI.H:261
void flip()
Flip the face in-place.
Definition: triFaceI.H:152
pointHit nearestPointClassify(const point &p, const UList< point > &points, label &nearType, label &nearLabel) const
Return nearest point to face and classify it:
Definition: triFaceI.H:357
pointHit nearestPoint(const point &p, const UList< point > &points) const
Return nearest point to face.
Definition: triFaceI.H:347
bool valid() const
Return true if the vertices are unique and non-negative.
Definition: triFaceI.H:114
label nEdges() const noexcept
Return number of edges == 3.
Definition: triFaceI.H:379
triFace reverseFace() const
Return face with reverse direction.
Definition: triFaceI.H:229
int edgeDirection(const Foam::edge &e) const
Test the edge direction on the face.
Definition: triFaceI.H:457
Foam::edge rcEdge(const label edgei) const
Return i-th face edge in reverse walk order.
Definition: triFaceI.H:401
triFace()
Default construct, with invalid point labels (-1)
Definition: triFaceI.H:65
tensor inertia(const UList< point > &points, const point &refPt=vector::zero, scalar density=1.0) const
Return the inertia tensor, with optional reference.
Definition: triFaceI.H:294
label which(const label pointLabel) const
Find local index on face for the point label, same as find()
Definition: triFaceI.H:236
static int compare(const triFace &a, const triFace &b)
Compare triFaces.
Definition: triFaceI.H:36
label nextLabel(const label i) const
Next vertex on face.
Definition: triFaceI.H:248
label thisLabel(const label i) const
Definition: triFaceI.H:242
edgeList rcEdges() const
Return list of edges in reverse walk order.
Definition: triFaceI.H:440
label collapse()
'Collapse' face by marking duplicate point labels.
Definition: triFaceI.H:125
face triFaceFace() const
Return triangle as a face.
Definition: triFaceI.H:170
edgeList edges() const
Return list of edges in forward walk order.
Definition: triFaceI.H:423
label nTriangles() const noexcept
Number of triangles after splitting == 1.
Definition: triFaceI.H:223
label prevLabel(const label i) const
Previous vertex on face.
Definition: triFaceI.H:254
triPointRef tri(const UList< point > &points) const
Return the triangle.
Definition: triFaceI.H:176
A triangle primitive used to calculate face normals and swept volumes.
Definition: triangle.H:80
vector areaNormal() const
The area normal - with magnitude equal to area of triangle.
Definition: triangleI.H:112
vector unitNormal() const
Return unit normal.
Definition: triangleI.H:119
volScalarField & p
const pointField & points
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
bool operator!=(const eddy &a, const eddy &b)
Definition: eddy.H:239
tmp< faMatrix< Type > > operator==(const faMatrix< Type > &, const faMatrix< Type > &)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
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
const direction noexcept
Definition: Scalar.H:223
Vector< scalar > vector
Definition: vector.H:61
volScalarField & b
Definition: createFields.H:27
volScalarField & e
Definition: createFields.H:11