pEqn.H
Go to the documentation of this file.
1{
2 volScalarField rAU(1.0/UEqn.A());
3 volVectorField HbyA("HbyA", U);
4 HbyA = rAU*UEqn.H();
5
6
7 // Define coefficients and pseudo-velocities for RCM interpolation
8 // M[U] = AU - H = -grad(p)
9 // U = H/A - 1/A grad(p)
10 // H/A = U + 1/A grad(p)
11 surfaceScalarField rhorAUf
12 (
13 "rhorAUf",
14 fvc::interpolate(rho)/fvc::interpolate(UEqn.A())
15 );
16
17 surfaceVectorField rhoHbyAf
18 (
19 "rhoHbyAf",
20 fvc::interpolate(rho)*fvc::interpolate(U)
21 + rhorAUf*fvc::interpolate(fvc::grad(p))
22 );
23
24 #include "resetBoundaries.H"
25
26 if (pimple.nCorrPISO() <= 1)
27 {
28 tUEqn.clear();
29 }
30
31 if (pimple.transonic())
32 {
33 FatalError
34 << "\nTransonic option not available for " << args.executable()
35 << exit(FatalError);
36 }
37 else
38 {
39 // Rhie & Chow interpolation (part 1)
40 surfaceScalarField phiHbyA
41 (
42 "phiHbyA",
43 (
44 (rhoHbyAf & mesh.Sf())
45 + rhorAUf*fvc::interpolate(rho)*fvc::ddtCorr(U, phiByRho)
46 + fvc::interpolate(rho)
47 * fvc::alphaCorr(U, phiByRho, pimple.finalInnerIter())
48 )
49 );
50
51 MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
52
53 // Non-orthogonal pressure corrector loop
54 while (pimple.correctNonOrthogonal())
55 {
56 // Pressure corrector
57 fvScalarMatrix pEqn
58 (
59 fvm::ddt(psi, p)
60 + fvc::div(phiHbyA)
61 - fvm::laplacian(rhorAUf, p)
62 ==
63 fvOptions(psi, p, rho.name())
64 );
65
66 pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
67
68 // Rhie & Chow interpolation (part 2)
69 if (pimple.finalNonOrthogonalIter())
70 {
71 phi = phiHbyA + pEqn.flux();
72 }
73 }
74 }
75
76 phiByRho = phi/fvc::interpolate(rho);
77
78 #include "rhoEqn.H"
79 #include "compressibleContinuityErrs.H"
80
81 // Explicitly relax pressure for momentum corrector
82 p.relax();
83
84 U = HbyA - rAU*fvc::grad(p);
85 U.correctBoundaryConditions();
86 fvOptions.correct(U);
87}
88
89rho = thermo.rho();
fv::options & fvOptions
surfaceScalarField & phi
IOMRFZoneList & MRF
pimpleControl & pimple
const word & executable() const noexcept
Name of executable without the path.
Definition: argListI.H:51
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
phiByRho
Definition: pEqn.H:76
surfaceVectorField rhoHbyAf("rhoHbyAf", fvc::interpolate(rho) *fvc::interpolate(U)+rhorAUf *fvc::interpolate(fvc::grad(p)))
dynamicFvMesh & mesh
tmp< volScalarField > rAU
Definition: initCorrectPhi.H:1
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
Foam::argList args(argc, argv)