tetIndices.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) 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
27Class
28 Foam::tetIndices
29
30Description
31 Storage and named access for the indices of a tet which is part of
32 the decomposition of a cell.
33
34 Tets are designated by
35 - cell (of course)
36 - face on cell
37 - three points on face (faceBasePt, facePtA, facePtB)
38 When constructing from a mesh and index in the face (tetPtI):
39 - faceBasePt is the mesh.tetBasePtIs() base point
40 - facePtA is tetPtI away from faceBasePt
41 - facePtB is next one after/before facePtA
42 e.g.:
43
44 +---+
45 |2 /|
46 | / |
47 |/ 1| <- tetPt (so 1 for first triangle, 2 for second)
48 +---+
49 ^
50 faceBasePt
51
52SourceFiles
53 tetIndicesI.H
54 tetIndices.C
55
56\*---------------------------------------------------------------------------*/
57
58#ifndef Foam_tetIndices_H
59#define Foam_tetIndices_H
60
61#include "label.H"
62#include "tetPointRef.H"
63#include "triPointRef.H"
64#include "polyMesh.H"
65#include "triFace.H"
66#include "face.H"
67
68// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
69
70namespace Foam
71{
72
73// Forward Declarations
74class tetIndices;
75
76Istream& operator>>(Istream&, tetIndices&);
77Ostream& operator<<(Ostream&, const tetIndices&);
78
79
80/*---------------------------------------------------------------------------*\
81 Class tetIndices Declaration
82\*---------------------------------------------------------------------------*/
84class tetIndices
85{
86 // Private Data
87
88 //- Cell that this is a decomposed tet of
89 label celli_;
90
91 //- Face that holds this decomposed tet
92 label facei_;
93
94 //- Point on the face, *relative to the base point*,
95 //- which characterises this tet on the face.
96 label tetPti_;
97
98
99 // Private Static Data
100
101 //- Maximum number of bad tet warnings
102 static constexpr int maxNWarnings = 100;
103
104 //- Current number of bad tet warnings.
105 //- Warnings stop when this reaches the maximum number.
106 static int nWarnings_;
107
108
109public:
110
111 // Constructors
112
113 //- Default construct, with invalid labels (-1)
114 inline constexpr tetIndices() noexcept;
115
116 //- Construct from components
117 inline constexpr tetIndices
118 (
119 label celli,
120 label facei,
121 label tetPointi
122 ) noexcept;
123
124
125 //- Destructor
126 ~tetIndices() = default;
127
128
129 // Member Functions
130
131 // Access
132
133 //- Return the cell index
134 label cell() const noexcept { return celli_; }
135
136 //- Non-const access to the cell index
137 label& cell() noexcept { return celli_; }
138
139 //- Return the face index
140 label face() const noexcept { return facei_; }
141
142 //- Non-const access to the face index
143 label& face() noexcept { return facei_; }
144
145 //- Return the characterising tet point index
146 label tetPt() const noexcept { return tetPti_; }
147
148 //- Non-const access to the characterising tet point index
149 label& tetPt() noexcept { return tetPti_; }
150
151
152 // Searching
153
154 //- Return the indices corresponding to the tri on the face for
155 // this tet. The normal of the tri points out of the cell
156 inline triFace faceTriIs
157 (
158 const polyMesh& mesh,
159 const bool warn = true
160 ) const;
161
162 //- Return the local indices corresponding to the tri on the face
163 //- for this tet. The normal of the tri points out of the cell
164 inline triFace triIs
165 (
166 const polyMesh& mesh,
167 const bool warn = true
168 ) const;
169
170 //- Return the geometry corresponding to this tet
171 inline tetPointRef tet(const polyMesh& mesh) const;
172
173 //- Return the geometry corresponding to the tri on the face for
174 //- this tet. The normal of the tri points out of the cell
175 inline triPointRef faceTri(const polyMesh& mesh) const;
176
177 //- Return the geometry corresponding to the tri on the face for
178 //- this tet using the old positions
179 inline triPointRef oldFaceTri(const polyMesh& mesh) const;
180
181
182 // Other
183
184 //- Compare tetIndices for equality.
185 //- Compares cell, face, tetPt elements in order, stopping at the
186 //- first inequality.
187 //
188 // \returns negative/zero/positive from the last element compared
189 static inline label compare
190 (
191 const tetIndices& a,
192 const tetIndices& b
193 ) noexcept;
194
195
196 // IOstream Operators
199 friend Ostream& operator<<(Ostream&, const tetIndices&);
200};
201
202
203// * * * * * * * * * * * * * * * * * Traits * * * * * * * * * * * * * * * * //
204
205//- Contiguous data for tetIndices
206template<> struct is_contiguous<tetIndices> : std::true_type {};
207
208//- Contiguous label data for tetIndices
209template<> struct is_contiguous_label<tetIndices> : std::true_type {};
210
211
212// * * * * * * * * * * * * * * * Global Operators * * * * * * * * * * * * * //
213
214inline bool operator==(const tetIndices& a, const tetIndices& b) noexcept;
215inline bool operator!=(const tetIndices& a, const tetIndices& b) noexcept;
216
218inline bool operator<(const tetIndices& a, const tetIndices& b) noexcept
219{
220 return (tetIndices::compare(a, b) < 0);
221}
223inline bool operator<=(const tetIndices& a, const tetIndices& b) noexcept
224{
225 return (tetIndices::compare(a, b) <= 0);
226}
228inline bool operator>(const tetIndices& a, const tetIndices& b) noexcept
229{
230 return (tetIndices::compare(a, b) > 0);
231}
233inline bool operator>=(const tetIndices& a, const tetIndices& b) noexcept
234{
235 return (tetIndices::compare(a, b) >= 0);
236}
237
238
239// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
240
241} // End namespace Foam
242
243// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
244
245#include "tetIndicesI.H"
246
247// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
248
249#endif
250
251// ************************************************************************* //
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
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
Storage and named access for the indices of a tet which is part of the decomposition of a cell.
Definition: tetIndices.H:84
label face() const noexcept
Return the face index.
Definition: tetIndices.H:139
triFace triIs(const polyMesh &mesh, const bool warn=true) const
Definition: tetIndicesI.H:95
~tetIndices()=default
Destructor.
static label compare(const tetIndices &a, const tetIndices &b) noexcept
Definition: tetIndicesI.H:181
label & face() noexcept
Non-const access to the face index.
Definition: tetIndices.H:142
friend Ostream & operator<<(Ostream &, const tetIndices &)
triFace faceTriIs(const polyMesh &mesh, const bool warn=true) const
Return the indices corresponding to the tri on the face for.
Definition: tetIndicesI.H:55
tetPointRef tet(const polyMesh &mesh) const
Return the geometry corresponding to this tet.
Definition: tetIndicesI.H:134
triPointRef faceTri(const polyMesh &mesh) const
Definition: tetIndicesI.H:149
label tetPt() const noexcept
Return the characterising tet point index.
Definition: tetIndices.H:145
label cell() const noexcept
Return the cell index.
Definition: tetIndices.H:133
constexpr tetIndices() noexcept
Default construct, with invalid labels (-1)
Definition: tetIndicesI.H:31
friend Istream & operator>>(Istream &, tetIndices &)
triPointRef oldFaceTri(const polyMesh &mesh) const
Definition: tetIndicesI.H:164
label & tetPt() noexcept
Non-const access to the characterising tet point index.
Definition: tetIndices.H:148
label & cell() noexcept
Non-const access to the cell index.
Definition: tetIndices.H:136
A tetrahedron primitive.
Definition: tetrahedron.H:87
A triangular face using a FixedList of labels corresponding to mesh vertices.
Definition: triFace.H:72
A triangle primitive used to calculate face normals and swept volumes.
Definition: triangle.H:80
dynamicFvMesh & mesh
Namespace for OpenFOAM.
bool operator!=(const eddy &a, const eddy &b)
Definition: eddy.H:239
bool operator<=(const IOstreamOption::versionNumber &a, const IOstreamOption::versionNumber &b) noexcept
Version A same or older than B.
bool operator>=(const IOstreamOption::versionNumber &a, const IOstreamOption::versionNumber &b) noexcept
Version A same or newer than B.
Ostream & operator<<(Ostream &, const boundaryPatch &p)
Write boundaryPatch as dictionary entries (without surrounding braces)
Definition: boundaryPatch.C:83
tmp< faMatrix< Type > > operator==(const faMatrix< Type > &, const faMatrix< Type > &)
Istream & operator>>(Istream &, directionInfo &)
const direction noexcept
Definition: Scalar.H:223
bool operator<(const IOstreamOption::versionNumber &a, const IOstreamOption::versionNumber &b) noexcept
Version A older than B.
bool operator>(const IOstreamOption::versionNumber &a, const IOstreamOption::versionNumber &b) noexcept
Version A newer than B.
volScalarField & b
Definition: createFields.H:27
A template class to specify if a data type is composed solely of Foam::label elements.
Definition: contiguous.H:86
A template class to specify that a data type can be considered as being contiguous in memory.
Definition: contiguous.H:78