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
compressible
rhoPimpleFoam
pEqn.H
Go to the documentation of this file.
1
rho
=
thermo
.rho();
2
3
volScalarField
rUA
= 1.0/
UEqn
().A();
4
U
= rUA*
UEqn
().H();
5
6
if
(
nCorr
<= 1)
7
{
8
UEqn
.clear();
9
}
10
11
if
(
transonic
)
12
{
13
surfaceScalarField
phid
14
(
15
"phid"
,
16
fvc::interpolate
(
psi
)
17
*(
18
(
fvc::interpolate
(
U
) &
mesh
.Sf())
19
+
fvc::ddtPhiCorr
(rUA,
rho
,
U
,
phi
)
20
)
21
);
22
23
for
(
int
nonOrth=0; nonOrth<=
nNonOrthCorr
; nonOrth++)
24
{
25
fvScalarMatrix
pEqn
26
(
27
fvm::ddt
(
psi
,
p
)
28
+
fvm::div
(phid,
p
)
29
-
fvm::laplacian
(
rho
*rUA,
p
)
30
);
31
32
if
33
(
34
oCorr ==
nOuterCorr
-1
35
&& corr ==
nCorr
-1
36
&& nonOrth ==
nNonOrthCorr
37
)
38
{
39
pEqn.solve(
mesh
.solver(
"pFinal"
));
40
}
41
else
42
{
43
pEqn.solve();
44
}
45
46
if
(nonOrth ==
nNonOrthCorr
)
47
{
48
phi
== pEqn.flux();
49
}
50
}
51
}
52
else
53
{
54
phi
=
55
fvc::interpolate
(
rho
)*
56
(
57
(
fvc::interpolate
(
U
) &
mesh
.Sf())
58
+
fvc::ddtPhiCorr
(rUA,
rho
,
U
,
phi
)
59
);
60
61
for
(
int
nonOrth=0; nonOrth<=
nNonOrthCorr
; nonOrth++)
62
{
63
// Pressure corrector
64
fvScalarMatrix
pEqn
65
(
66
fvm::ddt
(
psi
,
p
)
67
+
fvc::div
(
phi
)
68
-
fvm::laplacian
(
rho
*rUA,
p
)
69
);
70
71
if
72
(
73
oCorr ==
nOuterCorr
-1
74
&& corr ==
nCorr
-1
75
&& nonOrth ==
nNonOrthCorr
76
)
77
{
78
pEqn.solve(
mesh
.solver(
"pFinal"
));
79
}
80
else
81
{
82
pEqn.solve();
83
}
84
85
if
(nonOrth ==
nNonOrthCorr
)
86
{
87
phi
+= pEqn.flux();
88
}
89
}
90
}
91
92
#include <
finiteVolume/rhoEqn.H
>
93
#include <
finiteVolume/compressibleContinuityErrs.H
>
94
95
//if (oCorr != nOuterCorr-1)
96
{
97
// Explicitly relax pressure for momentum corrector
98
p
.relax();
99
100
rho
=
thermo
.rho();
101
rho
.relax();
102
Info
<<
"rho max/min : "
<<
max
(
rho
).value()
103
<<
" "
<<
min
(
rho
).value() <<
endl
;
104
}
105
106
U
-= rUA*
fvc::grad
(
p
);
107
U
.correctBoundaryConditions();
108
109
DpDt
=
fvc::DDt
(
surfaceScalarField
(
"phiU"
,
phi
/
fvc::interpolate
(
rho
)),
p
);
110
111
bound
(
p
,
pMin
);
112
113
// For closed-volume cases adjust the pressure and density levels
114
// to obey overall mass continuity
115
/*
116
if (closedVolume)
117
{
118
p += (initialMass - fvc::domainIntegrate(psi*p))
119
/fvc::domainIntegrate(psi);
120
}
121
*/
122
123
// ************************ vim: set sw=4 sts=4 et: ************************ //