ENH: autoLayerDriver: better layering information message
[OpenFOAM-2.0.x.git] / src / OpenFOAM / meshes / meshShapes / face / face.H
blobc1c9716c03f4f55ee1f62e232f94b168d20c8c1e
1 /*---------------------------------------------------------------------------*\
2   =========                 |
3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
4    \\    /   O peration     |
5     \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
6      \\/     M anipulation  |
7 -------------------------------------------------------------------------------
8 License
9     This file is part of OpenFOAM.
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.
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.
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/>.
24 Class
25     Foam::face
27 Description
28     A face is a list of labels corresponding to mesh vertices.
30 SeeAlso
31     Foam::triFace
33 SourceFiles
34     faceI.H
35     face.C
36     faceIntersection.C
37     faceContactSphere.C
38     faceAreaInContact.C
39     faceTemplates.C
41 \*---------------------------------------------------------------------------*/
43 #ifndef face_H
44 #define face_H
46 #include "pointField.H"
47 #include "labelList.H"
48 #include "edgeList.H"
49 #include "vectorField.H"
50 #include "faceListFwd.H"
51 #include "intersection.H"
52 #include "pointHit.H"
54 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
56 namespace Foam
59 // Forward declaration of friend functions and operators
61 class face;
62 class triFace;
64 template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
65 class DynamicList;
67 inline bool operator==(const face& a, const face& b);
68 inline bool operator!=(const face& a, const face& b);
69 inline Istream& operator>>(Istream&, face&);
71 /*---------------------------------------------------------------------------*\
72                            Class face Declaration
73 \*---------------------------------------------------------------------------*/
75 class face
77     public labelList
79     // Private Member Functions
81         //- Edge to the right of face vertex i
82         inline label right(const label i) const;
84         //- Edge to the left of face vertex i
85         inline label left(const label i) const;
87         //- Construct list of edge vectors for face
88         tmp<vectorField> calcEdges
89         (
90             const pointField& points
91         ) const;
93         //- Cos between neighbouring edges
94         scalar edgeCos
95         (
96             const vectorField& edges,
97             const label index
98         ) const;
100         //- Find index of largest internal angle on face
101         label mostConcaveAngle
102         (
103             const pointField& points,
104             const vectorField& edges,
105             scalar& edgeCos
106         ) const;
108         //- Enumeration listing the modes for split()
109         enum splitMode
110         {
111             COUNTTRIANGLE,  // count if split into triangles
112             COUNTQUAD,      // count if split into triangles&quads
113             SPLITTRIANGLE,  // split into triangles
114             SPLITQUAD       // split into triangles&quads
115         };
117         //- Split face into triangles or triangles&quads.
118         //  Stores results quadFaces[quadI], triFaces[triI]
119         //  Returns number of new faces created
120         label split
121         (
122             const splitMode mode,
123             const pointField& points,
124             label& triI,
125             label& quadI,
126             faceList& triFaces,
127             faceList& quadFaces
128         ) const;
131 public:
133     //- Return types for classify
134     enum proxType
135     {
136         NONE,
137         POINT,  // Close to point
138         EDGE    // Close to edge
139     };
141     // Static data members
143         static const char* const typeName;
146     // Constructors
148         //- Construct null
149         inline face();
151         //- Construct given size
152         explicit inline face(label);
154         //- Construct from list of labels
155         explicit inline face(const labelUList&);
157         //- Construct from list of labels
158         explicit inline face(const labelList&);
160         //- Construct by transferring the parameter contents
161         explicit inline face(const Xfer<labelList>&);
163         //- Copy construct from triFace
164         face(const triFace&);
166         //- Construct from Istream
167         inline face(Istream&);
170     // Member Functions
172         //- Collapse face by removing duplicate point labels
173         //  return the collapsed size
174         label collapse();
176         //- Flip the face in-place.
177         //  The starting points of the original and flipped face are identical.
178         void flip();
180         //- Return the points corresponding to this face
181         inline pointField points(const pointField&) const;
183         //- Centre point of face
184         point centre(const pointField&) const;
186         //- Calculate average value at centroid of face
187         template<class Type>
188         Type average(const pointField&, const Field<Type>&) const;
190         //- Magnitude of face area
191         inline scalar mag(const pointField&) const;
193         //- Vector normal; magnitude is equal to area of face
194         vector normal(const pointField&) const;
196         //- Return face with reverse direction
197         //  The starting points of the original and reverse face are identical.
198         face reverseFace() const;
200         //- Navigation through face vertices
202             //- Which vertex on face (face index given a global index)
203             //  returns -1 if not found
204             label which(const label globalIndex) const;
206             //- Next vertex on face
207             inline label nextLabel(const label i) const;
209             //- Previous vertex on face
210             inline label prevLabel(const label i) const;
213         //- Return the volume swept out by the face when its points move
214         scalar sweptVol
215         (
216             const pointField& oldPoints,
217             const pointField& newPoints
218         ) const;
220         //- Return the inertia tensor, with optional reference
221         //  point and density specification
222         tensor inertia
223         (
224             const pointField&,
225             const point& refPt = vector::zero,
226             scalar density = 1.0
227         ) const;
229         //- Return potential intersection with face with a ray starting
230         //  at p, direction n (does not need to be normalized)
231         //  Does face-centre decomposition and returns triangle intersection
232         //  point closest to p. Face-centre is calculated from point average.
233         //  For a hit, the distance is signed.  Positive number
234         //  represents the point in front of triangle
235         //  In case of miss the point is the nearest point on the face
236         //  and the distance is the distance between the intersection point
237         //  and the original point.
238         //  The half-ray or full-ray intersection and the contact
239         //  sphere adjustment of the projection vector is set by the
240         //  intersection parameters
241         pointHit ray
242         (
243             const point& p,
244             const vector& n,
245             const pointField&,
246             const intersection::algorithm alg = intersection::FULL_RAY,
247             const intersection::direction dir = intersection::VECTOR
248         ) const;
250         //- Fast intersection with a ray.
251         //  Does face-centre decomposition and returns triangle intersection
252         //  point closest to p. See triangle::intersection for details.
253         pointHit intersection
254         (
255             const point& p,
256             const vector& q,
257             const point& ctr,
258             const pointField&,
259             const intersection::algorithm alg,
260             const scalar tol = 0.0
261         ) const;
263         //- Return nearest point to face
264         pointHit nearestPoint
265         (
266             const point& p,
267             const pointField&
268         ) const;
270         //- Return nearest point to face and classify it:
271         //  + near point (nearType=POINT, nearLabel=0, 1, 2)
272         //  + near edge (nearType=EDGE, nearLabel=0, 1, 2)
273         //    Note: edges are counted from starting vertex so
274         //    e.g. edge n is from f[n] to f[0], where the face has n + 1
275         //    points
276         pointHit nearestPointClassify
277         (
278             const point& p,
279             const pointField&,
280             label& nearType,
281             label& nearLabel
282         ) const;
284         //- Return contact sphere diameter
285         scalar contactSphereDiameter
286         (
287             const point& p,
288             const vector& n,
289             const pointField&
290         ) const;
292         //- Return area in contact, given the displacement in vertices
293         scalar areaInContact
294         (
295             const pointField&,
296             const scalarField& v
297         ) const;
299         //- Return number of edges
300         inline label nEdges() const;
302         //- Return edges in face point ordering,
303         //  i.e. edges()[0] is edge between [0] and [1]
304         edgeList edges() const;
306         //- Return n-th face edge
307         inline edge faceEdge(const label n) const;
309         //- Return the edge direction on the face
310         //  Returns:
311         //  -  0: edge not found on the face
312         //  - +1: forward (counter-clockwise) on the face
313         //  - -1: reverse (clockwise) on the face
314         int edgeDirection(const edge&) const;
316         // Face splitting utilities
318             //- Number of triangles after splitting
319             inline label nTriangles() const;
321             //- Number of triangles after splitting
322             label nTriangles(const pointField& points) const;
324             //- Split into triangles using existing points.
325             //  Result in triFaces[triI..triI+nTri].
326             //  Splits intelligently to maximize triangle quality.
327             //  Returns number of faces created.
328             label triangles
329             (
330                 const pointField& points,
331                 label& triI,
332                 faceList& triFaces
333             ) const;
335             //- Split into triangles using existing points.
336             //  Append to DynamicList.
337             //  Returns number of faces created.
338             template<unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
339             label triangles
340             (
341                 const pointField& points,
342                 DynamicList<face, SizeInc, SizeMult, SizeDiv>& triFaces
343             ) const;
345             //- Number of triangles and quads after splitting
346             //  Returns the sum of both
347             label nTrianglesQuads
348             (
349                 const pointField& points,
350                 label& nTris,
351                 label& nQuads
352             ) const;
354             //- Split into triangles and quads.
355             //  Results in triFaces (starting at triI) and quadFaces
356             //  (starting at quadI).
357             //  Returns number of new faces created.
358             label trianglesQuads
359             (
360                 const pointField& points,
361                 label& triI,
362                 label& quadI,
363                 faceList& triFaces,
364                 faceList& quadFaces
365             ) const;
367         //- compare faces
368         //   0: different
369         //  +1: identical
370         //  -1: same face, but different orientation
371         static int compare(const face&, const face&);
374     // Friend Operators
376         friend bool operator==(const face& a, const face& b);
377         friend bool operator!=(const face& a, const face& b);
380     // Istream Operator
382         friend Istream& operator>>(Istream&, face&);
386 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
388 } // End namespace Foam
390 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
392 #include "faceI.H"
394 #ifdef NoRepository
395 #   include "faceTemplates.C"
396 #endif
398 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
400 #endif
402 // ************************************************************************* //