pEqn.H
Go to the documentation of this file.
1 if (!mesh.steady() && !pimple.simpleRho())
2 {
3  rho = thermo.rho();
4 }
5 
6 volScalarField rAU("rAU", 1.0/UEqn.A());
9 if (pimple.nCorrPiso() <= 1)
10 {
11  tUEqn.clear();
12 }
13 
14 surfaceScalarField phig(-rhorAUf*ghf*fvc::snGrad(rho)*mesh.magSf());
15 
17 (
18  "phiHbyA",
19  fvc::flux(rho*HbyA)
20  + MRF.zeroFilter(rhorAUf*fvc::ddtCorr(rho, U, phi))
21 );
22 
23 MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
24 
25 const bool closedVolume = mesh.steady() && adjustPhi(phiHbyA, U, p_rgh);
26 const bool adjustMass = closedVolume && !thermo.incompressible();
27 
29 
30 // Update the pressure BCs to ensure flux consistency
31 constrainPressure(p_rgh, rho, U, phiHbyA, rhorAUf, MRF);
32 
33 {
34  fvScalarMatrix p_rghEqnComp
35  (
37  ==
38  fvOptions(psi, p_rgh, rho.name())
39  );
40 
41  if (pimple.transonic())
42  {
44  (
45  "phid",
47  );
48 
49  phiHbyA -= fvc::interpolate(psi*p_rgh)*phiHbyA/fvc::interpolate(rho);
50 
51  p_rghEqnComp += fvm::div(phid, p_rgh);
52  }
53 
54  // Thermodynamic density needs to be updated by psi*d(p) after the
55  // pressure solution
56  tmp<volScalarField> psip0(mesh.steady() ? tmp<volScalarField>() : psi*p);
57 
58  while (pimple.correctNonOrthogonal())
59  {
60  fvScalarMatrix p_rghEqnIncomp
61  (
63  - fvm::laplacian(rhorAUf, p_rgh)
64  );
65 
66  fvScalarMatrix p_rghEqn(p_rghEqnComp + p_rghEqnIncomp);
67 
68  p_rghEqn.setReference
69  (
70  pressureControl.refCell(),
71  pressureControl.refValue()
72  );
73 
74  p_rghEqn.solve();
75 
76  if (pimple.finalNonOrthogonalIter())
77  {
78  // Calculate the conservative fluxes
79  phi = phiHbyA + p_rghEqn.flux();
80 
81  // Explicitly relax pressure for momentum corrector
82  p_rgh.relax();
83 
84  // Correct the momentum source with the pressure gradient flux
85  // calculated from the relaxed pressure
86  U = HbyA
87  + rAU*fvc::reconstruct((phig + p_rghEqnIncomp.flux())/rhorAUf);
88  U.correctBoundaryConditions();
89  fvOptions.correct(U);
90  K = 0.5*magSqr(U);
91  }
92  }
93 
94  p = p_rgh + rho*gh;
95 
96  pressureControl.limit(p);
97 
98  // Thermodynamic density update
99  if (!mesh.steady())
100  {
101  thermo.correctRho(psi*p - psip0);
102  }
103 }
104 
105 // Update pressure time derivative if needed
106 if (thermo.dpdt())
107 {
108  dpdt = fvc::ddt(p);
109 }
110 
111 // Solve continuity
112 if (!mesh.steady())
113 {
114  #include "rhoEqn.H"
116 }
117 else
118 {
120 }
121 
122 
123 // For closed-volume compressible cases adjust the pressure level
124 // to obey overall mass continuity
125 if (adjustMass)
126 {
127  p += (initialMass - fvc::domainIntegrate(thermo.rho()))
129  p_rgh = p - rho*gh;
130  p.correctBoundaryConditions();
131 }
132 
133 // Density updates
134 if (adjustMass || mesh.steady() || pimple.simpleRho())
135 {
136  rho = thermo.rho();
137 }
138 
139 if (mesh.steady() && !pimple.transonic())
140 {
141  rho.relax();
142 }
143 
144 Info<< "Min/max rho:" << min(rho).value() << ' '
145  << max(rho).value() << endl;
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
Return the correction form of the given matrix.
rAU
Definition: pEqn.H:1
phiHbyA
Definition: pEqn.H:22
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
fv::options & fvOptions
rho
Definition: pEqn.H:1
pimpleNoLoopControl & pimple
dimensionedScalar initialMass
Definition: createFields.H:57
IOMRFZoneList & MRF
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
surfaceScalarField phig("phig", -rhorAUf *ghf *fvc::snGrad(rho) *mesh.magSf())
tmp< GeometricField< typename flux< Type >::type, fvsPatchField, surfaceMesh > > ddtCorr(const GeometricField< Type, fvPatchField, volMesh > &U, const GeometricField< Type, fvsPatchField, surfaceMesh > &Uf)
Definition: fvcDdt.C:170
p_rgh
Definition: pEqn.H:140
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
tmp< fvVectorMatrix > tUEqn(fvm::ddt(rho, U)+fvm::div(phi, U)+MRF.DDt(rho, U)+turbulence->divDevRhoReff(U)==fvOptions(rho, U))
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: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< 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
volScalarField & dpdt
dynamicFvMesh & mesh
pressureControl & pressureControl
tmp< volVectorField > constrainHbyA(const tmp< volVectorField > &tHbyA, const volVectorField &U, const volScalarField &p)
Definition: constrainHbyA.C:33
dimensioned< scalar > magSqr(const dimensioned< Type > &)
const surfaceScalarField & ghf
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
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.
adjustPhi(phiHbyA, Urel, p)
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> reconstruct(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
U
Definition: pEqn.H:72
volVectorField & HbyA
Definition: pEqn.H:13
fvVectorMatrix & UEqn
Definition: UEqn.H:13
Calculates and prints the continuity errors.
const volScalarField & gh
phi
Definition: pEqn.H:18
messageStream Info
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
const volScalarField psip0(psi *p)