All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
pEqn.H
Go to the documentation of this file.
1 if (simple.consistent())
2 {
3  rho = thermo.rho();
4 }
5 
6 const volScalarField rAU("rAU", 1.0/UEqn.A());
7 const surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho*rAU));
8 
9 tmp<volScalarField> rAtU
10 (
11  simple.consistent()
12  ? volScalarField::New("rAtU", 1.0/(1.0/rAU - UEqn.H1()))
13  : tmp<volScalarField>(nullptr)
14 );
15 tmp<surfaceScalarField> rhorAtUf
16 (
17  simple.consistent()
19  : tmp<surfaceScalarField>(nullptr)
20 );
21 
22 const volScalarField& rAAtU = simple.consistent() ? rAtU() : rAU;
23 const surfaceScalarField& rhorAAtUf =
24  simple.consistent() ? rhorAtUf() : rhorAUf;
25 
27 
28 tUEqn.clear();
29 
31 (
32  "phiHbyA",
34 );
35 
36 MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
37 
38 bool closedVolume = false;
39 
40 // Update the pressure BCs to ensure flux consistency
41 constrainPressure(p, rho, U, phiHbyA, rhorAAtUf, MRF);
42 
43 if (simple.transonic())
44 {
46  (
47  "phid",
49  );
50 
52 
53  if (simple.consistent())
54  {
55  phiHbyA += (rhorAAtUf - rhorAUf)*fvc::snGrad(p)*mesh.magSf();
56  HbyA += (rAAtU - rAU)*fvc::grad(p);
57  }
58 
59  while (simple.correctNonOrthogonal())
60  {
61  fvScalarMatrix pEqn
62  (
64  + fvm::div(phid, p)
65  - fvm::laplacian(rhorAAtUf, p)
66  ==
67  fvOptions(psi, p, rho.name())
68  );
69 
70  // Relax the pressure equation to ensure diagonal-dominance
71  pEqn.relax();
72 
73  pEqn.setReference
74  (
75  pressureControl.refCell(),
76  pressureControl.refValue()
77  );
78 
79  pEqn.solve();
80 
81  if (simple.finalNonOrthogonalIter())
82  {
83  phi = phiHbyA + pEqn.flux();
84  }
85  }
86 }
87 else
88 {
89  closedVolume = adjustPhi(phiHbyA, U, p);
90 
91  if (simple.consistent())
92  {
93  phiHbyA += (rhorAAtUf - rhorAUf)*fvc::snGrad(p)*mesh.magSf();
94  HbyA += (rAAtU - rAU)*fvc::grad(p);
95  }
96 
97  while (simple.correctNonOrthogonal())
98  {
99  fvScalarMatrix pEqn
100  (
102  - fvm::laplacian(rhorAAtUf, p)
103  ==
104  fvOptions(psi, p, rho.name())
105  );
106 
107  pEqn.setReference
108  (
109  pressureControl.refCell(),
110  pressureControl.refValue()
111  );
112 
113  pEqn.solve();
114 
115  if (simple.finalNonOrthogonalIter())
116  {
117  phi = phiHbyA + pEqn.flux();
118  }
119  }
120 }
121 
123 
124 // Explicitly relax pressure for momentum corrector
125 p.relax();
126 
127 U = HbyA - rAAtU*fvc::grad(p);
128 U.correctBoundaryConditions();
129 fvOptions.correct(U);
130 
131 pressureControl.limit(p);
132 
133 // For closed-volume cases adjust the pressure and density levels
134 // to obey overall mass continuity
135 if (closedVolume && !thermo.incompressible())
136 {
139  p.correctBoundaryConditions();
140 }
141 
142 rho = thermo.rho();
143 
144 if (!simple.transonic())
145 {
146  rho.relax();
147 }
tmp< surfaceScalarField > rhorAtUf(pimple.consistent() ? surfaceScalarField::New("rhoRAtUf", fvc::interpolate(rho *rAtU())) :tmp< surfaceScalarField >(nullptr))
rAU
Definition: pEqn.H:1
tmp< volScalarField > rAtU(pimple.consistent() ? volScalarField::New("rAtU", 1.0/(1.0/rAU - UEqn.H1())) :tmp< volScalarField >(nullptr))
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
fv::options & fvOptions
rho
Definition: pEqn.H:1
dimensionedScalar initialMass
Definition: createFields.H:69
p
Definition: pEqn.H:50
IOMRFZoneList & MRF
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
Calculates and prints the continuity errors.
dimensioned< Type > domainIntegrate(const GeometricField< Type, fvPatchField, volMesh > &vf)
rhoReactionThermo & thermo
Definition: createFields.H:28
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:58
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
phi
Definition: pEqn.H:104
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
dynamicFvMesh & mesh
autoPtr< BasicCompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleMomentumTransportModel::transportModel &transport)
pressureControl & pressureControl
tmp< volVectorField > constrainHbyA(const tmp< volVectorField > &tHbyA, const volVectorField &U, const volScalarField &p)
Definition: constrainHbyA.C:34
tmp< fvVectorMatrix > tUEqn(fvm::ddt(rho, U)+fvm::div(phi, U)+MRF.DDt(rho, U)+turbulence->divDevTau(U)==fvOptions(rho, U))
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
phiHbyA
Definition: pEqn.H:32
adjustPhi(phiHbyA, Urel, p)
U
Definition: pEqn.H:72
volVectorField & HbyA
Definition: pEqn.H:13
fvVectorMatrix & UEqn
Definition: UEqn.H:11
constrainPressure(p_rgh, rho, U, phiHbyA, rhorAUf, MRF)
const volScalarField & psi
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
surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho *rAU))
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcSnGrad.C:45
simpleControl simple(mesh)