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
applications
solvers
heatTransfer
buoyantSimpleFoam
createFields.H
Go to the documentation of this file.
1
Info
<<
"Reading thermophysical properties\n"
<<
endl
;
2
3
autoPtr<basicPsiThermo>
pThermo
4
(
5
basicPsiThermo::New(
mesh
)
6
);
7
basicPsiThermo&
thermo
=
pThermo
();
8
9
volScalarField
rho
10
(
11
IOobject
12
(
13
"rho"
,
14
runTime.timeName(),
15
mesh
,
16
IOobject::NO_READ,
17
IOobject::NO_WRITE
18
),
19
thermo.rho()
20
);
21
22
volScalarField
&
p
= thermo.p();
23
volScalarField
&
h
= thermo.h();
24
const
volScalarField
&
psi
= thermo.psi();
25
26
Info
<<
"Reading field U\n"
<<
endl
;
27
volVectorField
U
28
(
29
IOobject
30
(
31
"U"
,
32
runTime.timeName(),
33
mesh
,
34
IOobject::MUST_READ,
35
IOobject::AUTO_WRITE
36
),
37
mesh
38
);
39
40
#include <
finiteVolume/compressibleCreatePhi.H
>
41
42
Info
<<
"Creating turbulence model\n"
<<
endl
;
43
autoPtr<compressible::RASModel>
turbulence
44
(
45
compressible::RASModel::New
46
(
47
rho
,
48
U
,
49
phi
,
50
thermo
51
)
52
);
53
54
55
Info
<<
"Calculating field g.h\n"
<<
endl
;
56
volScalarField
gh
(
"gh"
, g &
mesh
.C());
57
surfaceScalarField
ghf
(
"ghf"
, g &
mesh
.Cf());
58
59
Info
<<
"Reading field p_rgh\n"
<<
endl
;
60
volScalarField
p_rgh
61
(
62
IOobject
63
(
64
"p_rgh"
,
65
runTime.timeName(),
66
mesh
,
67
IOobject::MUST_READ,
68
IOobject::AUTO_WRITE
69
),
70
mesh
71
);
72
73
// Force p_rgh to be consistent with p
74
p_rgh
= p -
rho
*
gh
;
75
76
77
label
pRefCell
= 0;
78
scalar
pRefValue
= 0.0;
79
setRefCell
80
(
81
p,
82
p_rgh
,
83
mesh
.solutionDict().subDict(
"SIMPLE"
),
84
pRefCell
,
85
pRefValue
86
);
87
88
dimensionedScalar
initialMass
=
fvc::domainIntegrate
(
rho
);
89
dimensionedScalar
totalVolume =
sum
(
mesh
.V());