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
finiteVolume
fvMesh
fvMeshGeometry.C
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-2011 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
\*---------------------------------------------------------------------------*/
25
26
#include "
fvMesh.H
"
27
#include <
OpenFOAM/Time.H
>
28
#include <
finiteVolume/volFields.H
>
29
#include <
finiteVolume/surfaceFields.H
>
30
#include <
finiteVolume/slicedVolFields.H
>
31
#include <
finiteVolume/slicedSurfaceFields.H
>
32
#include <
OpenFOAM/SubField.H
>
33
#include <
finiteVolume/cyclicFvPatchFields.H
>
34
35
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36
37
namespace
Foam
38
{
39
40
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
41
42
void
fvMesh::makeSf()
const
43
{
44
if
(
debug
)
45
{
46
Info
<<
"void fvMesh::makeSf() : "
47
<<
"assembling face areas"
48
<<
endl
;
49
}
50
51
// It is an error to attempt to recalculate
52
// if the pointer is already set
53
if
(SfPtr_)
54
{
55
FatalErrorIn
(
"fvMesh::makeSf()"
)
56
<<
"face areas already exist"
57
<<
abort
(
FatalError
);
58
}
59
60
SfPtr_ =
new
slicedSurfaceVectorField
61
(
62
IOobject
63
(
64
"S"
,
65
pointsInstance
(),
66
meshSubDir
,
67
*
this
,
68
IOobject::NO_READ
,
69
IOobject::NO_WRITE
,
70
false
71
),
72
*
this
,
73
dimArea
,
74
faceAreas
()
75
);
76
}
77
78
79
void
fvMesh::makeMagSf()
const
80
{
81
if
(
debug
)
82
{
83
Info
<<
"void fvMesh::makeMagSf() : "
84
<<
"assembling mag face areas"
85
<<
endl
;
86
}
87
88
// It is an error to attempt to recalculate
89
// if the pointer is already set
90
if
(magSfPtr_)
91
{
92
FatalErrorIn
(
"void fvMesh::makeMagSf()"
)
93
<<
"mag face areas already exist"
94
<<
abort
(
FatalError
);
95
}
96
97
// Note: Added stabilisation for faces with exactly zero area.
98
// These should be caught on mesh checking but at least this stops
99
// the code from producing Nans.
100
magSfPtr_ =
new
surfaceScalarField
101
(
102
IOobject
103
(
104
"magSf"
,
105
pointsInstance
(),
106
meshSubDir
,
107
*
this
,
108
IOobject::NO_READ
,
109
IOobject::NO_WRITE
,
110
false
111
),
112
mag
(
Sf
()) +
dimensionedScalar
(
"vs"
,
dimArea
, VSMALL)
113
);
114
}
115
116
117
void
fvMesh::makeC()
const
118
{
119
if
(
debug
)
120
{
121
Info
<<
"void fvMesh::makeC() : "
122
<<
"assembling cell centres"
123
<<
endl
;
124
}
125
126
// It is an error to attempt to recalculate
127
// if the pointer is already set
128
if
(CPtr_)
129
{
130
FatalErrorIn
(
"fvMesh::makeC()"
)
131
<<
"cell centres already exist"
132
<<
abort
(
FatalError
);
133
}
134
135
CPtr_ =
new
slicedVolVectorField
136
(
137
IOobject
138
(
139
"C"
,
140
pointsInstance
(),
141
meshSubDir
,
142
*
this
,
143
IOobject::NO_READ
,
144
IOobject::NO_WRITE
,
145
false
146
),
147
*
this
,
148
dimLength
,
149
cellCentres
(),
150
faceCentres
()
151
);
152
153
154
// Need to correct for cyclics transformation since absolute quantity.
155
// Ok on processor patches since hold opposite cell centre (no
156
// transformation)
157
slicedVolVectorField
&
C
= *CPtr_;
158
159
forAll
(C.boundaryField(),
patchi
)
160
{
161
if
(isA<cyclicFvPatchVectorField>(C.boundaryField()[
patchi
]))
162
{
163
// Note: cyclic is not slice but proper field
164
C.boundaryField()[
patchi
] ==
static_cast<
const
vectorField
&
>
165
(
166
static_cast<
const
List<vector>&
>
167
(
168
boundary_[
patchi
].patchSlice(
faceCentres
())
169
)
170
);
171
}
172
}
173
}
174
175
176
void
fvMesh::makeCf()
const
177
{
178
if
(
debug
)
179
{
180
Info
<<
"void fvMesh::makeCf() : "
181
<<
"assembling face centres"
182
<<
endl
;
183
}
184
185
// It is an error to attempt to recalculate
186
// if the pointer is already set
187
if
(CfPtr_)
188
{
189
FatalErrorIn
(
"fvMesh::makeCf()"
)
190
<<
"face centres already exist"
191
<<
abort
(
FatalError
);
192
}
193
194
CfPtr_ =
new
slicedSurfaceVectorField
195
(
196
IOobject
197
(
198
"Cf"
,
199
pointsInstance
(),
200
meshSubDir
,
201
*
this
,
202
IOobject::NO_READ
,
203
IOobject::NO_WRITE
,
204
false
205
),
206
*
this
,
207
dimLength
,
208
faceCentres
()
209
);
210
}
211
212
213
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
214
215
const
volScalarField::DimensionedInternalField
&
fvMesh::V
()
const
216
{
217
if
(!VPtr_)
218
{
219
VPtr_ =
new
slicedVolScalarField::DimensionedInternalField
220
(
221
IOobject
222
(
223
"V"
,
224
time
().
timeName
(),
225
*
this
,
226
IOobject::NO_READ
,
227
IOobject::NO_WRITE
228
),
229
*
this
,
230
dimVolume
,
231
cellVolumes
()
232
);
233
}
234
235
return
*
static_cast<
slicedVolScalarField::DimensionedInternalField
*
>
(VPtr_);
236
}
237
238
239
const
volScalarField::DimensionedInternalField
&
fvMesh::V0
()
const
240
{
241
if
(!V0Ptr_)
242
{
243
FatalErrorIn
(
"fvMesh::V0() const"
)
244
<<
"V0 is not available"
245
<<
abort
(
FatalError
);
246
}
247
248
return
*V0Ptr_;
249
}
250
251
252
volScalarField::DimensionedInternalField
&
fvMesh::setV0
()
253
{
254
if
(!V0Ptr_)
255
{
256
FatalErrorIn
(
"fvMesh::setV0()"
)
257
<<
"V0 is not available"
258
<<
abort
(
FatalError
);
259
}
260
261
return
*V0Ptr_;
262
}
263
264
265
const
volScalarField::DimensionedInternalField
&
fvMesh::V00
()
const
266
{
267
if
(!V00Ptr_)
268
{
269
V00Ptr_ =
new
DimensionedField<scalar, volMesh>
270
(
271
IOobject
272
(
273
"V00"
,
274
time
().
timeName
(),
275
*
this
,
276
IOobject::NO_READ
,
277
IOobject::NO_WRITE
278
),
279
V0
()
280
);
281
282
// If V00 is used then V0 should be stored for restart
283
V0Ptr_->
writeOpt
() =
IOobject::AUTO_WRITE
;
284
}
285
286
return
*V00Ptr_;
287
}
288
289
290
tmp<volScalarField::DimensionedInternalField>
fvMesh::Vsc
()
const
291
{
292
if
(
moving
() &&
time
().subCycling())
293
{
294
const
TimeState
& ts =
time
();
295
const
TimeState
& ts0 =
time
().
prevTimeState
();
296
297
scalar tFrac =
298
(
299
ts.
value
() - (ts0.
value
() - ts0.
deltaTValue
())
300
)/ts0.
deltaTValue
();
301
302
if
(tFrac < (1 - SMALL))
303
{
304
return
V0
() + tFrac*(
V
() -
V0
());
305
}
306
else
307
{
308
return
V
();
309
}
310
}
311
else
312
{
313
return
V
();
314
}
315
}
316
317
318
tmp<volScalarField::DimensionedInternalField>
fvMesh::Vsc0
()
const
319
{
320
if
(
moving
() &&
time
().subCycling())
321
{
322
const
TimeState
& ts =
time
();
323
const
TimeState
& ts0 =
time
().
prevTimeState
();
324
325
scalar t0Frac =
326
(
327
(ts.
value
() - ts.
deltaTValue
())
328
- (ts0.
value
() - ts0.
deltaTValue
())
329
)/ts0.
deltaTValue
();
330
331
if
(t0Frac > SMALL)
332
{
333
return
V0
() + t0Frac*(
V
() -
V0
());
334
}
335
else
336
{
337
return
V0
();
338
}
339
}
340
else
341
{
342
return
V0
();
343
}
344
}
345
346
347
const
surfaceVectorField
&
fvMesh::Sf
()
const
348
{
349
if
(!SfPtr_)
350
{
351
makeSf();
352
}
353
354
return
*SfPtr_;
355
}
356
357
358
const
surfaceScalarField
&
fvMesh::magSf
()
const
359
{
360
if
(!magSfPtr_)
361
{
362
makeMagSf();
363
}
364
365
return
*magSfPtr_;
366
}
367
368
369
const
volVectorField
&
fvMesh::C
()
const
370
{
371
if
(!CPtr_)
372
{
373
makeC();
374
}
375
376
return
*CPtr_;
377
}
378
379
380
const
surfaceVectorField
&
fvMesh::Cf
()
const
381
{
382
if
(!CfPtr_)
383
{
384
makeCf();
385
}
386
387
return
*CfPtr_;
388
}
389
390
391
const
surfaceScalarField
&
fvMesh::phi
()
const
392
{
393
if
(!phiPtr_)
394
{
395
FatalErrorIn
(
"fvMesh::phi()"
)
396
<<
"mesh flux field does not exists, is the mesh actually moving?"
397
<<
exit
(
FatalError
);
398
}
399
400
// Set zero current time
401
// mesh motion fluxes if the time has been incremented
402
if
(phiPtr_->
timeIndex
() !=
time
().
timeIndex
())
403
{
404
(*phiPtr_) =
dimensionedScalar
(
"0"
,
dimVolume
/
dimTime
, 0.0);
405
}
406
407
return
*phiPtr_;
408
}
409
410
411
surfaceScalarField
&
fvMesh::setPhi
()
412
{
413
if
(!phiPtr_)
414
{
415
FatalErrorIn
(
"fvMesh::setPhi()"
)
416
<<
"mesh flux field does not exists, is the mesh actually moving?"
417
<<
exit
(
FatalError
);
418
}
419
420
return
*phiPtr_;
421
}
422
423
424
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
425
426
}
// End namespace Foam
427
428
// ************************ vim: set sw=4 sts=4 et: ************************ //