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
rhoSimpleFoam
pEqn.H
Go to the documentation of this file.
1
rho
=
thermo
.rho();
2
rho
=
max
(
rho
,
rhoMin
);
3
rho
=
min
(
rho
,
rhoMax
);
4
rho
.relax();
5
6
volScalarField
rUA
= 1.0/
UEqn
().A();
7
U
= rUA*
UEqn
().H();
8
UEqn
.clear();
9
10
bool
closedVolume
=
false
;
11
12
if
(
transonic
)
13
{
14
surfaceScalarField
phid
15
(
16
"phid"
,
17
fvc::interpolate
(
psi
)*(
fvc::interpolate
(
U
) &
mesh
.Sf())
18
);
19
20
for
(
int
nonOrth=0; nonOrth<=
nNonOrthCorr
; nonOrth++)
21
{
22
fvScalarMatrix
pEqn
23
(
24
fvm::div
(phid,
p
)
25
-
fvm::laplacian
(
rho
*rUA,
p
)
26
);
27
28
// Relax the pressure equation to ensure diagonal-dominance
29
pEqn.relax(
mesh
.relaxationFactor(
"pEqn"
));
30
31
pEqn.setReference(
pRefCell
,
pRefValue
);
32
33
// retain the residual from the first iteration
34
if
(nonOrth == 0)
35
{
36
eqnResidual
= pEqn.solve().initialResidual();
37
maxResidual
=
max
(
eqnResidual
,
maxResidual
);
38
}
39
else
40
{
41
pEqn.solve();
42
}
43
44
if
(nonOrth ==
nNonOrthCorr
)
45
{
46
phi
== pEqn.flux();
47
}
48
}
49
}
50
else
51
{
52
phi
=
fvc::interpolate
(
rho
)*(
fvc::interpolate
(
U
) &
mesh
.Sf());
53
closedVolume =
adjustPhi
(
phi
,
U
,
p
);
54
55
for
(
int
nonOrth=0; nonOrth<=
nNonOrthCorr
; nonOrth++)
56
{
57
fvScalarMatrix
pEqn
58
(
59
fvm::laplacian
(
rho
*rUA,
p
) ==
fvc::div
(
phi
)
60
);
61
62
pEqn.setReference(
pRefCell
,
pRefValue
);
63
64
// Retain the residual from the first iteration
65
if
(nonOrth == 0)
66
{
67
eqnResidual
= pEqn.solve().initialResidual();
68
maxResidual
=
max
(
eqnResidual
,
maxResidual
);
69
}
70
else
71
{
72
pEqn.solve();
73
}
74
75
if
(nonOrth ==
nNonOrthCorr
)
76
{
77
phi
-= pEqn.flux();
78
}
79
}
80
}
81
82
83
#include <
finiteVolume/continuityErrs.H
>
84
85
// Explicitly relax pressure for momentum corrector
86
p
.relax();
87
88
U
-= rUA*
fvc::grad
(
p
);
89
U
.correctBoundaryConditions();
90
91
// For closed-volume cases adjust the pressure and density levels
92
// to obey overall mass continuity
93
if
(closedVolume)
94
{
95
p
+= (
initialMass
-
fvc::domainIntegrate
(
psi
*
p
))
96
/
fvc::domainIntegrate
(
psi
);
97
}
98
99
rho
=
thermo
.rho();
100
rho
=
max
(
rho
,
rhoMin
);
101
rho
=
min
(
rho
,
rhoMax
);
102
rho
.relax();
103
Info
<<
"rho max/min : "
<<
max
(
rho
).value() <<
" "
<<
min
(
rho
).value() <<
endl
;