tetCellI.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 | Copyright (C) 2011-2012 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 Description
25 
26 \*---------------------------------------------------------------------------*/
27 
28 #include "IOstreams.H"
29 
30 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
31 
33 {}
34 
35 
37 (
38  const label a,
39  const label b,
40  const label c,
41  const label d
42 )
43 {
44  operator[](0) = a;
45  operator[](1) = b;
46  operator[](2) = c;
47  operator[](3) = d;
48 }
49 
50 
52 :
53  FixedList<label, 4>(lst)
54 {}
55 
56 
58 :
59  FixedList<label, 4>(is)
60 {}
61 
62 
63 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
64 
65 inline Foam::triFace Foam::tetCell::face(const label faceI) const
66 {
67  // Warning. Ordering of faces needs to be the same for a tetrahedron
68  // class, a tetrahedron cell shape model and a tetCell
69  static const label a[] = {1, 0, 0, 0};
70  static const label b[] = {2, 3, 1, 2};
71  static const label c[] = {3, 2, 3, 1};
72 
73 # ifdef FULLDEBUG
74  if (faceI >= 4)
75  {
76  FatalErrorIn("tetCell::face(const label faceI) const")
77  << "index out of range 0 -> 3. faceI = " << faceI
78  << abort(FatalError);
79  }
80 # endif
81 
82  return triFace
83  (
84  operator[](a[faceI]),
85  operator[](b[faceI]),
86  operator[](c[faceI])
87  );
88 }
89 
90 
91 inline Foam::label Foam::tetCell::edgeFace(const label edgeI) const
92 {
93  // Warning. Ordering of faces needs to be the same for a tetrahedron
94  // class, a tetrahedron cell shape model and a tetCell
95  //static const label edgeFaces[6] = {2, 1, 1, 0, 0, 0};
96  static const label edgeFaces[6] = {2, 3, 1, 0, 0, 1};
97 
98 # ifdef FULLDEBUG
99  if (edgeI >= 6)
100  {
102  (
103  "tetCell::edgeFace(const label edgeI)"
104  "const"
105  ) << "edge index out of range 0 -> 5. edgeI = " << edgeI
106  << abort(FatalError);
107  }
108 # endif
109 
110  return edgeFaces[edgeI];
111 }
112 
113 
115 (
116  const label edgeI,
117  const label faceI
118 ) const
119 {
120  // Warning. Ordering of faces needs to be the same for a tetrahedron
121  // class, a tetrahedron cell shape model and a tetCell
122  static const label adjacentFace[6][4] =
123  {
124  {-1, -1, 3, 2},
125  {-1, 3, -1, 1},
126  {-1, 2, 1, -1},
127  {2, -1, 0, -1},
128  {3, -1, -1, 0},
129  {1, 0, -1, -1}
130  };
131 
132 # ifdef FULLDEBUG
133  if (faceI >= 4)
134  {
136  (
137  "tetCell::edgeAdjacentFace(const label edgeI, const label faceI)"
138  "const"
139  ) << "face index out of range 0 -> 3. faceI = " << faceI
140  << abort(FatalError);
141  }
142 
143  if (edgeI >= 6)
144  {
146  (
147  "tetCell::edgeAdjacentFace(const label edgeI, const label faceI)"
148  "const"
149  ) << "edge index out of range 0 -> 5. edgeI = " << edgeI
150  << abort(FatalError);
151  }
152 # endif
153 
154  return adjacentFace[edgeI][faceI];
155 }
156 
157 
158 inline Foam::edge Foam::tetCell::tetEdge(const label edgeI) const
159 {
160  // Warning. Ordering of edges needs to be the same for a tetrahedron
161  // class, a tetrahedron cell shape model and a tetCell
162  //
163  static const label start[] = {0, 0, 0, 3, 1, 3};
164  static const label end[] = {1, 2, 3, 1, 2, 2};
165 
166 # ifdef FULLDEBUG
167  if (edgeI >= 6)
168  {
169  FatalErrorIn("tetCell::tetEdge(const label edgeI) const")
170  << "index out of range 0 -> 5. edgeI = " << edgeI
171  << abort(FatalError);
172  }
173 # endif
174 
175  return edge(operator[](start[edgeI]), operator[](end[edgeI]));
176 }
177 
178 
180 {
181  return tetPointRef
182  (
183  points[operator[](0)],
184  points[operator[](1)],
185  points[operator[](2)],
186  points[operator[](3)]
187  );
188 }
189 
190 
191 // ************************************************************************* //
edge tetEdge(const label edgeI) const
Return i-th edge.
Definition: tetCellI.H:158
const pointField & points
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
Definition: createFields.H:28
tetCell()
Construct null.
Definition: tetCellI.H:32
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
face triFace(3)
label edgeAdjacentFace(const label edgeI, const label faceI) const
Return face adjacent to the given face sharing the same edge.
Definition: tetCellI.H:115
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
tetrahedron< point, const point & > tetPointRef
Definition: tetrahedron.H:78
Useful combination of include files which define Sin, Sout and Serr and the use of IO streams general...
tetPointRef tet(const pointField &) const
Return the tetrahedron.
Definition: tetCellI.H:179
label & operator[](const label)
Return element of FixedList.
iterator end()
Return an iterator to end traversing the FixedList.
A tetrahedron primitive.
Definition: tetrahedron.H:62
errorManip< error > abort(error &err)
Definition: errorManip.H:131
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
error FatalError
triFace face(const label faceI) const
Return i-th face.
Definition: tetCellI.H:65
A triangular face using a FixedList of labels corresponding to mesh vertices.
Definition: triFace.H:68
const dimensionedScalar c
Speed of light in a vacuum.
label edgeFace(const label edgeI) const
Return first face adjacent to the given edge.
Definition: tetCellI.H:91