All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
polyMeshCheck.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-2024 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::meshCheck
26 
27 Description
28 
29 SourceFiles
30  polyMeshCheck.C
31 
32 \*---------------------------------------------------------------------------*/
33 
34 #ifndef polyMeshCheck_H
35 #define polyMeshCheck_H
36 
37 #include "polyMesh.H"
38 #include "primitiveMeshCheck.H"
39 
40 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
41 
42 namespace Foam
43 {
44 
45 /*---------------------------------------------------------------------------*\
46  Namespace meshCheck Declaration
47 \*---------------------------------------------------------------------------*/
48 
49 namespace meshCheck
50 {
51  //- Generate orthogonality field. (1 for fully orthogonal, < 1 for
52  // non-orthogonal)
53  tmp<scalarField> faceOrthogonality
54  (
55  const polyMesh& mesh,
56  const vectorField& fAreas,
57  const vectorField& cellCtrs
58  );
59 
60  //- Generate skewness field
61  tmp<scalarField> faceSkewness
62  (
63  const polyMesh& mesh,
64  const pointField& points,
65  const vectorField& fCtrs,
66  const vectorField& fAreas,
67  const vectorField& cellCtrs
68  );
69 
70  //- Generate interpolation factors field
71  tmp<scalarField> faceWeights
72  (
73  const polyMesh& mesh,
74  const vectorField& fCtrs,
75  const vectorField& fAreas,
76  const vectorField& cellCtrs
77  );
78 
79  //- Generate volume ratio field
80  tmp<scalarField> volRatio
81  (
82  const polyMesh& mesh,
83  const scalarField& vol
84  );
85 
86  //- Check face orthogonality
88  (
89  const polyMesh& mesh,
90  const scalar nonOrthThreshold,
91  const bool report = false,
92  labelHashSet* setPtr = nullptr
93  );
94 
95  //- Check face skewness
97  (
98  const polyMesh& mesh,
99  const scalar skewThreshold,
100  const bool report = false,
101  labelHashSet* setPtr = nullptr
102  );
103 
104  //- Check edge alignment for 1D/2D cases
105  bool checkEdgeAlignment
106  (
107  const polyMesh& mesh,
108  const bool report,
109  const Vector<label>& directions,
110  labelHashSet* setPtr
111  );
112 
114  (
115  const polyMesh& mesh,
116  const bool report,
117  labelHashSet* setPtr
118  );
119 
120  //- Check for face weights
121  bool checkFaceWeight
122  (
123  const polyMesh& mesh,
124  const bool report,
125  const scalar minWeight = 0.05,
126  labelHashSet* setPtr = nullptr
127  );
128 
129  //- Check for neighbouring cell volumes
130  bool checkVolRatio
131  (
132  const polyMesh& mesh,
133  const bool report,
134  const scalar minRatio = 0.01,
135  labelHashSet* setPtr = nullptr
136  );
137 
138  //- Check face orthogonality
140  (
141  const bool report,
142  const scalar orthWarn,
143  const polyMesh&,
144  const vectorField& cellCentres,
145  const vectorField& faceAreas,
146  const labelList& checkFaces,
147  const List<labelPair>& baffles,
148  labelHashSet* setPtr
149  );
150 
151  //- Check face pyramid volumes
152  bool checkFacePyramids
153  (
154  const bool report,
155  const scalar minPyrVol,
156  const polyMesh&,
157  const vectorField& cellCentres,
158  const pointField& p,
159  const labelList& checkFaces,
160  const List<labelPair>& baffles,
161  labelHashSet*
162  );
163 
164  //- Check face tetrahedra volumes
165  bool checkFaceTets
166  (
167  const bool report,
168  const scalar minPyrVol,
169  const polyMesh&,
170  const vectorField& cellCentres,
171  const vectorField& faceCentres,
172  const pointField& p,
173  const labelList& checkFaces,
174  const List<labelPair>& baffles,
175  labelHashSet*
176  );
177 
178  //- Check face skewness
179  bool checkFaceSkewness
180  (
181  const bool report,
182  const scalar internalSkew,
183  const scalar boundarySkew,
184  const polyMesh& mesh,
185  const pointField& points,
186  const vectorField& cellCentres,
187  const vectorField& faceCentres,
188  const vectorField& faceAreas,
189  const labelList& checkFaces,
190  const List<labelPair>& baffles,
191  labelHashSet* setPtr
192  );
193 
194  //- Check interpolation weights (0.5 for regular mesh)
195  bool checkFaceWeights
196  (
197  const bool report,
198  const scalar warnWeight,
199  const polyMesh& mesh,
200  const vectorField& cellCentres,
201  const vectorField& faceCentres,
202  const vectorField& faceAreas,
203  const labelList& checkFaces,
204  const List<labelPair>& baffles,
205  labelHashSet* setPtr
206  );
207 
208  //- Cell volume ratio of neighbouring cells (1 for regular mesh)
209  bool checkVolRatio
210  (
211  const bool report,
212  const scalar warnRatio,
213  const polyMesh& mesh,
214  const scalarField& cellVolumes,
215  const labelList& checkFaces,
216  const List<labelPair>& baffles,
217  labelHashSet* setPtr
218  );
219 
220  //- Check convexity of angles in a face. See primitiveMesh for explanation.
221  bool checkFaceAngles
222  (
223  const bool report,
224  const scalar maxConcave,
225  const polyMesh& mesh,
226  const vectorField& faceAreas,
227  const pointField& p,
228  const labelList& checkFaces,
229  labelHashSet* setPtr
230  );
231 
232  //- Check the difference between normals of individual face-triangles (from
233  // face-centre decomposition) and the cell-cell centre vector
234  bool checkFaceTwist
235  (
236  const bool report,
237  const scalar minTwist,
238  const polyMesh&,
239  const vectorField& cellCentres,
240  const vectorField& faceAreas,
241  const vectorField& faceCentres,
242  const pointField& p,
243  const labelList& checkFaces,
244  labelHashSet* setPtr
245  );
246 
247  //- Check consecutive face-triangle (from face-centre decomposition) normals
248  bool checkTriangleTwist
249  (
250  const bool report,
251  const scalar minTwist,
252  const polyMesh&,
253  const vectorField& faceAreas,
254  const vectorField& faceCentres,
255  const pointField& p,
256  const labelList& checkFaces,
257  labelHashSet* setPtr
258  );
259 
260  //- Check for face areas v.s. sum of face-triangle (from face-centre
261  // decomposition) areas
262  bool checkFaceFlatness
263  (
264  const bool report,
265  const scalar minFlatness,
266  const polyMesh&,
267  const vectorField& faceAreas,
268  const vectorField& faceCentres,
269  const pointField& p,
270  const labelList& checkFaces,
271  labelHashSet* setPtr
272  );
273 
274  //- Check for small faces
275  bool checkFaceArea
276  (
277  const bool report,
278  const scalar minArea,
279  const polyMesh&,
280  const vectorField& faceAreas,
281  const labelList& checkFaces,
282  labelHashSet* setPtr
283  );
284 
285  //- Check the area of internal faces v.s. boundary faces
287  (
288  const bool report,
289  const scalar minDet,
290  const polyMesh&,
291  const vectorField& faceAreas,
292  const labelList& checkFaces,
293  labelHashSet* setPtr
294  );
295 }
296 
297 
298 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
299 
300 } // End namespace Foam
301 
302 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
303 
304 #endif
305 
306 // ************************************************************************* //
const pointField & points
bool checkFaceWeight(const polyMesh &mesh, const bool report, const scalar minWeight=0.05, labelHashSet *setPtr=nullptr)
Check for face weights.
bool checkFacePyramids(const bool report, const scalar minPyrVol, const polyMesh &, const vectorField &cellCentres, const pointField &p, const labelList &checkFaces, const List< labelPair > &baffles, labelHashSet *)
Check face pyramid volumes.
tmp< scalarField > faceOrthogonality(const polyMesh &mesh, const vectorField &fAreas, const vectorField &cellCtrs)
Generate orthogonality field. (1 for fully orthogonal, < 1 for.
Definition: polyMeshCheck.C:34
bool checkEdgeAlignment(const polyMesh &mesh, const bool report, const Vector< label > &directions, labelHashSet *setPtr)
Check edge alignment for 1D/2D cases.
tmp< scalarField > faceSkewness(const polyMesh &mesh, const pointField &points, const vectorField &fCtrs, const vectorField &fAreas, const vectorField &cellCtrs)
Generate skewness field.
Definition: polyMeshCheck.C:89
bool checkTriangleTwist(const bool report, const scalar minTwist, const polyMesh &, const vectorField &faceAreas, const vectorField &faceCentres, const pointField &p, const labelList &checkFaces, labelHashSet *setPtr)
Check consecutive face-triangle (from face-centre decomposition) normals.
bool checkFaceTwist(const bool report, const scalar minTwist, const polyMesh &, const vectorField &cellCentres, const vectorField &faceAreas, const vectorField &faceCentres, const pointField &p, const labelList &checkFaces, labelHashSet *setPtr)
Check the difference between normals of individual face-triangles (from.
bool checkFaceArea(const bool report, const scalar minArea, const polyMesh &, const vectorField &faceAreas, const labelList &checkFaces, labelHashSet *setPtr)
Check for small faces.
bool checkFaceWeights(const bool report, const scalar warnWeight, const polyMesh &mesh, const vectorField &cellCentres, const vectorField &faceCentres, const vectorField &faceAreas, const labelList &checkFaces, const List< labelPair > &baffles, labelHashSet *setPtr)
Check interpolation weights (0.5 for regular mesh)
bool checkCellDeterminant(const polyMesh &mesh, const bool report, labelHashSet *setPtr)
bool checkFaceAngles(const bool report, const scalar maxConcave, const polyMesh &mesh, const vectorField &faceAreas, const pointField &p, const labelList &checkFaces, labelHashSet *setPtr)
Check convexity of angles in a face. See primitiveMesh for explanation.
bool checkVolRatio(const polyMesh &mesh, const bool report, const scalar minRatio=0.01, labelHashSet *setPtr=nullptr)
Check for neighbouring cell volumes.
tmp< scalarField > faceWeights(const polyMesh &mesh, const vectorField &fCtrs, const vectorField &fAreas, const vectorField &cellCtrs)
Generate interpolation factors field.
tmp< scalarField > volRatio(const polyMesh &mesh, const scalarField &vol)
Generate volume ratio field.
bool checkFaceTets(const bool report, const scalar minPyrVol, const polyMesh &, const vectorField &cellCentres, const vectorField &faceCentres, const pointField &p, const labelList &checkFaces, const List< labelPair > &baffles, labelHashSet *)
Check face tetrahedra volumes.
bool checkFaceFlatness(const bool report, const scalar minFlatness, const polyMesh &, const vectorField &faceAreas, const vectorField &faceCentres, const pointField &p, const labelList &checkFaces, labelHashSet *setPtr)
Check for face areas v.s. sum of face-triangle (from face-centre.
bool checkFaceOrthogonality(const polyMesh &mesh, const scalar nonOrthThreshold, const bool report=false, labelHashSet *setPtr=nullptr)
Check face orthogonality.
bool checkFaceSkewness(const polyMesh &mesh, const scalar skewThreshold, const bool report=false, labelHashSet *setPtr=nullptr)
Check face skewness.
Namespace for OpenFOAM.
List< label > labelList
A List of labels.
Definition: labelList.H:56
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Field< vector > vectorField
Specialisation of Field<T> for vector.
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys.
Definition: HashSet.H:213
volScalarField & p