pcEqn.H
Go to the documentation of this file.
1rho = thermo.rho();
2
3volScalarField rAU(1.0/UEqn.A());
4volScalarField rAtU(1.0/(1.0/rAU - UEqn.H1()));
5volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p));
6tUEqn.clear();
7
8bool closedVolume = false;
9
10surfaceScalarField phiHbyA("phiHbyA", fvc::interpolate(rho)*fvc::flux(HbyA));
11MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
12
13volScalarField rhorAtU("rhorAtU", rho*rAtU);
14
15// Update the pressure BCs to ensure flux consistency
17
18if (simple.transonic())
19{
20 surfaceScalarField phid
21 (
22 "phid",
23 (fvc::interpolate(psi)/fvc::interpolate(rho))*phiHbyA
24 );
25
26 phiHbyA +=
27 fvc::interpolate(rho*(rAtU - rAU))*fvc::snGrad(p)*mesh.magSf()
28 - fvc::interpolate(psi*p)*phiHbyA/fvc::interpolate(rho);
29
30 HbyA -= (rAU - rAtU)*fvc::grad(p);
31
32 while (simple.correctNonOrthogonal())
33 {
34 fvScalarMatrix pEqn
35 (
36 fvc::div(phiHbyA)
37 + fvm::div(phid, p)
38 - fvm::laplacian(rhorAtU, p)
39 ==
40 fvOptions(psi, p, rho.name())
41 );
42
43 // Relax the pressure equation to maintain diagonal dominance
44 pEqn.relax();
45
46 pEqn.setReference
47 (
48 pressureControl.refCell(),
49 pressureControl.refValue()
50 );
51
52 pEqn.solve();
53
54 if (simple.finalNonOrthogonalIter())
55 {
56 phi = phiHbyA + pEqn.flux();
57 }
58 }
59}
60else
61{
63
64 phiHbyA += fvc::interpolate(rho*(rAtU - rAU))*fvc::snGrad(p)*mesh.magSf();
65 HbyA -= (rAU - rAtU)*fvc::grad(p);
66
67 while (simple.correctNonOrthogonal())
68 {
69 fvScalarMatrix pEqn
70 (
71 fvc::div(phiHbyA)
72 - fvm::laplacian(rhorAtU, p)
73 ==
74 fvOptions(psi, p, rho.name())
75 );
76
77 pEqn.setReference
78 (
79 pressureControl.refCell(),
80 pressureControl.refValue()
81 );
82
83 pEqn.solve();
84
85 if (simple.finalNonOrthogonalIter())
86 {
87 phi = phiHbyA + pEqn.flux();
88 }
89 }
90}
91
92// The incompressible form of the continuity error check is appropriate for
93// steady-state compressible also.
95
96// Explicitly relax pressure for momentum corrector
97p.relax();
98
99U = HbyA - rAtU*fvc::grad(p);
100U.correctBoundaryConditions();
101fvOptions.correct(U);
102
104
105// For closed-volume cases adjust the pressure and density levels
106// to obey overall mass continuity
108{
109 p += (initialMass - fvc::domainIntegrate(psi*p))
110 /fvc::domainIntegrate(psi);
111}
112
114{
115 p.correctBoundaryConditions();
116}
117
118rho = thermo.rho();
119
120if (!simple.transonic())
121{
122 rho.relax();
123}
fv::options & fvOptions
surfaceScalarField & phi
IOMRFZoneList & MRF
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
volScalarField & p
const volScalarField & psi
tmp< fvVectorMatrix > tUEqn(fvm::ddt(rho, U)+fvm::div(phi, U)+MRF.DDt(rho, U)+turbulence->divDevRhoReff(U)==fvOptions(rho, U))
fvVectorMatrix & UEqn
Definition: UEqn.H:13
volScalarField rhorAtU("rhorAtU", rho *rAtU)
phiHbyA
Definition: pcEqn.H:73
U
Definition: pcEqn.H:107
rho
Definition: pcEqn.H:1
volScalarField rAtU(1.0/(1.0/rAU - UEqn.H1()))
HbyA
Definition: pcEqn.H:74
bool closedVolume
Definition: pcEqn.H:8
bool pLimited
Definition: pcEqn.H:103
surfaceScalarField phid("phid", fvc::interpolate(psi) *(fvc::flux(HbyA)+MRF.zeroFilter(rhorAUf *fvc::ddtCorr(rho, U, phi)/fvc::interpolate(rho))))
dynamicFvMesh & mesh
const pressureControl & pressureControl
adjustPhi(phiHbyA, U, p_rgh)
tmp< volScalarField > rAU
Definition: initCorrectPhi.H:1
const dictionary & simple
dimensionedScalar initialMass
Definition: createFields.H:57
Calculates and prints the continuity errors.