tetrahedron.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 Class
25  Foam::tetrahedron
26 
27 Description
28  A tetrahedron primitive.
29 
30  Ordering of edges needs to be the same for a tetrahedron
31  class, a tetrahedron cell shape model and a tetCell.
32 
33 SourceFiles
34  tetrahedronI.H
35  tetrahedron.C
36 
37 \*---------------------------------------------------------------------------*/
38 
39 #ifndef tetrahedron_H
40 #define tetrahedron_H
41 
42 #include "point.H"
43 #include "primitiveFieldsFwd.H"
44 #include "pointHit.H"
45 #include "Random.H"
46 #include "FixedList.H"
47 #include "UList.H"
48 #include "triPointRef.H"
49 #include "boundBox.H"
50 #include "barycentric.H"
51 
52 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
53 
54 namespace Foam
55 {
56 
57 class Istream;
58 class Ostream;
59 
60 // Forward declaration of friend functions and operators
61 
62 template<class Point, class PointRef> class tetrahedron;
63 
64 template<class Point, class PointRef>
65 inline Istream& operator>>
66 (
67  Istream&,
69 );
70 
71 template<class Point, class PointRef>
72 inline Ostream& operator<<
73 (
74  Ostream&,
76 );
77 
78 /*---------------------------------------------------------------------------*\
79  class tetrahedron Declaration
80 \*---------------------------------------------------------------------------*/
81 
82 template<class Point, class PointRef>
83 class tetrahedron
84 {
85  // Private Data
86 
87  PointRef a_, b_, c_, d_;
88 
89 
90 public:
91 
92  // Member constants
93 
94  enum
95  {
96  nVertices = 4, // Number of vertices in tetrahedron
97  nEdges = 6 // Number of edges in tetrahedron
98  };
99 
100 
101  // Constructors
102 
103  //- Construct from points
104  inline tetrahedron
105  (
106  const Point& a,
107  const Point& b,
108  const Point& c,
109  const Point& d
110  );
111 
112  //- Construct from four points in the list of points
113  inline tetrahedron
114  (
115  const UList<Point>&,
116  const FixedList<label, 4>& indices
117  );
118 
119  //- Construct from Istream
120  inline tetrahedron(Istream&);
121 
122 
123  // Member Functions
124 
125  // Access
126 
127  //- Return vertices
128  inline const Point& a() const;
129 
130  inline const Point& b() const;
131 
132  inline const Point& c() const;
133 
134  inline const Point& d() const;
135 
136  //- Return i-th face
137  inline triPointRef tri(const label facei) const;
138 
139  // Properties
140 
141  //- Return face normal
142  inline vector Sa() const;
143 
144  inline vector Sb() const;
145 
146  inline vector Sc() const;
147 
148  inline vector Sd() const;
149 
150  //- Return centre (centroid)
151  inline Point centre() const;
152 
153  //- Return volume
154  inline scalar mag() const;
155 
156  //- Return circum-centre
157  inline Point circumCentre() const;
158 
159  //- Return circum-radius
160  inline scalar circumRadius() const;
161 
162  //- Return quality: Ratio of tetrahedron and circum-sphere
163  // volume, scaled so that a regular tetrahedron has a
164  // quality of 1
165  inline scalar quality() const;
166 
167  //- Return a random point in the tetrahedron from a
168  // uniform distribution
169  inline Point randomPoint(Random& rndGen) const;
170 
171  //- Calculate the point from the given barycentric coordinates.
172  inline Point barycentricToPoint(const barycentric& bary) const;
173 
174  //- Calculate the barycentric coordinates from the given point
175  inline barycentric pointToBarycentric(const point& pt) const;
176 
177  //- Calculate the barycentric coordinates from the given point.
178  // Returns the determinant.
179  inline scalar pointToBarycentric
180  (
181  const point& pt,
182  barycentric& bary
183  ) const;
184 
185  //- Return nearest point to p on tetrahedron. Is p itself
186  // if inside.
187  inline pointHit nearestPoint(const point& p) const;
188 
189  //- Return true if point is inside tetrahedron
190  inline bool inside(const point& pt) const;
191 
192  //- Return (min)containment sphere, i.e. the smallest sphere with
193  // all points inside. Returns pointHit with:
194  // - hit : if sphere is equal to circumsphere
195  // (biggest sphere)
196  // - point : centre of sphere
197  // - distance : radius of sphere
198  // - eligiblemiss: false
199  // Tol (small compared to 1, e.g. 1e-9) is used to determine
200  // whether point is inside: mag(pt - ctr) < (1+tol)*radius.
201  pointHit containmentSphere(const scalar tol) const;
202 
203  //- Fill buffer with shape function products
204  void gradNiSquared(scalarField& buffer) const;
205 
206  void gradNiDotGradNj(scalarField& buffer) const;
207 
208  void gradNiGradNi(tensorField& buffer) const;
209 
210  void gradNiGradNj(tensorField& buffer) const;
211 
212  //- Calculate the bounding box
213  boundBox bounds() const;
214 
215 
216  // IOstream Operators
217 
218  friend Istream& operator>> <Point, PointRef>
219  (
220  Istream&,
221  tetrahedron&
222  );
223 
224  friend Ostream& operator<< <Point, PointRef>
225  (
226  Ostream&,
227  const tetrahedron&
228  );
229 };
230 
231 
232 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
233 
234 } // End namespace Foam
235 
236 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
237 
238 #include "tetrahedronI.H"
239 
240 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
241 
242 #ifdef NoRepository
243  #include "tetrahedron.C"
244 #endif
245 
246 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
247 
248 #endif
249 
250 // ************************************************************************* //
A tetrahedron primitive.
Definition: tetrahedron.H:61
A triangle primitive used to calculate face areas and swept volumes.
Definition: triangle.H:57
pointHit containmentSphere(const scalar tol) const
Return (min)containment sphere, i.e. the smallest sphere with.
Definition: tetrahedron.C:34
scalar circumRadius() const
Return circum-radius.
Definition: tetrahedronI.H:204
tetrahedron(const Point &a, const Point &b, const Point &c, const Point &d)
Construct from points.
Definition: tetrahedronI.H:34
Point circumCentre() const
Return circum-centre.
Definition: tetrahedronI.H:177
const Point & c() const
Definition: tetrahedronI.H:86
vector Sd() const
Definition: tetrahedronI.H:156
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
void gradNiSquared(scalarField &buffer) const
Fill buffer with shape function products.
Definition: tetrahedron.C:246
const Point & a() const
Return vertices.
Definition: tetrahedronI.H:72
A bounding box defined in terms of the points at its extremities.
Definition: boundBox.H:58
Point barycentricToPoint(const barycentric &bary) const
Calculate the point from the given barycentric coordinates.
Definition: tetrahedronI.H:254
Random rndGen(label(0))
const Point & d() const
Definition: tetrahedronI.H:93
Forward declarations of the specialisations of Field<T> for scalar, vector and tensor.
const Point & b() const
Definition: tetrahedronI.H:79
scalar quality() const
Return quality: Ratio of tetrahedron and circum-sphere.
Definition: tetrahedronI.H:230
void gradNiGradNi(tensorField &buffer) const
Definition: tetrahedron.C:297
vector Sb() const
Definition: tetrahedronI.H:142
scalar mag() const
Return volume.
Definition: tetrahedronI.H:170
vector Sa() const
Return face normal.
Definition: tetrahedronI.H:135
Point centre() const
Return centre (centroid)
Definition: tetrahedronI.H:163
bool inside(const point &pt) const
Return true if point is inside tetrahedron.
Definition: tetrahedronI.H:413
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
Random number generator.
Definition: Random.H:57
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
triPointRef tri(const label facei) const
Return i-th face.
Definition: tetrahedronI.H:101
boundBox bounds() const
Calculate the bounding box.
Definition: tetrahedron.C:340
barycentric pointToBarycentric(const point &pt) const
Calculate the barycentric coordinates from the given point.
Definition: tetrahedronI.H:264
void gradNiGradNj(tensorField &buffer) const
Definition: tetrahedron.C:315
pointHit nearestPoint(const point &p) const
Return nearest point to p on tetrahedron. Is p itself.
Definition: tetrahedronI.H:319
void gradNiDotGradNj(scalarField &buffer) const
Definition: tetrahedron.C:267
Point randomPoint(Random &rndGen) const
Return a random point in the tetrahedron from a.
Definition: tetrahedronI.H:244
volScalarField & p
vector Sc() const
Definition: tetrahedronI.H:149
Namespace for OpenFOAM.