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