pEqn.H
Go to the documentation of this file.
1if (!pimple.SIMPLErho())
2{
3 rho = thermo.rho();
4}
5// Thermodynamic density needs to be updated by psi*d(p) after the
6// pressure solution
7const volScalarField psip0(psi*p);
8
9volScalarField rAU("rAU", 1.0/UEqn.A());
10mesh.interpolate(rAU);
11
12surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho*rAU));
13volVectorField HbyA("HbyA", U);
14
16
17if (pimple.nCorrPISO() <= 1)
18{
19 tUEqn.clear();
20}
21
22surfaceScalarField phiHbyA
23(
24 "phiHbyA",
25 fvc::interpolate(rho)*fvc::flux(HbyA)
26);
27
29{
30 surfaceScalarField faceMaskOld
31 (
32 localMin<scalar>(mesh).interpolate(cellMask.oldTime())
33 );
34
35 phiHbyA +=
36 faceMaskOld*MRF.zeroFilter(rhorAUf*fvc::ddtCorr(rho, U, phi, rhoUf));
37}
38
39fvc::makeRelative(phiHbyA, rho, U);
40MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
41
42// Update the pressure BCs to ensure flux consistency
44
45if (pimple.transonic())
46{
47 surfaceScalarField phid
48 (
49 "phid",
50 (fvc::interpolate(psi)/fvc::interpolate(rho))*phiHbyA
51 );
52
53 phiHbyA -= fvc::interpolate(psi*p)*phiHbyA/fvc::interpolate(rho);
54
55 fvScalarMatrix pDDtEqn
56 (
57 fvc::ddt(rho) + psi*correction(fvm::ddt(p))
58 + fvc::div(phiHbyA) + fvm::div(phid, p)
59 ==
60 fvOptions(psi, p, rho.name())
61 );
62
63 while (pimple.correctNonOrthogonal())
64 {
65 fvScalarMatrix pEqn(pDDtEqn - fvm::laplacian(rhorAUf, p));
66
67 // Relax the pressure equation to ensure diagonal-dominance
68 pEqn.relax();
69
70 pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
71
72 if (pimple.finalNonOrthogonalIter())
73 {
74 phi = phiHbyA + pEqn.flux();
75 }
76 }
77}
78else
79{
80 fvScalarMatrix pDDtEqn
81 (
82 fvc::ddt(rho) + psi*correction(fvm::ddt(p))
83 + fvc::div(phiHbyA)
84 ==
85 fvOptions(psi, p, rho.name())
86 );
87
88 while (pimple.correctNonOrthogonal())
89 {
90 fvScalarMatrix pEqn(pDDtEqn - fvm::laplacian(rhorAUf, p));
91
92 pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
93
94 if (pimple.finalNonOrthogonalIter())
95 {
96 phi = phiHbyA + pEqn.flux();
97 }
98 }
99}
100
101#include "rhoEqn.H"
103
104// Explicitly relax pressure for momentum corrector
105p.relax();
106
107volVectorField gradP(fvc::grad(p));
108//mesh.interpolate(gradP);
109U = cellMask*(HbyA - rAU*gradP);
110U.correctBoundaryConditions();
111fvOptions.correct(U);
112K = 0.5*magSqr(U);
113
115{
116 p.correctBoundaryConditions();
117}
118
119thermo.correctRho(psi*p - psip0, rhoMin, rhoMax) ;
120rho = thermo.rho();
121
122{
123 // Correct rhoUf if the mesh is moving
124 fvc::correctRhoUf(rhoUf, rho, U, phi);
125}
126
127if (thermo.dpdt())
128{
129 dpdt = fvc::ddt(p);
130
131 if (mesh.moving())
132 {
133 dpdt -= fvc::div(fvc::meshPhi(rho, U), p);
134 }
135}
136
137surfaceScalarField faceMask
138(
139 localMin<scalar>(mesh).interpolate(cellMask)
140);
141phi *= faceMask;
CGAL::Exact_predicates_exact_constructions_kernel K
fv::options & fvOptions
surfaceScalarField & phi
IOMRFZoneList & MRF
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 rhorAUf("rhorAUf", fvc::interpolate(rho *rAU))
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
phiHbyA
Definition: pcEqn.H:73
HbyA
Definition: pcEqn.H:74
const volScalarField psip0(psi *p)
surfaceScalarField faceMask(localMin< scalar >(mesh).interpolate(cellMask))
volVectorField gradP(fvc::grad(p))
mesh interpolate(rAU)
ddtCorr
Definition: readControls.H:9
fvScalarMatrix pDDtEqn(fvc::ddt(rho)+psi *correction(fvm::ddt(p))+fvc::div(phiHbyA)==fvOptions(psi, p, rho.name()))
surfaceScalarField phid("phid", fvc::interpolate(psi) *(fvc::flux(HbyA)+MRF.zeroFilter(rhorAUf *fvc::ddtCorr(rho, U, phi)/fvc::interpolate(rho))))
dynamicFvMesh & mesh
const dimensionedScalar rhoMin
const dimensionedScalar rhoMax
volScalarField & dpdt
const pressureControl & pressureControl
tmp< volScalarField > rAU
Definition: initCorrectPhi.H:1
Calculates and prints the continuity errors.