directionInfo.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-2022 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 Class
25  Foam::directionInfo
26 
27 Description
28  Holds direction in which to split cell (in fact a local coordinate axes).
29  Information is a label and a direction.
30 
31  The direction is the normal
32  direction to cut in. The label's meaning depends on whether the info
33  is on a cell or on a face:
34  - in cell: edge that is being cut. (determines for hex how cut is)
35  - in face: local face point that is being cut or -1.
36  -# (-1) : cut is tangential to plane
37  -# (>= 0): edge fp..fp+1 is cut
38 
39  (has to be facepoint, not vertex since vertex not valid across
40  processors whereas f[0] should correspond to f[0] on other side)
41 
42  The rule is that if the label is set (-1 or higher) it is used
43  (topological information only), otherwise the vector is used. This makes
44  sure that we use topological information as much as possible and so a
45  hex mesh is cut purely topologically. All other shapes are cut
46  geometrically.
47 
48 SourceFiles
49  directionInfoI.H
50  directionInfo.C
51 
52 \*---------------------------------------------------------------------------*/
53 
54 #ifndef directionInfo_H
55 #define directionInfo_H
56 
57 #include "point.H"
58 #include "labelList.H"
59 #include "tensor.H"
60 
61 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
62 
63 namespace Foam
64 {
65 class polyPatch;
66 class polyMesh;
67 class primitiveMesh;
68 class edge;
69 class face;
70 class polyMesh;
71 class transformer;
72 
73 
74 // Forward declaration of friend functions and operators
75 
76 class directionInfo;
77 
78 Istream& operator>>(Istream&, directionInfo&);
79 Ostream& operator<<(Ostream&, const directionInfo&);
80 
81 
82 /*---------------------------------------------------------------------------*\
83  Class directionInfo Declaration
84 \*---------------------------------------------------------------------------*/
85 
86 class directionInfo
87 {
88  // Private Data
89 
90  // Either mesh edge or face point
91  label index_;
92 
93  // Local n axis
94  vector n_;
95 
96 
97  // Private Member Functions
98 
99  //- Find edge among edgeLabels that uses v0 and v1
100  static label findEdge
101  (
102  const primitiveMesh& mesh,
103  const labelList& edgeLabels,
104  const label v1,
105  const label v0
106  );
107 
108  //- Return 'lowest' of a,b in face of size.
109  static label lowest
110  (
111  const label size,
112  const label a,
113  const label b
114  );
115 
116 public:
117 
118  // Static Functions
119 
120  //- Given edge on hex cell find corresponding edge on face. Is either
121  // index in face or -1 (cut tangential to face). Public since is
122  // needed to fill in seed faces in meshWave.
123  static label edgeToFaceIndex
124  (
125  const primitiveMesh& mesh,
126  const label celli,
127  const label facei,
128  const label edgeI
129  );
130 
131  // Constructors
132 
133  //- Construct null
134  inline directionInfo();
135 
136  //- Construct from components
137  inline directionInfo
138  (
139  const label,
140  const vector& n
141  );
142 
143 
144  // Member Functions
145 
146  // Access
148  inline label index() const
149  {
150  return index_;
151  }
153  inline const vector& n() const
154  {
155  return n_;
156  }
157 
158  // Needed by FaceCellWave
159 
160  //- Check whether origin has been changed at all or
161  // still contains original (invalid) value.
162  template<class TrackingData>
163  inline bool valid(TrackingData& td) const;
164 
165  //- Check for identical geometrical data. Used for cyclics checking.
166  template<class TrackingData>
167  inline bool sameGeometry
168  (
169  const polyMesh&,
170  const directionInfo&,
171  const scalar,
172  TrackingData& td
173  ) const;
174 
175  //- Transform across an interface
176  template<class TrackingData>
177  inline void transform
178  (
179  const polyPatch& patch,
180  const label patchFacei,
181  const transformer& transform,
182  TrackingData& td
183  );
184 
185  //- Influence of neighbouring face.
186  template<class TrackingData>
187  inline bool updateCell
188  (
189  const polyMesh&,
190  const label thisCelli,
191  const label neighbourFacei,
192  const directionInfo& neighbourInfo,
193  const scalar tol,
194  TrackingData& td
195  );
196 
197  //- Influence of neighbouring cell.
198  template<class TrackingData>
199  inline bool updateFace
200  (
201  const polyMesh&,
202  const label thisFacei,
203  const label neighbourCelli,
204  const directionInfo& neighbourInfo,
205  const scalar tol,
206  TrackingData& td
207  );
208 
209  //- Influence of different value on same face.
210  template<class TrackingData>
211  inline bool updateFace
212  (
213  const polyMesh&,
214  const label thisFacei,
215  const directionInfo& neighbourInfo,
216  const scalar tol,
217  TrackingData& td
218  );
219 
220  //- Same (like operator==)
221  template<class TrackingData>
222  inline bool equal(const directionInfo&, TrackingData& td) const;
223 
224  // Member Operators
225 
226  // Needed for List IO
227  inline bool operator==(const directionInfo&) const;
228 
229  inline bool operator!=(const directionInfo&) const;
230 
231 
232  // IOstream Operators
233 
234  friend Ostream& operator<<(Ostream&, const directionInfo&);
236 };
237 
238 
239 //- Data associated with directionInfo type are contiguous
240 template<>
241 inline bool contiguous<directionInfo>()
242 {
243  return true;
244 }
245 
246 
247 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
248 
249 } // End namespace Foam
250 
251 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
252 
253 #include "directionInfoI.H"
254 
255 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
256 
257 #endif
258 
259 // ************************************************************************* //
bool sameGeometry(const polyMesh &, const directionInfo &, const scalar, TrackingData &td) const
Check for identical geometrical data. Used for cyclics checking.
Vector-tensor class used to perform translations, rotations and scaling operations in 3D space...
Definition: transformer.H:83
bool operator!=(const directionInfo &) const
bool equal(const directionInfo &, TrackingData &td) const
Same (like operator==)
Cell-face mesh analysis engine.
Definition: primitiveMesh.H:74
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
bool operator==(const directionInfo &) const
const dimensionedScalar b
Wien displacement law constant: default SI units: [m K].
Definition: createFields.H:27
friend Ostream & operator<<(Ostream &, const directionInfo &)
bool updateCell(const polyMesh &, const label thisCelli, const label neighbourFacei, const directionInfo &neighbourInfo, const scalar tol, TrackingData &td)
Influence of neighbouring face.
bool contiguous< directionInfo >()
Data associated with directionInfo type are contiguous.
fvMesh & mesh
void transform(const polyPatch &patch, const label patchFacei, const transformer &transform, TrackingData &td)
Transform across an interface.
bool valid(TrackingData &td) const
Check whether origin has been changed at all or.
Istream & operator>>(Istream &, directionInfo &)
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
friend Istream & operator>>(Istream &, directionInfo &)
label index() const
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
Holds direction in which to split cell (in fact a local coordinate axes). Information is a label and ...
Definition: directionInfo.H:85
directionInfo()
Construct null.
const vector & n() const
bool updateFace(const polyMesh &, const label thisFacei, const label neighbourCelli, const directionInfo &neighbourInfo, const scalar tol, TrackingData &td)
Influence of neighbouring cell.
Ostream & operator<<(Ostream &, const ensightPart &)
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:76
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
Namespace for OpenFOAM.