motionSmootherAlgoCheck.C
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-2022 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 \*---------------------------------------------------------------------------*/
25 
26 #include "motionSmootherAlgo.H"
27 #include "polyMeshCheck.H"
28 #include "IOmanip.H"
29 
30 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
31 
33 (
34  const bool report,
35  const polyMesh& mesh,
36  const dictionary& dict,
37  const labelList& checkFaces,
38  labelHashSet& wrongFaces
39 )
40 {
41  List<labelPair> emptyBaffles;
42  return checkMesh
43  (
44  report,
45  mesh,
46  dict,
47  checkFaces,
48  emptyBaffles,
49  wrongFaces
50  );
51 }
52 
54 (
55  const bool report,
56  const polyMesh& mesh,
57  const dictionary& dict,
58  const labelList& checkFaces,
59  const List<labelPair>& baffles,
60  labelHashSet& wrongFaces
61 )
62 {
63  const scalar maxNonOrtho
64  (
65  dict.lookup<scalar>("maxNonOrtho", true)
66  );
67  const scalar minVol
68  (
69  dict.lookup<scalar>("minVol", true)
70  );
71  const scalar minTetQuality
72  (
73  dict.lookup<scalar>("minTetQuality", true)
74  );
75  const scalar maxConcave
76  (
77  dict.lookup<scalar>("maxConcave", true)
78  );
79  const scalar maxIntSkew
80  (
81  dict.lookup<scalar>("maxInternalSkewness", true)
82  );
83  const scalar maxBounSkew
84  (
85  dict.lookup<scalar>("maxBoundarySkewness", true)
86  );
87  const scalar minWeight
88  (
89  dict.lookup<scalar>("minFaceWeight", true)
90  );
91  const scalar minVolRatio
92  (
93  dict.lookup<scalar>("minVolRatio", true)
94  );
95  const scalar minTwist
96  (
97  dict.lookup<scalar>("minTwist", true)
98  );
99  scalar minFaceFlatness = -1.0;
100  dict.readIfPresent("minFaceFlatness", minFaceFlatness, true);
101  const scalar minDet
102  (
103  dict.lookup<scalar>("minDeterminant", true)
104  );
105  label nWrongFaces = 0;
106 
107  Info<< "Checking faces in error :" << endl;
108  // Pout.setf(ios_base::left);
109 
110  if (maxNonOrtho < 180.0-small)
111  {
113  (
114  report,
115  maxNonOrtho,
116  mesh,
117  mesh.cellCentres(),
118  mesh.faceAreas(),
119  checkFaces,
120  baffles,
121  &wrongFaces
122  );
123 
124  label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
125 
126  Info<< " non-orthogonality > "
127  << setw(3) << maxNonOrtho
128  << " degrees : "
129  << nNewWrongFaces-nWrongFaces << endl;
130 
131  nWrongFaces = nNewWrongFaces;
132  }
133 
134  if (minVol > -great)
135  {
136  const scalar refVol = pow3(mesh.bounds().minDim());
137 
139  (
140  report,
141  minVol*refVol,
142  mesh,
143  mesh.cellCentres(),
144  mesh.points(),
145  checkFaces,
146  baffles,
147  &wrongFaces
148  );
149 
150  label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
151 
152  Info<< " faces with face pyramid volume < "
153  << setw(5) << minVol << " : "
154  << nNewWrongFaces-nWrongFaces << endl;
155 
156  nWrongFaces = nNewWrongFaces;
157  }
158 
159  if (minTetQuality > -great)
160  {
162  (
163  report,
164  minTetQuality,
165  mesh,
166  mesh.cellCentres(),
167  mesh.faceCentres(),
168  mesh.points(),
169  checkFaces,
170  baffles,
171  &wrongFaces
172  );
173 
174  label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
175 
176  Info<< " faces with face-decomposition tet quality < "
177  << setw(5) << minTetQuality << " : "
178  << nNewWrongFaces-nWrongFaces << endl;
179 
180  nWrongFaces = nNewWrongFaces;
181  }
182 
183  if (maxConcave < 180.0-small)
184  {
186  (
187  report,
188  maxConcave,
189  mesh,
190  mesh.faceAreas(),
191  mesh.points(),
192  checkFaces,
193  &wrongFaces
194  );
195 
196  label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
197 
198  Info<< " faces with concavity > "
199  << setw(3) << maxConcave
200  << " degrees : "
201  << nNewWrongFaces-nWrongFaces << endl;
202 
203  nWrongFaces = nNewWrongFaces;
204  }
205 
206  if (maxIntSkew > 0 || maxBounSkew > 0)
207  {
209  (
210  report,
211  maxIntSkew,
212  maxBounSkew,
213  mesh,
214  mesh.points(),
215  mesh.cellCentres(),
216  mesh.faceCentres(),
217  mesh.faceAreas(),
218  checkFaces,
219  baffles,
220  &wrongFaces
221  );
222 
223  label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
224 
225  Info<< " faces with skewness > "
226  << setw(3) << maxIntSkew
227  << " (internal) or " << setw(3) << maxBounSkew
228  << " (boundary) : " << nNewWrongFaces-nWrongFaces << endl;
229 
230  nWrongFaces = nNewWrongFaces;
231  }
232 
233  if (minWeight >= 0 && minWeight < 1)
234  {
236  (
237  report,
238  minWeight,
239  mesh,
240  mesh.cellCentres(),
241  mesh.faceCentres(),
242  mesh.faceAreas(),
243  checkFaces,
244  baffles,
245  &wrongFaces
246  );
247 
248  label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
249 
250  Info<< " faces with interpolation weights (0..1) < "
251  << setw(5) << minWeight
252  << " : "
253  << nNewWrongFaces-nWrongFaces << endl;
254 
255  nWrongFaces = nNewWrongFaces;
256  }
257 
258  if (minVolRatio >= 0)
259  {
261  (
262  report,
263  minVolRatio,
264  mesh,
265  mesh.cellVolumes(),
266  checkFaces,
267  baffles,
268  &wrongFaces
269  );
270 
271  label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
272 
273  Info<< " faces with volume ratio of neighbour cells < "
274  << setw(5) << minVolRatio
275  << " : "
276  << nNewWrongFaces-nWrongFaces << endl;
277 
278  nWrongFaces = nNewWrongFaces;
279  }
280 
281  if (minTwist > -1)
282  {
283  // Pout<< "Checking face twist: dot product of face normal "
284  // << "with face triangle normals" << endl;
286  (
287  report,
288  minTwist,
289  mesh,
290  mesh.cellCentres(),
291  mesh.faceAreas(),
292  mesh.faceCentres(),
293  mesh.points(),
294  checkFaces,
295  &wrongFaces
296  );
297 
298  label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
299 
300  Info<< " faces with face twist < "
301  << setw(5) << minTwist
302  << " : "
303  << nNewWrongFaces-nWrongFaces << endl;
304 
305  nWrongFaces = nNewWrongFaces;
306  }
307 
308  if (minFaceFlatness > -small)
309  {
311  (
312  report,
313  minFaceFlatness,
314  mesh,
315  mesh.faceAreas(),
316  mesh.faceCentres(),
317  mesh.points(),
318  checkFaces,
319  &wrongFaces
320  );
321 
322  label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
323 
324  Info<< " faces with flatness < "
325  << setw(5) << minFaceFlatness
326  << " : "
327  << nNewWrongFaces-nWrongFaces << endl;
328 
329  nWrongFaces = nNewWrongFaces;
330  }
331 
332  if (minDet > -1)
333  {
335  (
336  report,
337  minDet,
338  mesh,
339  mesh.faceAreas(),
340  checkFaces,
341  &wrongFaces
342  );
343 
344  label nNewWrongFaces = returnReduce(wrongFaces.size(), sumOp<label>());
345 
346  Info<< " faces on cells with determinant < "
347  << setw(5) << minDet << " : "
348  << nNewWrongFaces-nWrongFaces << endl;
349 
350  nWrongFaces = nNewWrongFaces;
351  }
352 
353  // Pout.setf(ios_base::right);
354 
355  return nWrongFaces > 0;
356 }
357 
358 
360 (
361  const bool report,
362  const polyMesh& mesh,
363  const dictionary& dict,
364  labelHashSet& wrongFaces
365 )
366 {
367  return checkMesh
368  (
369  report,
370  mesh,
371  dict,
372  identity(mesh.nFaces()),
373  wrongFaces
374  );
375 }
376 
377 
378 // ************************************************************************* //
static bool checkMesh(const bool report, const polyMesh &mesh, const dictionary &dict, labelHashSet &wrongFaces)
Check mesh with mesh settings in dict. Collects incorrect faces.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
label nFaces() const
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
Definition: ListOps.C:104
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)
Check face non-orthogonality.
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)
Check face skewness.
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1211
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.
PolyMesh checking routines. Checks various criteria for a mesh and supplied geometry, with the mesh only used for topology. Coupled patch aware (i.e., coupled faces are treated as internal).
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
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)
const vectorField & cellCentres() const
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.
Istream and Ostream manipulators taking arguments.
bool checkFaceAngles(const bool report, const scalar maxDeg, 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.
const vectorField & faceCentres() const
dimensionedScalar pow3(const dimensionedScalar &ds)
const boundBox & bounds() const
Return mesh bounding box.
Definition: polyMesh.H:459
bool checkCellDeterminant(const bool report, const scalar minDet, const polyMesh &, const vectorField &faceAreas, const labelList &checkFaces, labelHashSet *setPtr)
Check the area of internal faces v.s. boundary faces.
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.
const vectorField & faceAreas() const
messageStream Info
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
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)
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:76
Omanip< int > setw(const int i)
Definition: IOmanip.H:199
scalar minDim() const
Smallest length/height/width dimension.
Definition: boundBoxI.H:102
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)
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:864
const scalarField & cellVolumes() const