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
combustion
rhoReactingFoam
pEqn.H
Go to the documentation of this file.
1
{
2
rho
=
thermo
.rho();
3
4
volScalarField
rUA
= 1.0/
UEqn
.A();
5
U
= rUA*
UEqn
.H();
6
7
if
(
transonic
)
8
{
9
surfaceScalarField
phiv =
10
(
fvc::interpolate
(
U
) &
mesh
.Sf())
11
+
fvc::ddtPhiCorr
(rUA,
rho
,
U
,
phi
);
12
13
phi
=
fvc::interpolate
(
rho
)*phiv;
14
15
surfaceScalarField
phid
16
(
17
"phid"
,
18
fvc::interpolate
(
thermo
.psi())*phiv
19
);
20
21
fvScalarMatrix
pDDtEqn
22
(
23
fvc::ddt
(
rho
) +
fvc::div
(
phi
)
24
+
correction
(
fvm::ddt
(
psi
,
p
) +
fvm::div
(
phid
,
p
))
25
);
26
27
// Thermodynamic density needs to be updated by psi*d(p) after the
28
// pressure solution - done in 2 parts. Part 1:
29
thermo
.rho() -=
psi
*
p
;
30
31
for
(
int
nonOrth=0; nonOrth<=
nNonOrthCorr
; nonOrth++)
32
{
33
fvScalarMatrix
pEqn
34
(
35
pDDtEqn -
fvm::laplacian
(
rho
*rUA, p)
36
);
37
38
if
(ocorr ==
nOuterCorr
&& corr ==
nCorr
&& nonOrth ==
nNonOrthCorr
)
39
{
40
pEqn.solve(
mesh
.solver(p.name() +
"Final"
));
41
}
42
else
43
{
44
pEqn.solve();
45
}
46
47
if
(nonOrth ==
nNonOrthCorr
)
48
{
49
phi
+= pEqn.flux();
50
}
51
}
52
53
// Second part of thermodynamic density update
54
thermo
.rho() +=
psi
*
p
;
55
}
56
else
57
{
58
phi
=
59
fvc::interpolate
(
rho
)
60
*(
61
(
fvc::interpolate
(
U
) &
mesh
.Sf())
62
+
fvc::ddtPhiCorr
(rUA,
rho
,
U
,
phi
)
63
);
64
65
fvScalarMatrix
pDDtEqn
66
(
67
fvc::ddt
(
rho
) +
psi
*
correction
(
fvm::ddt
(p))
68
+
fvc::div
(
phi
)
69
);
70
71
// Thermodynamic density needs to be updated by psi*d(p) after the
72
// pressure solution - done in 2 parts. Part 1:
73
thermo
.rho() -=
psi
*
p
;
74
75
for
(
int
nonOrth=0; nonOrth<=
nNonOrthCorr
; nonOrth++)
76
{
77
fvScalarMatrix
pEqn
78
(
79
pDDtEqn -
fvm::laplacian
(
rho
*rUA, p)
80
);
81
82
if
(ocorr ==
nOuterCorr
&& corr ==
nCorr
&& nonOrth ==
nNonOrthCorr
)
83
{
84
pEqn.solve(
mesh
.solver(p.name() +
"Final"
));
85
}
86
else
87
{
88
pEqn.solve();
89
}
90
91
if
(nonOrth ==
nNonOrthCorr
)
92
{
93
phi
+= pEqn.flux();
94
}
95
}
96
97
// Second part of thermodynamic density update
98
thermo
.rho() +=
psi
*
p
;
99
}
100
101
#include <
finiteVolume/rhoEqn.H
>
102
#include <
finiteVolume/compressibleContinuityErrs.H
>
103
104
U
-= rUA*
fvc::grad
(p);
105
U
.correctBoundaryConditions();
106
107
DpDt
=
fvc::DDt
(
surfaceScalarField
(
"phiU"
,
phi
/
fvc::interpolate
(
rho
)), p);
108
}