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
chtMultiRegionFoam
fluid
pEqn.H
Go to the documentation of this file.
1
{
2
bool
closedVolume
=
p_rgh
.needReference();
3
dimensionedScalar
compressibility =
fvc::domainIntegrate
(
psi
);
4
bool
compressible = (compressibility.value() > SMALL);
5
6
rho
=
thermo
.rho();
7
8
volScalarField
rUA
= 1.0/
UEqn
().A();
9
surfaceScalarField
rhorUAf
(
"(rho*(1|A(U)))"
,
fvc::interpolate
(
rho
*rUA));
10
11
U
= rUA*
UEqn
().H();
12
13
surfaceScalarField
phiU
14
(
15
fvc::interpolate
(
rho
)
16
*(
17
(
fvc::interpolate
(
U
) &
mesh
.Sf())
18
+
fvc::ddtPhiCorr
(rUA,
rho
,
U
,
phi
)
19
)
20
);
21
22
phi
=
phiU
-
rhorUAf
*
ghf
*
fvc::snGrad
(
rho
)*
mesh
.magSf();
23
24
{
25
fvScalarMatrix
p_rghDDtEqn
26
(
27
fvc::ddt
(
rho
) +
psi
*
correction
(
fvm::ddt
(
p_rgh
))
28
+
fvc::div
(
phi
)
29
);
30
31
// Thermodynamic density needs to be updated by psi*d(p) after the
32
// pressure solution - done in 2 parts. Part 1:
33
thermo
.rho() -=
psi
*
p_rgh
;
34
35
for
(
int
nonOrth=0; nonOrth<=
nNonOrthCorr
; nonOrth++)
36
{
37
fvScalarMatrix
p_rghEqn
38
(
39
p_rghDDtEqn
40
-
fvm::laplacian
(
rhorUAf
, p_rgh)
41
);
42
43
p_rghEqn.solve
44
(
45
mesh
.solver
46
(
47
p_rgh.select
48
(
49
(
50
oCorr ==
nOuterCorr
-1
51
&& corr ==
nCorr
-1
52
&& nonOrth ==
nNonOrthCorr
53
)
54
)
55
)
56
);
57
58
if
(nonOrth ==
nNonOrthCorr
)
59
{
60
phi
+= p_rghEqn.flux();
61
}
62
}
63
64
// Second part of thermodynamic density update
65
thermo
.rho() +=
psi
*
p_rgh
;
66
}
67
68
// Correct velocity field
69
U
+= rUA*
fvc::reconstruct
((
phi
-
phiU
)/
rhorUAf
);
70
U
.correctBoundaryConditions();
71
72
p
= p_rgh +
rho
*
gh
;
73
74
// Update pressure substantive derivative
75
DpDt
=
fvc::DDt
(
surfaceScalarField
(
"phiU"
,
phi
/
fvc::interpolate
(
rho
)),
p
);
76
77
// Solve continuity
78
#include <
finiteVolume/rhoEqn.H
>
79
80
// Update continuity errors
81
#include "
compressibleContinuityErrors.H
"
82
83
// For closed-volume cases adjust the pressure and density levels
84
// to obey overall mass continuity
85
if
(closedVolume && compressible)
86
{
87
p
+= (
initialMass
-
fvc::domainIntegrate
(
thermo
.rho()))
88
/compressibility;
89
rho
=
thermo
.rho();
90
p_rgh =
p
-
rho
*
gh
;
91
}
92
93
// Update thermal conductivity
94
K
=
thermoFluid
[i].Cp()*
turb
.alphaEff();
95
}
96
97
// ************************ vim: set sw=4 sts=4 et: ************************ //