polyMeshGeometry.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-2018 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::polyMeshGeometry
26 
27 Description
28  Updateable mesh geometry and checking routines.
29 
30  - non-ortho done across coupled faces.
31  - faceWeight (delta factors) done across coupled faces.
32 
33 SourceFiles
34  polyMeshGeometry.C
35 
36 \*---------------------------------------------------------------------------*/
37 
38 #ifndef polyMeshGeometry_H
39 #define polyMeshGeometry_H
40 
41 #include "pointFields.H"
42 #include "HashSet.H"
43 
44 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
45 
46 namespace Foam
47 {
48 
49 /*---------------------------------------------------------------------------*\
50  Class polyMeshGeometry Declaration
51 \*---------------------------------------------------------------------------*/
52 
53 class polyMeshGeometry
54 {
55  //- Reference to polyMesh.
56  const polyMesh& mesh_;
57 
58  //- Uptodate copy of face areas
59  vectorField faceAreas_;
60 
61  //- Uptodate copy of face centres
62  vectorField faceCentres_;
63 
64  //- Uptodate copy of cell centres
65  vectorField cellCentres_;
66 
67  //- Uptodate copy of cell volumes
68  scalarField cellVolumes_;
69 
70 
71  // Private Member Functions
72 
73  //- Update face areas and centres on selected faces.
74  void updateFaceCentresAndAreas
75  (
76  const pointField& p,
77  const labelList& changedFaces
78  );
79 
80  //- Update cell volumes and centres on selected cells. Requires
81  // cells and faces to be consistent set.
82  void updateCellCentresAndVols
83  (
84  const labelList& changedCells,
85  const labelList& changedFaces
86  );
87 
88  //- Detect&report non-ortho error for single face.
89  static scalar checkNonOrtho
90  (
91  const polyMesh& mesh,
92  const bool report,
93  const scalar severeNonorthogonalityThreshold,
94  const label facei,
95  const vector& s, // face area vector
96  const vector& d, // cc-cc vector
97 
98  label& severeNonOrth,
99  label& errorNonOrth,
100  labelHashSet* setPtr
101  );
102 
103  //- Calculate skewness given two cell centres and one face centre.
104  static scalar calcSkewness
105  (
106  const point& ownCc,
107  const point& neiCc,
108  const point& fc
109  );
110 
111  //- Detect&report incorrect face-triangle orientation
112  static bool checkFaceTet
113  (
114  const polyMesh&,
115  const bool report,
116  const scalar minTetQuality,
117  const pointField& p,
118  const label facei,
119  const point& fc, // face centre
120  const point& cc, // cell centre
121  labelHashSet* setPtr
122  );
123 
124 
125 public:
126 
127  ClassName("polyMeshGeometry");
128 
129  // Constructors
130 
131  //- Construct from mesh
132  polyMeshGeometry(const polyMesh&);
133 
134 
135  // Member Functions
136 
137  // Access
139  const polyMesh& mesh() const
140  {
141  return mesh_;
142  }
144  const vectorField& faceAreas() const
145  {
146  return faceAreas_;
147  }
148  const vectorField& faceCentres() const
149  {
150  return faceCentres_;
151  }
152  const vectorField& cellCentres() const
153  {
154  return cellCentres_;
155  }
156  const scalarField& cellVolumes() const
157  {
158  return cellVolumes_;
159  }
160 
161  // Edit
162 
163  //- Take over properties from mesh
164  void correct();
165 
166  //- Recalculate on selected faces. Recalculates cell properties
167  // on owner and neighbour of these cells.
168  void correct
169  (
170  const pointField& p,
171  const labelList& changedFaces
172  );
173 
174  //- Helper function: get affected cells from faces
175  static labelList affectedCells
176  (
177  const polyMesh&,
178  const labelList& changedFaces
179  );
180 
181 
182  // Checking of selected faces with supplied geometry (mesh only used for
183  // topology). Coupled aware (coupled faces treated as internal ones)
184 
185  //- See primitiveMesh
186  static bool checkFaceDotProduct
187  (
188  const bool report,
189  const scalar orthWarn,
190  const polyMesh&,
191  const vectorField& cellCentres,
192  const vectorField& faceAreas,
193  const labelList& checkFaces,
194  const List<labelPair>& baffles,
195  labelHashSet* setPtr
196  );
197 
198  //- See primitiveMesh
199  static bool checkFacePyramids
200  (
201  const bool report,
202  const scalar minPyrVol,
203  const polyMesh&,
204  const vectorField& cellCentres,
205  const pointField& p,
206  const labelList& checkFaces,
207  const List<labelPair>& baffles,
208  labelHashSet*
209  );
210 
211  //- See primitiveMesh
212  static bool checkFaceTets
213  (
214  const bool report,
215  const scalar minPyrVol,
216  const polyMesh&,
217  const vectorField& cellCentres,
218  const vectorField& faceCentres,
219  const pointField& p,
220  const labelList& checkFaces,
221  const List<labelPair>& baffles,
222  labelHashSet*
223  );
224 
225  //- See primitiveMesh
226  static bool checkFaceSkewness
227  (
228  const bool report,
229  const scalar internalSkew,
230  const scalar boundarySkew,
231  const polyMesh& mesh,
232  const pointField& points,
233  const vectorField& cellCentres,
234  const vectorField& faceCentres,
235  const vectorField& faceAreas,
236  const labelList& checkFaces,
237  const List<labelPair>& baffles,
238  labelHashSet* setPtr
239  );
240 
241  //- Interpolation weights (0.5 for regular mesh)
242  static bool checkFaceWeights
243  (
244  const bool report,
245  const scalar warnWeight,
246  const polyMesh& mesh,
247  const vectorField& cellCentres,
248  const vectorField& faceCentres,
249  const vectorField& faceAreas,
250  const labelList& checkFaces,
251  const List<labelPair>& baffles,
252  labelHashSet* setPtr
253  );
254 
255  //- Cell volume ratio of neighbouring cells (1 for regular mesh)
256  static bool checkVolRatio
257  (
258  const bool report,
259  const scalar warnRatio,
260  const polyMesh& mesh,
261  const scalarField& cellVolumes,
262  const labelList& checkFaces,
263  const List<labelPair>& baffles,
264  labelHashSet* setPtr
265  );
266 
267  //- See primitiveMesh
268  static bool checkFaceAngles
269  (
270  const bool report,
271  const scalar maxDeg,
272  const polyMesh& mesh,
273  const vectorField& faceAreas,
274  const pointField& p,
275  const labelList& checkFaces,
276  labelHashSet* setPtr
277  );
278 
279  //- Triangle (from face-centre decomposition) normal v.s.
280  // average face normal
281  static bool checkFaceTwist
282  (
283  const bool report,
284  const scalar minTwist,
285  const polyMesh&,
286  const vectorField& cellCentres,
287  const vectorField& faceAreas,
288  const vectorField& faceCentres,
289  const pointField& p,
290  const labelList& checkFaces,
291  labelHashSet* setPtr
292  );
293 
294  //- Consecutive triangle (from face-centre decomposition) normals
295  static bool checkTriangleTwist
296  (
297  const bool report,
298  const scalar minTwist,
299  const polyMesh&,
300  const vectorField& faceAreas,
301  const vectorField& faceCentres,
302  const pointField& p,
303  const labelList& checkFaces,
304  labelHashSet* setPtr
305  );
306 
307  //- Area of faces v.s. sum of triangle areas
308  static bool checkFaceFlatness
309  (
310  const bool report,
311  const scalar minFlatness,
312  const polyMesh&,
313  const vectorField& faceAreas,
314  const vectorField& faceCentres,
315  const pointField& p,
316  const labelList& checkFaces,
317  labelHashSet* setPtr
318  );
319 
320  //- Small faces
321  static bool checkFaceArea
322  (
323  const bool report,
324  const scalar minArea,
325  const polyMesh&,
326  const vectorField& faceAreas,
327  const labelList& checkFaces,
328  labelHashSet* setPtr
329  );
330 
331  //- Area of internal faces v.s. boundary faces
332  static bool checkCellDeterminant
333  (
334  const bool report,
335  const scalar minDet,
336  const polyMesh&,
337  const vectorField& faceAreas,
338  const labelList& checkFaces,
339  const labelList& affectedCells,
340  labelHashSet* setPtr
341  );
342 
343 
344  // Checking of selected faces with local geometry. Uses above static
345  // functions. Parallel aware.
346 
348  (
349  const bool report,
350  const scalar orthWarn,
351  const labelList& checkFaces,
352  const List<labelPair>& baffles,
353  labelHashSet* setPtr
354  ) const;
355 
356  bool checkFacePyramids
357  (
358  const bool report,
359  const scalar minPyrVol,
360  const pointField& p,
361  const labelList& checkFaces,
362  const List<labelPair>& baffles,
363  labelHashSet* setPtr
364  ) const;
365 
366  bool checkFaceTets
367  (
368  const bool report,
369  const scalar minTetQuality,
370  const pointField& p,
371  const labelList& checkFaces,
372  const List<labelPair>& baffles,
373  labelHashSet* setPtr
374  ) const;
375 
376  bool checkFaceWeights
377  (
378  const bool report,
379  const scalar warnWeight,
380  const labelList& checkFaces,
381  const List<labelPair>& baffles,
382  labelHashSet* setPtr
383  ) const;
384 
385  bool checkVolRatio
386  (
387  const bool report,
388  const scalar warnRatio,
389  const labelList& checkFaces,
390  const List<labelPair>& baffles,
391  labelHashSet* setPtr
392  ) const;
393 
394  bool checkFaceAngles
395  (
396  const bool report,
397  const scalar maxDeg,
398  const pointField& p,
399  const labelList& checkFaces,
400  labelHashSet* setPtr
401  ) const;
402 
403  bool checkFaceTwist
404  (
405  const bool report,
406  const scalar minTwist,
407  const pointField& p,
408  const labelList& checkFaces,
409  labelHashSet* setPtr
410  ) const;
411 
412  bool checkTriangleTwist
413  (
414  const bool report,
415  const scalar minTwist,
416  const pointField& p,
417  const labelList& checkFaces,
418  labelHashSet* setPtr
419  ) const;
420 
421  bool checkFaceFlatness
422  (
423  const bool report,
424  const scalar minFlatness,
425  const pointField& p,
426  const labelList& checkFaces,
427  labelHashSet* setPtr
428  ) const;
429 
430  bool checkFaceArea
431  (
432  const bool report,
433  const scalar minArea,
434  const labelList& checkFaces,
435  labelHashSet* setPtr
436  ) const;
437 
439  (
440  const bool report,
441  const scalar warnDet,
442  const labelList& checkFaces,
443  const labelList& affectedCells,
444  labelHashSet* setPtr
445  ) const;
446 };
447 
448 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
449 
450 } // End namespace Foam
451 
452 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
453 
454 #endif
455 
456 // ************************************************************************* //
static bool checkFacePyramids(const bool report, const scalar minPyrVol, const polyMesh &, const vectorField &cellCentres, const pointField &p, const labelList &checkFaces, const List< labelPair > &baffles, labelHashSet *)
See primitiveMesh.
static 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 *)
See primitiveMesh.
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
static bool checkTriangleTwist(const bool report, const scalar minTwist, const polyMesh &, const vectorField &faceAreas, const vectorField &faceCentres, const pointField &p, const labelList &checkFaces, labelHashSet *setPtr)
Consecutive triangle (from face-centre decomposition) normals.
const vectorField & faceAreas() const
static 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)
Triangle (from face-centre decomposition) normal v.s.
const vectorField & cellCentres() const
static bool checkVolRatio(const bool report, const scalar warnRatio, const polyMesh &mesh, const scalarField &cellVolumes, const labelList &checkFaces, const List< labelPair > &baffles, labelHashSet *setPtr)
Cell volume ratio of neighbouring cells (1 for regular mesh)
static bool checkFaceFlatness(const bool report, const scalar minFlatness, const polyMesh &, const vectorField &faceAreas, const vectorField &faceCentres, const pointField &p, const labelList &checkFaces, labelHashSet *setPtr)
Area of faces v.s. sum of triangle areas.
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
const pointField & points
Updateable mesh geometry and checking routines.
const scalarField & cellVolumes() const
static bool checkFaceArea(const bool report, const scalar minArea, const polyMesh &, const vectorField &faceAreas, const labelList &checkFaces, labelHashSet *setPtr)
Small faces.
const vectorField & faceCentres() const
static bool checkFaceAngles(const bool report, const scalar maxDeg, const polyMesh &mesh, const vectorField &faceAreas, const pointField &p, const labelList &checkFaces, labelHashSet *setPtr)
See primitiveMesh.
static bool checkFaceDotProduct(const bool report, const scalar orthWarn, const polyMesh &, const vectorField &cellCentres, const vectorField &faceAreas, const labelList &checkFaces, const List< labelPair > &baffles, labelHashSet *setPtr)
See primitiveMesh.
void correct()
Take over properties from mesh.
static labelList affectedCells(const polyMesh &, const labelList &changedFaces)
Helper function: get affected cells from faces.
static 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)
Interpolation weights (0.5 for regular mesh)
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
static bool checkFaceSkewness(const bool report, const scalar internalSkew, const scalar boundarySkew, const polyMesh &mesh, const pointField &points, const vectorField &cellCentres, const vectorField &faceCentres, const vectorField &faceAreas, const labelList &checkFaces, const List< labelPair > &baffles, labelHashSet *setPtr)
See primitiveMesh.
volScalarField & p
static bool checkCellDeterminant(const bool report, const scalar minDet, const polyMesh &, const vectorField &faceAreas, const labelList &checkFaces, const labelList &affectedCells, labelHashSet *setPtr)
Area of internal faces v.s. boundary faces.
const polyMesh & mesh() const
ClassName("polyMeshGeometry")
Namespace for OpenFOAM.
polyMeshGeometry(const polyMesh &)
Construct from mesh.