pEqn.H
Go to the documentation of this file.
1{
2 if (pimple.nCorrPIMPLE() == 1)
3 {
4 p =
5 (
6 rho
7 - alphal*rhol0
8 - ((alphav*psiv + alphal*psil) - psi)*pSat
9 )/psi;
10 }
11
12 surfaceScalarField rhof("rhof", fvc::interpolate(rho));
13
14 volScalarField rAU(1.0/UEqn.A());
15 surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho*rAU));
16 volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p));
17
18 phi = fvc::flux(HbyA)
19 + rhorAUf*fvc::ddtCorr(U, phi);
20
21 surfaceScalarField phiGradp(rhorAUf*mesh.magSf()*fvc::snGrad(p));
22
23 phi -= phiGradp/rhof;
24
25 while (pimple.correctNonOrthogonal())
26 {
27 fvScalarMatrix pEqn
28 (
29 fvm::ddt(psi, p)
30 - (rhol0 + (psil - psiv)*pSat)*fvc::ddt(alphav) - pSat*fvc::ddt(psi)
31 + fvc::div(phi, rho)
32 + fvc::div(phiGradp)
33 - fvm::laplacian(rhorAUf, p)
34 );
35
36 pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
37
38 if (pimple.finalNonOrthogonalIter())
39 {
40 phi += (phiGradp + pEqn.flux())/rhof;
41 }
42 }
43
44 Info<< "Predicted p max-min : " << max(p).value()
45 << " " << min(p).value() << endl;
46
47 rho == max
48 (
49 psi*p
50 + alphal*rhol0
51 + ((alphav*psiv + alphal*psil) - psi)*pSat,
52 rhoMin
53 );
54
55 #include "alphavPsi.H"
56
57 p =
58 (
59 rho
60 - alphal*rhol0
61 - ((alphav*psiv + alphal*psil) - psi)*pSat
62 )/psi;
63
64 p.correctBoundaryConditions();
65
66 Info<< "Phase-change corrected p max-min : " << max(p).value()
67 << " " << min(p).value() << endl;
68
69 // Correct velocity
70
71 U = HbyA - rAU*fvc::grad(p);
72
73 // Remove the swirl component of velocity for "wedge" cases
74 if (pimple.dict().found("removeSwirl"))
75 {
76 label swirlCmpt(pimple.dict().get<label>("removeSwirl"));
77
78 Info<< "Removing swirl component-" << swirlCmpt << " of U" << endl;
79 U.field().replace(swirlCmpt, Zero);
80 }
81
82 U.correctBoundaryConditions();
83
84 Info<< "max(U) " << max(mag(U)).value() << endl;
85}
Y[inertIndex] max(0.0)
alphal
Definition: alphavPsi.H:12
surfaceScalarField & phi
pimpleControl & pimple
U
Definition: pEqn.H:72
volScalarField & p
const volScalarField & psi
surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho *rAU))
fvVectorMatrix & UEqn
Definition: UEqn.H:13
HbyA
Definition: pcEqn.H:74
dynamicFvMesh & mesh
const dimensionedScalar rhoMin
tmp< volScalarField > rAU
Definition: initCorrectPhi.H:1
surfaceScalarField rhof(fvc::interpolate(rho, "div(phi,rho)"))
surfaceScalarField phiGradp(rhorAUf *mesh.magSf() *fvc::snGrad(p))
messageStream Info
Information stream (stdout output on master, null elsewhere)
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33