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 + MRF.zeroFilter(fvc::interpolate(rho*rAU)*fvc::ddtCorr(U, phi))
10 );
11 MRF.makeRelative(phiHbyA);
12
13 surfaceScalarField phig
14 (
15 (
16 mixture.surfaceTensionForce()
17 - ghf*fvc::snGrad(rho)
18 )*rAUf*mesh.magSf()
19 );
20
22
23 // Update the pressure BCs to ensure flux consistency
25
26 tmp<fvScalarMatrix> p_rghEqnComp1;
27 tmp<fvScalarMatrix> p_rghEqnComp2;
28
29 if (pimple.transonic())
30 {
31 #include "rhofs.H"
32
33 surfaceScalarField phid1("phid1", fvc::interpolate(psi1)*phi);
34 surfaceScalarField phid2("phid2", fvc::interpolate(psi2)*phi);
35
37 pos(alpha1)
38 *(
39 (
40 fvc::ddt(alpha1, rho1) + fvc::div(alphaPhi1*rho1f)
41 - (fvOptions(alpha1, mixture.thermo1().rho())&rho1)
42 )/rho1
43 - fvc::ddt(alpha1) - fvc::div(alphaPhi1)
44 + (alpha1/rho1)
45 *correction
46 (
47 psi1*fvm::ddt(p_rgh)
48 + fvm::div(phid1, p_rgh) - fvm::Sp(fvc::div(phid1), p_rgh)
49 )
50 );
51 p_rghEqnComp1.ref().relax();
52
54 pos(alpha2)
55 *(
56 (
57 fvc::ddt(alpha2, rho2) + fvc::div(alphaPhi2*rho2f)
58 - (fvOptions(alpha2, mixture.thermo2().rho())&rho2)
59 )/rho2
60 - fvc::ddt(alpha2) - fvc::div(alphaPhi2)
61 + (alpha2/rho2)
62 *correction
63 (
64 psi2*fvm::ddt(p_rgh)
65 + fvm::div(phid2, p_rgh) - fvm::Sp(fvc::div(phid2), p_rgh)
66 )
67 );
68 p_rghEqnComp2.ref().relax();
69 }
70 else
71 {
72 #include "rhofs.H"
73
75 pos(alpha1)
76 *(
77 (
78 fvc::ddt(alpha1, rho1) + fvc::div(alphaPhi1*rho1f)
79 - (fvOptions(alpha1, mixture.thermo1().rho())&rho1)
80 )/rho1
81 - fvc::ddt(alpha1) - fvc::div(alphaPhi1)
82 + (alpha1*psi1/rho1)*correction(fvm::ddt(p_rgh))
83 ) - surfaceFilm.Srho()/rho1;
84
86 pos(alpha2)
87 *(
88 (
89 fvc::ddt(alpha2, rho2) + fvc::div(alphaPhi2*rho2f)
90 - (fvOptions(alpha2, mixture.thermo2().rho())&rho2)
91 )/rho2
92 - fvc::ddt(alpha2) - fvc::div(alphaPhi2)
93 + (alpha2*psi2/rho2)*correction(fvm::ddt(p_rgh))
94 );
95 }
96
97 // Cache p_rgh prior to solve for density update
98 volScalarField p_rgh_0(p_rgh);
99
100 while (pimple.correctNonOrthogonal())
101 {
102 fvScalarMatrix p_rghEqnIncomp
103 (
104 fvc::div(phiHbyA)
105 - fvm::laplacian(rAUf, p_rgh)
106 );
107
108 solve
109 (
110 p_rghEqnComp1() + p_rghEqnComp2() + p_rghEqnIncomp,
111 mesh.solver(p_rgh.select(pimple.finalInnerIter()))
112 );
113
114 if (pimple.finalNonOrthogonalIter())
115 {
116 p = max(p_rgh + (alpha1*rho1 + alpha2*rho2)*gh, pMin);
117 p_rgh = p - (alpha1*rho1 + alpha2*rho2)*gh;
118
119 dgdt =
120 (
123 );
124
125 phi = phiHbyA + p_rghEqnIncomp.flux();
126
127 U = HbyA
128 + rAU*fvc::reconstruct((phig + p_rghEqnIncomp.flux())/rAUf);
129 U.correctBoundaryConditions();
130 fvOptions.correct(U);
131 }
132 }
133
134 // Update densities from change in p_rgh
135 mixture.thermo1().correctRho(psi1*(p_rgh - p_rgh_0));
136 mixture.thermo2().correctRho(psi2*(p_rgh - p_rgh_0));
137
139
140 // Correct p_rgh for consistency with p and the updated densities
142 p_rgh.correctBoundaryConditions();
143
144 K = 0.5*magSqr(U);
145}
CGAL::Exact_predicates_exact_constructions_kernel K
Y[inertIndex] max(0.0)
volScalarField & p_rgh
fv::options & fvOptions
surfaceScalarField & phi
const surfaceScalarField & ghf
IOMRFZoneList & MRF
const volScalarField & gh
pimpleControl & pimple
const volScalarField & alpha1
const volScalarField & psi2
volScalarField & rho2
const volScalarField & alpha2
const surfaceScalarField & alphaPhi1
const volScalarField & psi1
const dimensionedScalar & pMin
volScalarField & rho1
const surfaceScalarField & alphaPhi2
U
Definition: pEqn.H:72
volScalarField & p
regionModels::surfaceFilmModel & surfaceFilm
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
surfaceScalarField rho2f(fvc::interpolate(rho2))
surfaceScalarField rho1f(fvc::interpolate(rho1))
dynamicFvMesh & mesh
surfaceScalarField rAUf("rAUf", fvc::interpolate(rAU))
tmp< volScalarField > rAU
Definition: initCorrectPhi.H:1
tmp< fvScalarMatrix > p_rghEqnComp1
Definition: pEqn.H:29
volScalarField p_rgh_0(p_rgh)
tmp< fvScalarMatrix > p_rghEqnComp2
Definition: pEqn.H:30
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
CEqn solve()
Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n"<< endl;autoPtr< temperaturePhaseChangeTwoPhaseMixture > mixture
Definition: createFields.H:39