pEqn.H
Go to the documentation of this file.
1 {
2  volScalarField rAU("rAU", 1.0/UEqn.A());
6  (
7  "phiHbyA",
10  );
11 
13  (
14  (
15  interface.surfaceTensionForce()
17  )*rAUf*mesh.magSf()
18  );
19 
21 
22  // Update the pressure BCs to ensure flux consistency
24 
25  tmp<fvScalarMatrix> p_rghEqnComp1;
26  tmp<fvScalarMatrix> p_rghEqnComp2;
27 
28  if (pimple.transonic())
29  {
31  surfaceScalarField phid2("phid2", fvc::interpolate(psi2)*phi);
32 
33  p_rghEqnComp1 =
34  fvc::ddt(rho1) + fvc::div(phi, rho1) - fvc::Sp(fvc::div(phi), rho1)
35  + correction
36  (
38  + fvm::div(phid1, p_rgh) - fvm::Sp(fvc::div(phid1), p_rgh)
39  );
40  deleteDemandDrivenData(p_rghEqnComp1.ref().faceFluxCorrectionPtr());
41  p_rghEqnComp1.ref().relax();
42 
43  p_rghEqnComp2 =
44  fvc::ddt(rho2) + fvc::div(phi, rho2) - fvc::Sp(fvc::div(phi), rho2)
45  + correction
46  (
48  + fvm::div(phid2, p_rgh) - fvm::Sp(fvc::div(phid2), p_rgh)
49  );
50  deleteDemandDrivenData(p_rghEqnComp2.ref().faceFluxCorrectionPtr());
51  p_rghEqnComp2.ref().relax();
52  }
53  else
54  {
55  p_rghEqnComp1 =
58 
59  p_rghEqnComp2 =
61  + fvc::div(phi, rho2) - fvc::Sp(fvc::div(phi), rho2);
62  }
63 
64  // Cache p_rgh prior to solve for density update
66 
67  while (pimple.correctNonOrthogonal())
68  {
69  fvScalarMatrix p_rghEqnIncomp
70  (
73  );
74 
75  solve
76  (
77  (
78  (max(alpha1, scalar(0))/rho1)*p_rghEqnComp1()
79  + (max(alpha2, scalar(0))/rho2)*p_rghEqnComp2()
80  )
81  + p_rghEqnIncomp,
82  mesh.solver(p_rgh.select(pimple.finalInnerIter()))
83  );
84 
85  if (pimple.finalNonOrthogonalIter())
86  {
87  p = max(p_rgh + (alpha1*rho1 + alpha2*rho2)*gh, pMin);
88  p_rgh = p - (alpha1*rho1 + alpha2*rho2)*gh;
89 
90  dgdt =
91  (
92  pos(alpha2)*(p_rghEqnComp2 & p_rgh)/rho2
93  - pos(alpha1)*(p_rghEqnComp1 & p_rgh)/rho1
94  );
95 
96  phi = phiHbyA + p_rghEqnIncomp.flux();
97 
98  U = HbyA
99  + rAU*fvc::reconstruct((phig + p_rghEqnIncomp.flux())/rAUf);
100  U.correctBoundaryConditions();
101  }
102  }
103 
104  // Update densities from change in p_rgh
106  rho2 += psi2*(p_rgh - p_rgh_0);
107 
109 
110  // Correct p_rgh for consistency with p and the updated densities
111  p_rgh = p - rho*gh;
112  p_rgh.correctBoundaryConditions();
113 
114  K = 0.5*magSqr(U);
115 
116  Info<< "max(U) " << max(mag(U)).value() << endl;
117  Info<< "min(p_rgh) " << min(p_rgh).value() << endl;
118 }
tmp< fvScalarMatrix > p_rghEqnComp2
Definition: pEqn.H:29
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
Return the correction form of the given matrix.
psi2
Definition: TEqns.H:35
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
U
Definition: pEqn.H:83
p
Definition: pEqn.H:50
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
interfaceProperties interface(alpha1, U, mixture())
phiHbyA
Definition: pEqn.H:20
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
tmp< GeometricField< typename flux< Type >::type, fvsPatchField, surfaceMesh > > ddtCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)
Definition: fvcDdt.C:155
const dictionary & pimple
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
tmp< GeometricField< Type, fvPatchField, volMesh > > Sp(const volScalarField &sp, const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvcSup.C:67
tmp< surfaceScalarField > interpolate(const RhoType &rho)
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
CGAL::Exact_predicates_exact_constructions_kernel K
tmp< fvScalarMatrix > p_rghEqnComp1
Definition: pEqn.H:28
rho1
Definition: pEqn.H:114
const surfaceScalarField & ghf
tmp< GeometricField< Type, fvPatchField, volMesh > > ddt(const dimensioned< Type > dt, const fvMesh &mesh)
Definition: fvcDdt.C:45
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
dimensionedScalar pos(const dimensionedScalar &ds)
dynamicFvMesh & mesh
rhoEqn solve()
rho2
Definition: pEqn.H:115
surfaceScalarField phig("phig",-rhorAUf *ghf *fvc::snGrad(rho)*mesh.magSf())
tmp< volVectorField > constrainHbyA(const tmp< volVectorField > &tHbyA, const volVectorField &U, const volScalarField &p)
Definition: constrainHbyA.C:33
alpha2
Definition: alphaEqn.H:112
p_rgh
Definition: pEqn.H:120
dimensioned< scalar > magSqr(const dimensioned< Type > &)
constrainPressure(p, rho, U, phiHbyA, rhorAUf, MRF)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
volScalarField & alpha1
dimensionedScalar pMin("pMin", dimPressure, fluid)
surfaceScalarField rAUf("rAUf", fvc::interpolate(rAU))
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> reconstruct(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
psi1
Definition: TEqns.H:34
volVectorField & HbyA
Definition: pEqn.H:13
fvVectorMatrix & UEqn
Definition: UEqn.H:13
const volScalarField & gh
volScalarField p_rgh_0(p_rgh)
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
phi
Definition: pEqn.H:18
tmp< surfaceScalarField > flux(const volVectorField &vvf)
Return the face-flux field obtained from the given volVectorField.
Definition: fvcFlux.C:32
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
rho
Definition: pEqn.H:1
void deleteDemandDrivenData(DataPtr &dataPtr)
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcSnGrad.C:45
volScalarField rAU(1.0/UEqn.A())