pEqn.H
Go to the documentation of this file.
1dimensionedScalar compressibility = fvc::domainIntegrate(psi);
2bool compressible = (compressibility.value() > SMALL);
3
4rho = thermo.rho();
5
6// Thermodynamic density needs to be updated by psi*d(p) after the
7// pressure solution
8const volScalarField psip0(psi*p);
9
10volScalarField rAU(1.0/UEqn.A());
11mesh.interpolate(rAU);
12
13surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho*rAU));
14volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p_rgh));
15
16surfaceScalarField phig(-rhorAUf*ghf*fvc::snGrad(rho)*mesh.magSf());
17
18surfaceScalarField phiHbyA
19(
20 "phiHbyA",
21 fvc::flux(rho*HbyA) + phig
22);
23
25{
26 surfaceScalarField faceMaskOld
27 (
28 localMin<scalar>(mesh).interpolate(cellMask.oldTime())
29 );
30
31 phiHbyA +=
32 faceMaskOld*MRF.zeroFilter(rhorAUf*fvc::ddtCorr(rho, U, phi));
33}
34
35MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
36
37// Update the pressure BCs to ensure flux consistency
39
40fvScalarMatrix p_rghDDtEqn
41(
42 fvc::ddt(rho) + psi*correction(fvm::ddt(p_rgh))
43 + fvc::div(phiHbyA)
44 ==
45 fvOptions(psi, p_rgh, rho.name())
46);
47
48while (pimple.correctNonOrthogonal())
49{
50 fvScalarMatrix p_rghEqn
51 (
53 - fvm::laplacian(rhorAUf, p_rgh)
54 );
55
56 p_rghEqn.solve(mesh.solver(p_rgh.select(pimple.finalInnerIter())));
57
58 if (pimple.finalNonOrthogonalIter())
59 {
60 // Calculate the conservative fluxes
61 phi = phiHbyA + p_rghEqn.flux();
62
63 // Explicitly relax pressure for momentum corrector
64 p_rgh.relax();
65
66 // Correct the momentum source with the pressure gradient flux
67 // calculated from the relaxed pressure
68 U =
69 cellMask*
70 (
71 HbyA + rAU*fvc::reconstruct((phig + p_rghEqn.flux())/rhorAUf)
72 );
73 U.correctBoundaryConditions();
74 fvOptions.correct(U);
75 K = 0.5*magSqr(U);
76 }
77}
78
80
81#include "rhoEqn.H"
82#include "compressibleContinuityErrs.H"
83
84if (p_rgh.needReference())
85{
86 if (!compressible)
87 {
88 p += dimensionedScalar
89 (
90 "p",
91 p.dimensions(),
92 pRefValue - getRefCellValue(p, pRefCell)
93 );
94 }
95 else
96 {
97 p += (initialMass - fvc::domainIntegrate(psi*p))
99 thermo.correctRho(psi*p - psip0);
100 rho = thermo.rho();
101 p_rgh = p - rho*gh;
102 }
103}
104else
105{
106 thermo.correctRho(psi*p - psip0);
107}
108
109rho = thermo.rho();
110
111{
112 fvc::correctRhoUf(rhoUf, rho, U, phi);
113}
114
115if (thermo.dpdt())
116{
117 dpdt = fvc::ddt(p);
118
119 if (mesh.moving())
120 {
121 dpdt -= fvc::div(fvc::meshPhi(rho, U), p);
122 }
123}
124
125surfaceScalarField faceMask
126(
127 localMin<scalar>(mesh).interpolate(cellMask)
128);
129phi *= faceMask;
CGAL::Exact_predicates_exact_constructions_kernel K
volScalarField & p_rgh
fv::options & fvOptions
surfaceScalarField & phi
const scalar pRefValue
const surfaceScalarField & ghf
const label pRefCell
IOMRFZoneList & MRF
const volScalarField & gh
pimpleControl & pimple
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
U
Definition: pEqn.H:72
rhoUf
Definition: pEqn.H:89
volScalarField & p
const volScalarField & psi
surfaceScalarField phig("phig", -rhorAUf *ghf *fvc::snGrad(rho) *mesh.magSf())
surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho *rAU))
fvVectorMatrix & UEqn
Definition: UEqn.H:13
phiHbyA
Definition: pcEqn.H:73
HbyA
Definition: pcEqn.H:74
const volScalarField psip0(psi *p)
fvScalarMatrix p_rghDDtEqn(fvc::ddt(rho)+psi *correction(fvm::ddt(p_rgh))+fvc::div(phiHbyA)==fvOptions(psi, p_rgh, rho.name()))
surfaceScalarField faceMask(localMin< scalar >(mesh).interpolate(cellMask))
mesh interpolate(rAU)
ddtCorr
Definition: readControls.H:9
dynamicFvMesh & mesh
volScalarField & dpdt
dimensionedScalar compressibility
Definition: pEqn.H:1
bool compressible
Definition: pEqn.H:2
tmp< volScalarField > rAU
Definition: initCorrectPhi.H:1
dimensionedScalar initialMass
Definition: createFields.H:57