FreeFOAM The Cross-Platform CFD Toolkit
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 | Copyright (C) 1991-2010 OpenCFD Ltd.
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 <OpenFOAM/point.H>
58 #include <OpenFOAM/labelList.H>
59 #include <OpenFOAM/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 
72 /*---------------------------------------------------------------------------*\
73  Class directionInfo Declaration
74 \*---------------------------------------------------------------------------*/
75 
77 {
78  // Private data
79 
80  // Either mesh edge or face point
81  label index_;
82 
83  // Local n axis
84  vector n_;
85 
86 
87  // Private Member Functions
88 
89  //- edge uses two labels
90  static bool equal(const edge& e, const label, const label);
91 
92  //- Calculate mid point of edge.
93  static point eMid(const primitiveMesh& mesh, const label edgeI);
94 
95  //- Find edge among edgeLabels that uses v0 and v1
96  static label findEdge
97  (
98  const primitiveMesh& mesh,
99  const labelList& edgeLabels,
100  const label v1,
101  const label v0
102  );
103 
104  //- Return 'lowest' of a,b in face of size.
105  static label lowest
106  (
107  const label size,
108  const label a,
109  const label b
110  );
111 
112 public:
113 
114  // Static Functions
115 
116  //- Given edge on hex cell find corresponding edge on face. Is either
117  // index in face or -1 (cut tangential to face). Public since is
118  // needed to fill in seed faces in meshWave.
119  static label edgeToFaceIndex
120  (
121  const primitiveMesh& mesh,
122  const label cellI,
123  const label faceI,
124  const label edgeI
125  );
126 
127  // Constructors
128 
129  //- Construct null
130  inline directionInfo();
131 
132  //- Construct from components
133  inline directionInfo
134  (
135  const label,
136  const vector& n
137  );
138 
139  //- Construct as copy
140  inline directionInfo(const directionInfo&);
141 
142 
143  // Member Functions
144 
145  // Access
146 
147  inline label index() const
148  {
149  return index_;
150  }
151 
152  inline const vector& n() const
153  {
154  return n_;
155  }
156 
157  // Needed by FaceCellWave
158 
159  //- Check whether origin has been changed at all or
160  // still contains original (invalid) value.
161  inline bool valid() const;
162 
163  //- Check for identical geometrical data. Used for cyclics checking.
164  inline bool sameGeometry
165  (
166  const polyMesh&,
167  const directionInfo&,
168  const scalar
169  ) const;
170 
171  //- Convert any absolute coordinates into relative to (patch)face
172  // centre
173  inline void leaveDomain
174  (
175  const polyMesh&,
176  const polyPatch&,
177  const label patchFaceI,
178  const point& faceCentre
179  );
180 
181  //- Reverse of leaveDomain
182  inline void enterDomain
183  (
184  const polyMesh&,
185  const polyPatch&,
186  const label patchFaceI,
187  const point& faceCentre
188  );
189 
190  //- Apply rotation matrix to any coordinates
191  inline void transform
192  (
193  const polyMesh&,
194  const tensor&
195  );
196 
197  //- Influence of neighbouring face.
198  bool updateCell
199  (
200  const polyMesh&,
201  const label thisCellI,
202  const label neighbourFaceI,
203  const directionInfo& neighbourInfo,
204  const scalar tol
205  );
206 
207  //- Influence of neighbouring cell.
208  bool updateFace
209  (
210  const polyMesh&,
211  const label thisFaceI,
212  const label neighbourCellI,
213  const directionInfo& neighbourInfo,
214  const scalar tol
215  );
216 
217  //- Influence of different value on same face.
218  bool updateFace
219  (
220  const polyMesh&,
221  const label thisFaceI,
222  const directionInfo& neighbourInfo,
223  const scalar tol
224  );
225 
226  // Member Operators
227 
228  // Needed for List IO
229  inline bool operator==(const directionInfo&) const;
230 
231  inline bool operator!=(const directionInfo&) const;
232 
233 
234  // IOstream Operators
235 
236  friend Ostream& operator<<(Ostream&, const directionInfo&);
238 };
239 
240 
241 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
242 
243 } // End namespace Foam
244 
245 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
246 
248 
249 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
250 
251 #endif
252 
253 // ************************ vim: set sw=4 sts=4 et: ************************ //