pEqn.H
Go to the documentation of this file.
1{
2 volScalarField rAU("rAU", 1.0/UEqn.A());
3 surfaceScalarField rAUf("rAUf", fvc::interpolate(rAU));
4 volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p_rgh));
5 surfaceScalarField phiHbyA
6 (
7 "phiHbyA",
8 fvc::flux(HbyA)
9 + fvc::interpolate(rho*rAU)*fvc::ddtCorr(U, phi)
10 );
12
13 surfaceScalarField phig
14 (
15 (
16 interface.surfaceTensionForce()
17 - ghf*fvc::snGrad(rho)
18 )*rAUf*mesh.magSf()
19 );
20
22
23 // Update the pressure BCs to ensure flux consistency
25
26 Pair<tmp<volScalarField>> vDotP = mixture->vDotP();
27 const volScalarField& vDotcP = vDotP[0]();
28 const volScalarField& vDotvP = vDotP[1]();
29
30 while (pimple.correctNonOrthogonal())
31 {
32 fvScalarMatrix p_rghEqn
33 (
34 fvc::div(phiHbyA) - fvm::laplacian(rAUf, p_rgh)
35 - (vDotvP - vDotcP)*(mixture->pSat() - rho*gh)
36 + fvm::Sp(vDotvP - vDotcP, p_rgh)
37 );
38
39 p_rghEqn.setReference(pRefCell, pRefValue);
40
41 p_rghEqn.solve(mesh.solver(p_rgh.select(pimple.finalInnerIter())));
42
43 if (pimple.finalNonOrthogonalIter())
44 {
45 phi = phiHbyA + p_rghEqn.flux();
46
47 U = HbyA + rAU*fvc::reconstruct((phig + p_rghEqn.flux())/rAUf);
48 U.correctBoundaryConditions();
49 fvOptions.correct(U);
50 }
51 }
52
53 p == p_rgh + rho*gh;
54
55 if (p_rgh.needReference())
56 {
57 p += dimensionedScalar
58 (
59 "p",
60 p.dimensions(),
61 pRefValue - getRefCellValue(p, pRefCell)
62 );
63 p_rgh = p - rho*gh;
64 }
65}
volScalarField & p_rgh
fv::options & fvOptions
surfaceScalarField & phi
const scalar pRefValue
const surfaceScalarField & ghf
const label pRefCell
const volScalarField & gh
pimpleControl & pimple
U
Definition: pEqn.H:72
volScalarField & p
surfaceScalarField phig("phig", -rhorAUf *ghf *fvc::snGrad(rho) *mesh.magSf())
fvVectorMatrix & UEqn
Definition: UEqn.H:13
phiHbyA
Definition: pcEqn.H:73
HbyA
Definition: pcEqn.H:74
dynamicFvMesh & mesh
surfaceScalarField rAUf("rAUf", fvc::interpolate(rAU))
adjustPhi(phiHbyA, U, p_rgh)
tmp< volScalarField > rAU
Definition: initCorrectPhi.H:1
const volScalarField & vDotcP
Definition: pEqn.H:33
Pair< tmp< volScalarField > > vDotP
Definition: pEqn.H:32
const volScalarField & vDotvP
Definition: pEqn.H:34
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
interfaceProperties interface(alpha1, U, thermo->transportPropertiesDict())
Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n"<< endl;autoPtr< temperaturePhaseChangeTwoPhaseMixture > mixture
Definition: createFields.H:39