tetIndicesI.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-2021 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 \*---------------------------------------------------------------------------*/
25 
26 #include "tetIndices.H"
27 #include "Time.H"
28 
29 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
30 
32 {
33  return celli_;
34 }
35 
36 
38 {
39  return celli_;
40 }
41 
42 
44 {
45  return facei_;
46 }
47 
48 
50 {
51  return facei_;
52 }
53 
54 
56 {
57  return tetPti_;
58 }
59 
60 
62 {
63  return tetPti_;
64 }
65 
66 
68 {
69  const Foam::face& f = mesh.faces()[face()];
70 
71  label faceBasePtI = mesh.tetBasePtIs()[face()];
72 
73  if (faceBasePtI < 0)
74  {
75  static labelHashSet badFaces;
76  static label badTimeIndex = -1;
77 
78  if (badTimeIndex != mesh.time().timeIndex())
79  {
80  badFaces.clear();
81  badTimeIndex = mesh.time().timeIndex();
82  }
83 
84  if (!badFaces[face()])
85  {
87  << "No base point for face " << face() << ", " << f
88  << ", produces a valid tet decomposition." << endl;
89 
90  badFaces.set(face());
91  }
92 
93  faceBasePtI = 0;
94  }
95 
96  label facePtI = (tetPt() + faceBasePtI) % f.size();
97  label faceOtherPtI = f.fcIndex(facePtI);
98 
99  if (mesh.faceOwner()[face()] != cell())
100  {
101  Swap(facePtI, faceOtherPtI);
102  }
103 
104  return triFace(f[faceBasePtI], f[facePtI], f[faceOtherPtI]);
105 }
106 
107 
109 {
110  const pointField& meshPoints = mesh.points();
111  const triFace tri = faceTriIs(mesh);
112 
113  return tetPointRef
114  (
115  mesh.cellCentres()[cell()],
116  meshPoints[tri[0]],
117  meshPoints[tri[1]],
118  meshPoints[tri[2]]
119  );
120 }
121 
122 
124 {
125  const pointField& meshPoints = mesh.points();
126  const triFace tri = faceTriIs(mesh);
127 
128  return triPointRef
129  (
130  meshPoints[tri[0]],
131  meshPoints[tri[1]],
132  meshPoints[tri[2]]
133  );
134 }
135 
136 
138 (
139  const polyMesh& mesh
140 ) const
141 {
142  const pointField& meshOldPoints = mesh.oldPoints();
143  const triFace tri = faceTriIs(mesh);
144 
145  return triPointRef
146  (
147  meshOldPoints[tri[0]],
148  meshOldPoints[tri[1]],
149  meshOldPoints[tri[2]]
150  );
151 }
152 
153 
154 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
155 
156 inline bool Foam::tetIndices::operator==(const Foam::tetIndices& rhs) const
157 {
158  return
159  cell() == rhs.cell()
160  && face() == rhs.face()
161  && tetPt() == rhs.tetPt();
162 }
163 
164 
165 inline bool Foam::tetIndices::operator!=(const Foam::tetIndices& rhs) const
166 {
167  return !(*this == rhs);
168 }
169 
170 
171 // ************************************************************************* //
bool set(const Key &key)
Same as insert (cannot overwrite nil content)
Definition: HashSet.H:121
A tetrahedron primitive.
Definition: tetrahedron.H:61
tetrahedron< point, const point & > tetPointRef
Definition: tetPointRef.H:44
label tetPt() const
Return the characterising tetPtI.
Definition: tetIndicesI.H:55
A triangle primitive used to calculate face areas and swept volumes.
Definition: triangle.H:57
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
const labelIOList & tetBasePtIs() const
Return the tetBasePtIs.
Definition: polyMesh.C:928
triPointRef faceTri(const polyMesh &mesh) const
Return the geometry corresponding to the tri on the face for.
Definition: tetIndicesI.H:123
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
label fcIndex(const label i) const
Return the forward circular index, i.e. the next index.
Definition: UListI.H:58
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
triFace faceTriIs(const polyMesh &mesh) const
Return the indices corresponding to the tri on the face for.
Definition: tetIndicesI.H:67
tetPointRef tet(const polyMesh &mesh) const
Return the geometry corresponding to this tet.
Definition: tetIndicesI.H:108
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1211
void Swap(T &a, T &b)
Definition: Swap.H:43
bool operator==(const tetIndices &) const
Definition: tetIndicesI.H:156
A triangular face using a FixedList of labels corresponding to mesh vertices.
Definition: triFace.H:68
label cell() const
Return the cell.
Definition: tetIndicesI.H:31
face triFace(3)
void clear()
Clear all entries from table.
Definition: HashTable.C:468
virtual const pointField & oldPoints() const
Return old points for mesh motion.
Definition: polyMesh.C:1249
virtual const labelList & faceOwner() const
Return face owner.
Definition: polyMesh.C:1237
Storage and named access for the indices of a tet which is part of the decomposition of a cell...
Definition: tetIndices.H:81
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1224
const vectorField & cellCentres() const
label face() const
Return the face.
Definition: tetIndicesI.H:43
const Time & time() const
Return time.
labelList f(nPoints)
bool operator!=(const tetIndices &) const
Definition: tetIndicesI.H:165
label timeIndex() const
Return current time index.
Definition: TimeStateI.H:28
#define WarningInFunction
Report a warning using Foam::Warning.
triangle< point, const point & > triPointRef
Definition: triPointRef.H:44
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:76
triPointRef oldFaceTri(const polyMesh &mesh) const
Return the geometry corresponding to the tri on the face for.
Definition: tetIndicesI.H:138