meshTools.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 Namespace
25  Foam::meshTools
26 
27 Description
28  Collection of static functions to do various simple mesh related things.
29 
30 SourceFiles
31  meshTools.C
32 
33 \*---------------------------------------------------------------------------*/
34 
35 #ifndef meshTools_H
36 #define meshTools_H
37 
38 #include "label.H"
39 #include "vector.H"
40 #include "triad.H"
41 #include "labelList.H"
42 #include "pointField.H"
43 #include "faceList.H"
44 #include "cellList.H"
45 #include "primitivePatch.H"
46 
47 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
48 
49 namespace Foam
50 {
51 
52 class primitiveMesh;
53 class polyMesh;
54 
55 /*---------------------------------------------------------------------------*\
56  Namespace meshTools Declaration
57 \*---------------------------------------------------------------------------*/
58 
59 namespace meshTools
60 {
61  // Bit identifiers for octants (p=plus, m=min e.g. plusXminYminZ)
62 
63  static const label mXmYmZ = 0;
64  static const label pXmYmZ = 1;
65  static const label mXpYmZ = 2;
66  static const label pXpYmZ = 3;
67 
68  static const label mXmYpZ = 4;
69  static const label pXmYpZ = 5;
70  static const label mXpYpZ = 6;
71  static const label pXpYpZ = 7;
72 
73  static const label mXmYmZMask = 1 << mXmYmZ;
74  static const label pXmYmZMask = 1 << pXmYmZ;
75  static const label mXpYmZMask = 1 << mXpYmZ;
76  static const label pXpYmZMask = 1 << pXpYmZ;
77 
78  static const label mXmYpZMask = 1 << mXmYpZ;
79  static const label pXmYpZMask = 1 << pXmYpZ;
80  static const label mXpYpZMask = 1 << mXpYpZ;
81  static const label pXpYpZMask = 1 << pXpYpZ;
82 
83 
84  // Normal handling
85 
86  //- Check if n is in same direction as normals of all faceLabels
87  bool visNormal
88  (
89  const vector& n,
90  const vectorField& faceNormals,
91  const labelList& faceLabels
92  );
93 
94  //- Calculate point normals on a 'box' mesh (all edges aligned with
95  // coordinate axes)
97 
98  //- Normalised edge vector
99  vector normEdgeVec(const primitiveMesh&, const label edgeI);
100 
101 
102  // OBJ writing
103 
104  //- Write obj representation of point
105  void writeOBJ
106  (
107  Ostream& os,
108  const point& pt
109  );
110 
111  //- Write obj representation of a triad. Requires the location of the
112  // triad to be supplied
113  void writeOBJ
114  (
115  Ostream& os,
116  const triad& t,
117  const point& pt
118  );
119 
120  //- Write obj representation of a line connecting two points
121  // Need to keep track of points that have been added. count starts at 0
122  void writeOBJ
123  (
124  Ostream& os,
125  const point& p1,
126  const point& p2,
127  label& count
128  );
129 
130  //- Write obj representation of a point p1 with a vector from p1 to p2
131  void writeOBJ
132  (
133  Ostream& os,
134  const point& p1,
135  const point& p2
136  );
137 
138  //- Write obj representation of faces subset
139  template<class FaceType>
140  void writeOBJ
141  (
142  Ostream& os,
143  const UList<FaceType>&,
144  const pointField&,
145  const labelList& faceLabels
146  );
147 
148  //- Write obj representation of faces
149  template<class FaceType>
150  void writeOBJ
151  (
152  Ostream& os,
153  const UList<FaceType>&,
154  const pointField&
155  );
156 
157  //- Write obj representation of cell subset
158  void writeOBJ
159  (
160  Ostream& os,
161  const cellList&,
162  const faceList&,
163  const pointField&,
164  const labelList& cellLabels
165  );
166 
167 
168  // Cell/face/edge walking
169 
170  //- Is edge used by cell
171  bool edgeOnCell
172  (
173  const primitiveMesh&,
174  const label celli,
175  const label edgeI
176  );
177 
178  //- Is edge used by face
179  bool edgeOnFace
180  (
181  const primitiveMesh&,
182  const label facei,
183  const label edgeI
184  );
185 
186  //- Is face used by cell
187  bool faceOnCell
188  (
189  const primitiveMesh&,
190  const label celli,
191  const label facei
192  );
193 
194  //- Return edge among candidates that uses the two vertices.
196  (
197  const edgeList& edges,
198  const labelList& candidates,
199  const label v0,
200  const label v1
201  );
202 
203  //- Return edge between two vertices. Returns -1 if no edge.
205  (
206  const primitiveMesh&,
207  const label v0,
208  const label v1
209  );
210 
211  //- Return edge shared by two faces. Throws error if no edge found.
213  (
214  const primitiveMesh&,
215  const label f0,
216  const label f1
217  );
218 
219  //- Return face shared by two cells. Throws error if none found.
221  (
222  const primitiveMesh&,
223  const label cell0,
224  const label cell1
225  );
226 
227  //- Get faces on cell using edgeI. Throws error if no two found.
228  void getEdgeFaces
229  (
230  const primitiveMesh&,
231  const label celli,
232  const label edgeI,
233  label& face0,
234  label& face1
235  );
236 
237  //- Return label of other edge (out of candidates edgeLabels)
238  // connected to vertex but not edgeI. Throws error if none found.
240  (
241  const primitiveMesh&,
242  const labelList& edgeLabels,
243  const label edgeI,
244  const label vertI
245  );
246 
247  //- Return face on cell using edgeI but not facei. Throws error
248  // if none found.
250  (
251  const primitiveMesh&,
252  const label celli,
253  const label facei,
254  const label edgeI
255  );
256 
257  //- Return cell on other side of face. Throws error
258  // if face not internal.
260  (
261  const primitiveMesh&,
262  const label celli,
263  const label facei
264  );
265 
266  //- Returns label of edge nEdges away from startEdge (in the direction
267  // of startVertI)
269  (
270  const primitiveMesh&,
271  const label facei,
272  const label startEdgeI,
273  const label startVertI,
274  const label nEdges
275  );
276 
277 
278  // Constraints on position
279 
280  //- Set the constrained components of position to mesh centre
282  (
283  const polyMesh& mesh,
284  point& pt
285  );
287  (
288  const polyMesh& mesh,
289  pointField& pt
290  );
291 
292  //- Set the constrained components of directions/velocity to zero
293  void constrainDirection
294  (
295  const polyMesh& mesh,
296  const Vector<label>& dirs,
297  vector& d
298  );
299  void constrainDirection
300  (
301  const polyMesh& mesh,
302  const Vector<label>& dirs,
303  vectorField& d
304  );
305 
306 
307  // Hex only functionality.
308 
309  //- Given edge on hex find other 'parallel', non-connected edges.
310  void getParallelEdges
311  (
312  const primitiveMesh&,
313  const label celli,
314  const label e0,
315  label&,
316  label&,
317  label&
318  );
319 
320  //- Given edge on hex find all 'parallel' (i.e. non-connected)
321  // edges and average direction of them
323  (
324  const primitiveMesh&,
325  const label celli,
326  const label edgeI
327  );
328 
329  //- Reverse of edgeToCutDir: given direction find edge bundle and
330  // return one of them.
332  (
333  const primitiveMesh&,
334  const label celli,
335  const vector& cutDir
336  );
337 
338 } // End namespace meshTools
339 
340 
341 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
342 
343 } // End namespace Foam
344 
345 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
346 
347 #ifdef NoRepository
348  #include "meshToolsTemplates.C"
349 #endif
350 
351 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
352 
353 #endif
354 
355 // ************************************************************************* //
label n
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: UList.H:74
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:80
Cell-face mesh analysis engine.
Definition: primitiveMesh.H:75
Representation of a 3D Cartesian coordinate system as a Vector of vectors.
Definition: triad.H:68
static const label mXmYpZMask
Definition: meshTools.H:78
void constrainDirection(const polyMesh &mesh, const Vector< label > &dirs, vector &d)
Set the constrained components of directions/velocity to zero.
Definition: meshTools.C:671
static const label mXpYmZ
Definition: meshTools.H:65
void constrainToMeshCentre(const polyMesh &mesh, point &pt)
Set the constrained components of position to mesh centre.
Definition: meshTools.C:613
static const label pXmYpZ
Definition: meshTools.H:69
void getEdgeFaces(const primitiveMesh &, const label celli, const label edgeI, label &face0, label &face1)
Get faces on cell using edgeI. Throws error if no two found.
Definition: meshTools.C:458
label getSharedEdge(const primitiveMesh &, const label f0, const label f1)
Return edge shared by two faces. Throws error if no edge found.
Definition: meshTools.C:386
static const label mXpYpZ
Definition: meshTools.H:70
vectorField calcBoxPointNormals(const primitivePatch &pp)
Calculate point normals on a 'box' mesh (all edges aligned with.
Definition: meshTools.C:52
static const label mXpYmZMask
Definition: meshTools.H:75
label findEdge(const edgeList &edges, const labelList &candidates, const label v0, const label v1)
Return edge among candidates that uses the two vertices.
Definition: meshTools.C:337
static const label pXpYmZ
Definition: meshTools.H:66
label walkFace(const primitiveMesh &, const label facei, const label startEdgeI, const label startVertI, const label nEdges)
Returns label of edge nEdges away from startEdge (in the direction.
Definition: meshTools.C:587
static const label pXpYpZ
Definition: meshTools.H:71
bool faceOnCell(const primitiveMesh &, const label celli, const label facei)
Is face used by cell.
Definition: meshTools.C:308
label otherCell(const primitiveMesh &, const label celli, const label facei)
Return cell on other side of face. Throws error.
Definition: meshTools.C:561
static const label pXpYmZMask
Definition: meshTools.H:76
static const label pXmYpZMask
Definition: meshTools.H:79
static const label pXmYmZ
Definition: meshTools.H:64
vector normEdgeVec(const primitiveMesh &, const label edgeI)
Normalised edge vector.
Definition: meshTools.C:189
label getSharedFace(const primitiveMesh &, const label cell0, const label cell1)
Return face shared by two cells. Throws error if none found.
Definition: meshTools.C:419
static const label pXmYmZMask
Definition: meshTools.H:74
bool visNormal(const vector &n, const vectorField &faceNormals, const labelList &faceLabels)
Check if n is in same direction as normals of all faceLabels.
Definition: meshTools.C:33
static const label mXmYmZ
Definition: meshTools.H:63
bool edgeOnCell(const primitiveMesh &, const label celli, const label edgeI)
Is edge used by cell.
Definition: meshTools.C:285
static const label mXmYmZMask
Definition: meshTools.H:73
label cutDirToEdge(const primitiveMesh &, const label celli, const vector &cutDir)
Reverse of edgeToCutDir: given direction find edge bundle and.
Definition: meshTools.C:795
static const label pXpYpZMask
Definition: meshTools.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
vector edgeToCutDir(const primitiveMesh &, const label celli, const label edgeI)
Given edge on hex find all 'parallel' (i.e. non-connected)
Definition: meshTools.C:747
void getParallelEdges(const primitiveMesh &, const label celli, const label e0, label &, label &, label &)
Given edge on hex find other 'parallel', non-connected edges.
Definition: meshTools.C:721
static const label mXpYpZMask
Definition: meshTools.H:80
bool edgeOnFace(const primitiveMesh &, const label facei, const label edgeI)
Is edge used by face.
Definition: meshTools.C:296
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
static const label mXmYpZ
Definition: meshTools.H:68
void writeOBJ(Ostream &os, const point &pt)
Write obj representation of point.
Definition: meshTools.C:203
Namespace for OpenFOAM.
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 count(const ListType &l, typename ListType::const_reference x)
Count the number of occurrences of a value in a list.