pEqn.H
Go to the documentation of this file.
1rho = thermo.rho();
2
3volScalarField rAU(1.0/UEqn.A());
4surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho*rAU));
5volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p));
6
7if (pimple.nCorrPISO() <= 1)
8{
9 tUEqn.clear();
10}
11
12if (pimple.transonic())
13{
14 surfaceScalarField phid
15 (
16 "phid",
17 fvc::interpolate(psi)
18 *(
19 fvc::flux(HbyA)
20 + MRF.zeroFilter
21 (
22 rhorAUf*fvc::ddtCorr(rho, U, phi)/fvc::interpolate(rho)
23 )
24 )
25 );
26
27 MRF.makeRelative(fvc::interpolate(psi), phid);
28
29 while (pimple.correctNonOrthogonal())
30 {
31 fvScalarMatrix pEqn
32 (
33 fvm::ddt(psi, p)
34 + fvm::div(phid, p)
35 - fvm::laplacian(rhorAUf, p)
36 ==
37 fvOptions(psi, p, rho.name())
38 );
39
40 pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
41
42 if (pimple.finalNonOrthogonalIter())
43 {
44 phi == pEqn.flux();
45 }
46 }
47}
48else
49{
50 surfaceScalarField phiHbyA
51 (
52 "phiHbyA",
53 (
54 fvc::flux(rho*HbyA)
55 + MRF.zeroFilter(rhorAUf*fvc::ddtCorr(rho, U, phi))
56 )
57 );
58
59 MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
60
61 // Update the pressure BCs to ensure flux consistency
63
64 while (pimple.correctNonOrthogonal())
65 {
66 fvScalarMatrix pEqn
67 (
68 fvm::ddt(psi, p)
69 + fvc::div(phiHbyA)
70 - fvm::laplacian(rhorAUf, p)
71 ==
72 fvOptions(psi, p, rho.name())
73 );
74
75 pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
76
77 if (pimple.finalNonOrthogonalIter())
78 {
79 phi = phiHbyA + pEqn.flux();
80 }
81 }
82}
83
84#include "rhoEqn.H"
85#include "compressibleContinuityErrs.H"
86
87// Explicitly relax pressure for momentum corrector
88p.relax();
89
90U = HbyA - rAU*fvc::grad(p);
91U.correctBoundaryConditions();
92fvOptions.correct(U);
93K = 0.5*magSqr(U);
94
95if (pressureControl.limit(p))
96{
97 p.correctBoundaryConditions();
98}
99
100rho = thermo.rho();
101
102if (thermo.dpdt())
103{
104 dpdt = fvc::ddt(p);
105}
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
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
surfaceScalarField phid("phid", fvc::interpolate(psi) *(fvc::flux(HbyA)+MRF.zeroFilter(rhorAUf *fvc::ddtCorr(rho, U, phi)/fvc::interpolate(rho))))
dynamicFvMesh & mesh
volScalarField & dpdt
const pressureControl & pressureControl
tmp< volScalarField > rAU
Definition: initCorrectPhi.H:1