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));
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 fvc::interpolate(rho*rAU)*fvc::ddtCorr(rho, U, phi)
23 /fvc::interpolate(rho)
24 )
25 )
26 );
27
28 MRF.makeRelative(fvc::interpolate(psi), phid);
29
30 surfaceScalarField phic
31 (
32 "phic",
33 fvc::interpolate(rho*(rAtU - rAU))*fvc::snGrad(p)*mesh.magSf()
34 );
35
36 HbyA -= (rAU - rAtU)*fvc::grad(p);
37
38 volScalarField rhorAtU("rhorAtU", rho*rAtU);
39
40 while (pimple.correctNonOrthogonal())
41 {
42 fvScalarMatrix pEqn
43 (
44 fvm::ddt(psi, p)
45 + fvm::div(phid, p)
46 + fvc::div(phic)
47 - fvm::laplacian(rhorAtU, p)
48 ==
49 fvOptions(psi, p, rho.name())
50 );
51
52 pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
53
54 if (pimple.finalNonOrthogonalIter())
55 {
56 phi == phic + pEqn.flux();
57 }
58 }
59}
60else
61{
62 surfaceScalarField phiHbyA
63 (
64 "phiHbyA",
65 (
66 fvc::flux(rho*HbyA)
67 + MRF.zeroFilter(fvc::interpolate(rho*rAU)*fvc::ddtCorr(rho, U, phi))
68 )
69 );
70
71 MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
72
73 phiHbyA += fvc::interpolate(rho*(rAtU - rAU))*fvc::snGrad(p)*mesh.magSf();
74 HbyA -= (rAU - rAtU)*fvc::grad(p);
75
76 volScalarField rhorAtU("rhorAtU", rho*rAtU);
77
78 // Update the pressure BCs to ensure flux consistency
80
81 while (pimple.correctNonOrthogonal())
82 {
83 fvScalarMatrix pEqn
84 (
85 fvm::ddt(psi, p)
86 + fvc::div(phiHbyA)
87 - fvm::laplacian(rhorAtU, p)
88 ==
89 fvOptions(psi, p, rho.name())
90 );
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"
102#include "compressibleContinuityErrs.H"
103
104// Explicitly relax pressure for momentum corrector
105p.relax();
106
107U = HbyA - rAtU*fvc::grad(p);
108U.correctBoundaryConditions();
109fvOptions.correct(U);
110K = 0.5*magSqr(U);
111
113{
114 p.correctBoundaryConditions();
115 rho = thermo.rho();
116}
117
118if (thermo.dpdt())
119{
120 dpdt = fvc::ddt(p);
121}
CGAL::Exact_predicates_exact_constructions_kernel K
surfaceScalarField phic(mixture.cAlpha() *mag(alphaPhic/mesh.magSf()))
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...
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
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