Home
Downloads
Documentation
Installation
User Guide
man-pages
API Documentation
README
Release Notes
Changes
License
Support
SourceForge Project
Main Page
Related Pages
Namespaces
Classes
Files
File List
File Members
All
Classes
Namespaces
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
Friends
Macros
Pages
src
dynamicMesh
motionSmoother
polyMeshGeometry
polyMeshGeometry.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::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 <
OpenFOAM/pointFields.H
>
42
#include <
OpenFOAM/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
public
:
112
113
ClassName
(
"polyMeshGeometry"
);
114
115
// Constructors
116
117
//- Construct from mesh
118
polyMeshGeometry
(
const
polyMesh
&);
119
120
121
// Member Functions
122
123
// Access
124
125
const
polyMesh
&
mesh
()
const
126
{
127
return
mesh_;
128
}
129
130
const
vectorField
&
faceAreas
()
const
131
{
132
return
faceAreas_;
133
}
134
const
vectorField
&
faceCentres
()
const
135
{
136
return
faceCentres_;
137
}
138
const
vectorField
&
cellCentres
()
const
139
{
140
return
cellCentres_;
141
}
142
const
scalarField
&
cellVolumes
()
const
143
{
144
return
cellVolumes_;
145
}
146
147
// Edit
148
149
//- Take over properties from mesh
150
void
correct
();
151
152
//- Recalculate on selected faces. Recalculates cell properties
153
// on owner and neighbour of these cells.
154
void
correct
155
(
156
const
pointField
& p,
157
const
labelList
& changedFaces
158
);
159
160
//- Helper function: get affected cells from faces
161
static
labelList
affectedCells
162
(
163
const
polyMesh
&,
164
const
labelList
& changedFaces
165
);
166
167
168
// Checking of selected faces with supplied geometry (mesh only used for
169
// topology). Coupled aware (coupled faces treated as internal ones)
170
171
//- See primitiveMesh
172
static
bool
checkFaceDotProduct
173
(
174
const
bool
report,
175
const
scalar orthWarn,
176
const
polyMesh
&,
177
const
vectorField
&
cellCentres
,
178
const
vectorField
&
faceAreas
,
179
const
labelList
& checkFaces,
180
const
List<labelPair>
& baffles,
181
labelHashSet
* setPtr
182
);
183
184
//- See primitiveMesh
185
static
bool
checkFacePyramids
186
(
187
const
bool
report,
188
const
scalar minPyrVol,
189
const
polyMesh
&,
190
const
vectorField
&
cellCentres
,
191
const
pointField
& p,
192
const
labelList
& checkFaces,
193
const
List<labelPair>
& baffles,
194
labelHashSet
*
195
);
196
197
//- See primitiveMesh
198
static
bool
checkFaceSkewness
199
(
200
const
bool
report,
201
const
scalar internalSkew,
202
const
scalar boundarySkew,
203
const
polyMesh
& mesh,
204
const
vectorField
&
cellCentres
,
205
const
vectorField
&
faceCentres
,
206
const
vectorField
&
faceAreas
,
207
const
labelList
& checkFaces,
208
const
List<labelPair>
& baffles,
209
labelHashSet
* setPtr
210
);
211
212
//- Interpolation weights (0.5 for regular mesh)
213
static
bool
checkFaceWeights
214
(
215
const
bool
report,
216
const
scalar warnWeight,
217
const
polyMesh
& mesh,
218
const
vectorField
&
cellCentres
,
219
const
vectorField
&
faceCentres
,
220
const
vectorField
&
faceAreas
,
221
const
labelList
& checkFaces,
222
const
List<labelPair>
& baffles,
223
labelHashSet
* setPtr
224
);
225
226
//- Cell volume ratio of neighbouring cells (1 for regular mesh)
227
static
bool
checkVolRatio
228
(
229
const
bool
report,
230
const
scalar warnRatio,
231
const
polyMesh
& mesh,
232
const
scalarField
&
cellVolumes
,
233
const
labelList
& checkFaces,
234
const
List<labelPair>
& baffles,
235
labelHashSet
* setPtr
236
);
237
238
//- See primitiveMesh
239
static
bool
checkFaceAngles
240
(
241
const
bool
report,
242
const
scalar maxDeg,
243
const
polyMesh
& mesh,
244
const
vectorField
&
faceAreas
,
245
const
pointField
& p,
246
const
labelList
& checkFaces,
247
labelHashSet
* setPtr
248
);
249
250
//- Triangle (from face-centre decomposition) normal v.s.
251
// average face normal
252
static
bool
checkFaceTwist
253
(
254
const
bool
report,
255
const
scalar minTwist,
256
const
polyMesh
&,
257
const
vectorField
&
cellCentres
,
258
const
vectorField
&
faceAreas
,
259
const
vectorField
&
faceCentres
,
260
const
pointField
& p,
261
const
labelList
& checkFaces,
262
labelHashSet
* setPtr
263
);
264
265
//- Consecutive triangle (from face-centre decomposition) normals
266
static
bool
checkTriangleTwist
267
(
268
const
bool
report,
269
const
scalar minTwist,
270
const
polyMesh
&,
271
const
vectorField
&
faceAreas
,
272
const
vectorField
&
faceCentres
,
273
const
pointField
& p,
274
const
labelList
& checkFaces,
275
labelHashSet
* setPtr
276
);
277
278
//- Small faces
279
static
bool
checkFaceArea
280
(
281
const
bool
report,
282
const
scalar minArea,
283
const
polyMesh
&,
284
const
vectorField
&
faceAreas
,
285
const
labelList
& checkFaces,
286
labelHashSet
* setPtr
287
);
288
289
static
bool
checkCellDeterminant
290
(
291
const
bool
report,
292
const
scalar minDet,
293
const
polyMesh
&,
294
const
vectorField
&
faceAreas
,
295
const
labelList
& checkFaces,
296
const
labelList
&
affectedCells
,
297
labelHashSet
* setPtr
298
);
299
300
301
// Checking of selected faces with local geometry. Uses above static
302
// functions. Parallel aware.
303
304
bool
checkFaceDotProduct
305
(
306
const
bool
report,
307
const
scalar orthWarn,
308
const
labelList
& checkFaces,
309
const
List<labelPair>
& baffles,
310
labelHashSet
* setPtr
311
)
const
;
312
313
bool
checkFacePyramids
314
(
315
const
bool
report,
316
const
scalar minPyrVol,
317
const
pointField
& p,
318
const
labelList
& checkFaces,
319
const
List<labelPair>
& baffles,
320
labelHashSet
* setPtr
321
)
const
;
322
323
bool
checkFaceSkewness
324
(
325
const
bool
report,
326
const
scalar internalSkew,
327
const
scalar boundarySkew,
328
const
labelList
& checkFaces,
329
const
List<labelPair>
& baffles,
330
labelHashSet
* setPtr
331
)
const
;
332
333
bool
checkFaceWeights
334
(
335
const
bool
report,
336
const
scalar warnWeight,
337
const
labelList
& checkFaces,
338
const
List<labelPair>
& baffles,
339
labelHashSet
* setPtr
340
)
const
;
341
342
bool
checkVolRatio
343
(
344
const
bool
report,
345
const
scalar warnRatio,
346
const
labelList
& checkFaces,
347
const
List<labelPair>
& baffles,
348
labelHashSet
* setPtr
349
)
const
;
350
351
bool
checkFaceAngles
352
(
353
const
bool
report,
354
const
scalar maxDeg,
355
const
pointField
& p,
356
const
labelList
& checkFaces,
357
labelHashSet
* setPtr
358
)
const
;
359
360
bool
checkFaceTwist
361
(
362
const
bool
report,
363
const
scalar minTwist,
364
const
pointField
& p,
365
const
labelList
& checkFaces,
366
labelHashSet
* setPtr
367
)
const
;
368
369
bool
checkTriangleTwist
370
(
371
const
bool
report,
372
const
scalar minTwist,
373
const
pointField
& p,
374
const
labelList
& checkFaces,
375
labelHashSet
* setPtr
376
)
const
;
377
378
bool
checkFaceArea
379
(
380
const
bool
report,
381
const
scalar minArea,
382
const
labelList
& checkFaces,
383
labelHashSet
* setPtr
384
)
const
;
385
386
bool
checkCellDeterminant
387
(
388
const
bool
report,
389
const
scalar warnDet,
390
const
labelList
& checkFaces,
391
const
labelList
&
affectedCells
,
392
labelHashSet
* setPtr
393
)
const
;
394
};
395
396
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
397
398
}
// End namespace Foam
399
400
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
401
402
#endif
403
404
// ************************ vim: set sw=4 sts=4 et: ************************ //