pEqn.H
Go to the documentation of this file.
1{
2 surfaceScalarField faceMask(localMin<scalar>(mesh).interpolate(cellMask));
3
4 volScalarField rAU(1.0/UEqn.A());
5 surfaceScalarField rhorAUf("rhorAUf", faceMask*fvc::interpolate(rho*rAU));
6 volVectorField HbyA("HbyA", U);
7 HbyA = constrainHbyA(cellMask*rAU*UEqn.H(), U, p);
8 tUEqn.clear();
9
10 bool closedVolume = false;
11
12 surfaceScalarField phiHbyA("phiHbyA", fvc::interpolate(rho)*fvc::flux(HbyA));
13 MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
14
15 // Update the pressure BCs to ensure flux consistency
17
18 if (simple.transonic())
19 {
20 surfaceScalarField phid
21 (
22 "phid",
23 (fvc::interpolate(psi)/fvc::interpolate(rho))*phiHbyA
24 );
25
26 phiHbyA -= fvc::interpolate(psi*p)*phiHbyA/fvc::interpolate(rho);
27
28 while (simple.correctNonOrthogonal())
29 {
30 fvScalarMatrix pEqn
31 (
32 fvc::div(phiHbyA)
33 + fvm::div(phid, p)
34 - fvm::laplacian(rhorAUf, p)
35 ==
36 fvOptions(psi, p, rho.name())
37 );
38
39 // Relax the pressure equation to ensure diagonal-dominance
40 pEqn.relax();
41
42 pEqn.setReference
43 (
44 pressureControl.refCell(),
45 pressureControl.refValue()
46 );
47
48 pEqn.solve();
49
50 if (simple.finalNonOrthogonalIter())
51 {
52 phi = phiHbyA + pEqn.flux();
53 }
54 }
55 }
56 else
57 {
59
61 {
62 oversetAdjustPhi(phiHbyA, U);
63 }
64
65 while (simple.correctNonOrthogonal())
66 {
67 fvScalarMatrix pEqn
68 (
69 fvc::div(phiHbyA)
70 - fvm::laplacian(rhorAUf, p)
71 ==
72 fvOptions(psi, p, rho.name())
73 );
74
75 pEqn.setReference
76 (
77 pressureControl.refCell(),
78 pressureControl.refValue()
79 );
80
81 pEqn.solve();
82
83 if (simple.finalNonOrthogonalIter())
84 {
85 phi = phiHbyA + pEqn.flux();
86 }
87 }
88 }
89
91
92 // Explicitly relax pressure for momentum corrector
93 p.relax();
94
95 volVectorField gradP(fvc::grad(p));
96
97 U = HbyA - rAU*cellMask*gradP;
98 U.correctBoundaryConditions();
99 fvOptions.correct(U);
100
101
103
104 // For closed-volume cases adjust the pressure and density levels
105 // to obey overall mass continuity
107 {
108 p += (initialMass - fvc::domainIntegrate(psi*p))
109 /fvc::domainIntegrate(psi);
110 }
111
113 {
114 p.correctBoundaryConditions();
115 }
116
117 rho = thermo.rho();
118
119 if (!simple.transonic())
120 {
121 rho.relax();
122 }
123}
fv::options & fvOptions
surfaceScalarField & phi
IOMRFZoneList & MRF
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 faceMask(localMin< scalar >(mesh).interpolate(cellMask))
volVectorField gradP(fvc::grad(p))
mesh interpolate(rAU)
bool closedVolume
Definition: pEqn.H:10
bool pLimited
Definition: pEqn.H:102
surfaceScalarField phid("phid", fvc::interpolate(psi) *(fvc::flux(HbyA)+MRF.zeroFilter(rhorAUf *fvc::ddtCorr(rho, U, phi)/fvc::interpolate(rho))))
dynamicFvMesh & mesh
const pressureControl & pressureControl
adjustPhi(phiHbyA, U, p_rgh)
tmp< volScalarField > rAU
Definition: initCorrectPhi.H:1
adjustFringe
Definition: readControls.H:16
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
const dictionary & simple
dimensionedScalar initialMass
Definition: createFields.H:57
Calculates and prints the continuity errors.