directionInfoI.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-2019 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 "polyMesh.H"
27 #include "meshTools.H"
28 #include "hexMatcher.H"
29 
30 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
31 
33 :
34  index_(-3),
35  n_(Zero)
36 {}
37 
38 
40 (
41  const label index,
42  const vector& n
43 )
44 :
45  index_(index),
46  n_(n)
47 {}
48 
49 
50 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
51 
52 template<class TrackingData>
53 inline bool Foam::directionInfo::valid(TrackingData& td) const
54 {
55  return index_ != -3;
56 }
57 
58 
59 template<class TrackingData>
61 (
62  const polyMesh&,
63  const directionInfo& w2,
64  const scalar tol,
65  TrackingData& td
66 )
67  const
68 {
69  return true;
70 }
71 
72 
73 template<class TrackingData>
75 (
76  const polyMesh&,
77  const polyPatch& patch,
78  const label patchFacei,
79  const point& faceCentre,
80  TrackingData& td
81 )
82 {}
83 
84 
85 template<class TrackingData>
87 (
88  const polyMesh&,
89  const polyPatch& patch,
90  const label patchFacei,
91  const point& faceCentre,
92  TrackingData& td
93 )
94 {
95  if (index_ >= 0)
96  {
97  const face& f = patch[patchFacei];
98 
99  index_ = (f.size() - index_) % f.size();
100  }
101 }
102 
103 
104 template<class TrackingData>
106 (
107  const polyMesh&,
108  const tensor& rotTensor,
109  TrackingData& td
110 )
111 {}
112 
113 
114 template<class TrackingData>
116 (
117  const polyMesh& mesh,
118  const label thisCelli,
119  const label neighbourFacei,
120  const directionInfo& neighbourInfo,
121  const scalar, // tol
122  TrackingData& td
123 )
124 {
125  if (index_ >= -2)
126  {
127  // Already determined.
128  return false;
129  }
130 
131  if (hexMatcher().isA(mesh, thisCelli))
132  {
133  const face& f = mesh.faces()[neighbourFacei];
134 
135  if (neighbourInfo.index() == -2)
136  {
137  // Geometric information from neighbour
138  index_ = -2;
139  }
140  else if (neighbourInfo.index() == -1)
141  {
142  // Cut tangential to face. Take any edge connected to face
143  // but not used in face.
144 
145  // Get first edge on face.
146  label edgeI = mesh.faceEdges()[neighbourFacei][0];
147 
148  const edge& e = mesh.edges()[edgeI];
149 
150  // Find face connected to face through edgeI and on same cell.
151  label facei =
153  (
154  mesh,
155  thisCelli,
156  neighbourFacei,
157  edgeI
158  );
159 
160  // Find edge on facei which is connected to e.start() but not edgeI.
161  index_ =
163  (
164  mesh,
165  mesh.faceEdges()[facei],
166  edgeI,
167  e.start()
168  );
169  }
170  else
171  {
172  // Index is a vertex on the face. Convert to mesh edge.
173 
174  // Get mesh edge between f[index_] and f[index_+1]
175  label v0 = f[neighbourInfo.index()];
176  label v1 = f[(neighbourInfo.index() + 1) % f.size()];
177 
178  index_ = findEdge(mesh, mesh.faceEdges()[neighbourFacei], v0, v1);
179  }
180  }
181  else
182  {
183  // Not a hex so mark this as geometric.
184  index_ = -2;
185  }
186 
187 
188  n_ = neighbourInfo.n();
189 
190  return true;
191 }
192 
193 
194 template<class TrackingData>
196 (
197  const polyMesh& mesh,
198  const label thisFacei,
199  const label neighbourCelli,
200  const directionInfo& neighbourInfo,
201  const scalar, // tol
202  TrackingData& td
203 )
204 {
205  // Handle special cases first
206 
207  if (index_ >= -2)
208  {
209  // Already determined
210  return false;
211  }
212 
213  // Handle normal cases where topological or geometrical info comes from
214  // neighbouring cell
215 
216  if (neighbourInfo.index() >= 0)
217  {
218  // Neighbour has topological direction (and hence is hex). Find cut
219  // edge on face.
220  index_ =
222  (
223  mesh,
224  neighbourCelli,
225  thisFacei,
226  neighbourInfo.index()
227  );
228  }
229  else
230  {
231  // Neighbour has geometric information. Use.
232  index_ = -2;
233  }
234 
235 
236  n_ = neighbourInfo.n();
237 
238  return true;
239 }
240 
241 
242 template<class TrackingData>
244 (
245  const polyMesh& mesh,
246  const label, // thisFacei
247  const directionInfo& neighbourInfo,
248  const scalar, // tol
249  TrackingData& td
250 )
251 {
252  if (index_ >= -2)
253  {
254  // Already visited.
255  return false;
256  }
257  else
258  {
259  index_ = neighbourInfo.index();
260 
261  n_ = neighbourInfo.n();
262 
263  return true;
264  }
265 }
266 
267 
268 template<class TrackingData>
269 inline bool Foam::directionInfo::equal
270 (
271  const directionInfo& rhs,
272  TrackingData& td
273 ) const
274 {
275  return operator==(rhs);
276 }
277 
278 
279 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
280 
281 inline bool Foam::directionInfo::operator==
282 (
283  const Foam::directionInfo& rhs
284 ) const
285 {
286  return index() == rhs.index() && n() == rhs.n();
287 }
288 
289 
290 inline bool Foam::directionInfo::operator!=
291 (
292  const Foam::directionInfo& rhs
293 ) const
294 {
295  return !(*this == rhs);
296 }
297 
298 
299 // ************************************************************************* //
bool isA(const Type &t)
Check if a dynamic_cast to typeid is possible.
Definition: typeInfo.H:134
bool sameGeometry(const polyMesh &, const directionInfo &, const scalar, TrackingData &td) const
Check for identical geometrical data. Used for cyclics checking.
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
label otherEdge(const primitiveMesh &, const labelList &edgeLabels, const label edgeI, const label vertI)
Return label of other edge (out of candidates edgeLabels)
Definition: meshTools.C:501
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
bool equal(const directionInfo &, TrackingData &td) const
Same (like operator==)
const labelListList & faceEdges() const
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
bool operator==(const directionInfo &) const
A cellMatcher for hex cells.
Definition: hexMatcher.H:51
void leaveDomain(const polyMesh &, const polyPatch &, const label patchFacei, const point &faceCentre, TrackingData &td)
Convert any absolute coordinates into relative to (patch)face.
bool updateCell(const polyMesh &, const label thisCelli, const label neighbourFacei, const directionInfo &neighbourInfo, const scalar tol, TrackingData &td)
Influence of neighbouring face.
bool valid(TrackingData &td) const
Check whether origin has been changed at all or.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
void enterDomain(const polyMesh &, const polyPatch &, const label patchFacei, const point &faceCentre, TrackingData &td)
Reverse of leaveDomain.
static label edgeToFaceIndex(const primitiveMesh &mesh, const label celli, const label facei, const label edgeI)
Given edge on hex cell find corresponding edge on face. Is either.
Definition: directionInfo.C:93
static const zero Zero
Definition: zero.H:97
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1156
label index() const
labelList f(nPoints)
Holds direction in which to split cell (in fact a local coordinate axes). Information is a label and ...
Definition: directionInfo.H:84
directionInfo()
Construct null.
const vector & n() const
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
void transform(const polyMesh &, const tensor &, TrackingData &td)
Apply rotation matrix to any coordinates.
bool updateFace(const polyMesh &, const label thisFacei, const label neighbourCelli, const directionInfo &neighbourInfo, const scalar tol, TrackingData &td)
Influence of neighbouring cell.
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:98
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
label start() const
Return start vertex label.
Definition: edgeI.H:81
label otherFace(const primitiveMesh &, const label celli, const label facei, const label edgeI)
Return face on cell using edgeI but not facei. Throws error.
Definition: meshTools.C:536