pcEqn.H
Go to the documentation of this file.
1 volScalarField rAU(1.0/UEqn().A());
2 volScalarField rAtU(1.0/(1.0/rAU - UEqn().H1()));
3 
4 volVectorField HbyA("HbyA", U);
5 HbyA = rAU*UEqn().H();
6 
7 UEqn.clear();
8 
9 bool closedVolume = false;
10 
11 if (simple.transonic())
12 {
14  (
15  "phid",
17  *(fvc::interpolate(HbyA) & mesh.Sf())
18  );
19 
20  MRF.makeRelative(fvc::interpolate(psi), phid);
21 
23  (
24  "phic",
26  );
27 
28  HbyA -= (rAU - rAtU)*fvc::grad(p);
29 
30  volScalarField rhorAtU("rhorAtU", rho*rAtU);
31 
32  while (simple.correctNonOrthogonal())
33  {
34  fvScalarMatrix pEqn
35  (
36  fvm::div(phid, p)
37  + fvc::div(phic)
38  - fvm::laplacian(rhorAtU, p)
39  ==
40  fvOptions(psi, p, rho.name())
41  );
42 
43  // Relax the pressure equation to maintain diagonal dominance
44  pEqn.relax();
45 
46  pEqn.setReference(pRefCell, pRefValue);
47 
48  pEqn.solve();
49 
50  if (simple.finalNonOrthogonalIter())
51  {
52  phi == phic + pEqn.flux();
53  }
54  }
55 }
56 else
57 {
59  (
60  "phiHbyA",
61  fvc::interpolate(rho*HbyA) & mesh.Sf()
62  );
63 
64  MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
65 
66  closedVolume = adjustPhi(phiHbyA, U, p);
67 
69  HbyA -= (rAU - rAtU)*fvc::grad(p);
70 
71  volScalarField rhorAtU("rhorAtU", rho*rAtU);
72 
73  while (simple.correctNonOrthogonal())
74  {
75  fvScalarMatrix pEqn
76  (
79  ==
80  fvOptions(psi, p, rho.name())
81  );
82 
83  pEqn.setReference(pRefCell, pRefValue);
84 
85  pEqn.solve();
86 
87  if (simple.finalNonOrthogonalIter())
88  {
89  phi = phiHbyA + pEqn.flux();
90  }
91  }
92 }
93 
94 // The incompressibe form of the continuity error check is appropriate for
95 // steady-state compressible also.
97 
98 // Explicitly relax pressure for momentum corrector
99 p.relax();
100 
101 U = HbyA - rAtU*fvc::grad(p);
102 U.correctBoundaryConditions();
103 fvOptions.correct(U);
104 
105 // For closed-volume cases adjust the pressure and density levels
106 // to obey overall mass continuity
107 if (closedVolume)
108 {
111 }
112 
113 // Recalculate density from the relaxed pressure
114 rho = thermo.rho();
115 rho = max(rho, rhoMin);
116 rho = min(rho, rhoMax);
117 
118 if (!simple.transonic())
119 {
120  rho.relax();
121 }
122 
123 Info<< "rho max/min : " << max(rho).value() << " " << min(rho).value() << endl;
rho
Definition: pcEqn.H:1
volScalarField rAU(1.0/UEqn().A())
volScalarField rAtU(1.0/(1.0/rAU-UEqn().H1()))
bool adjustPhi(surfaceScalarField &phi, const volVectorField &U, volScalarField &p)
Adjust the balance of fluxes to obey continuity.
Definition: adjustPhi.C:34
phiHbyA
Definition: pcEqn.H:74
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
messageStream Info
tmp< surfaceScalarField > interpolate(const RhoType &rho)
dynamicFvMesh & mesh
HbyA
Definition: pcEqn.H:9
fv::IOoptionList & fvOptions
const scalar pRefValue
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
IOMRFZoneList & MRF
fvVectorMatrix UEqn(fvm::ddt(rho, U)+fvm::div(phi, U)+turbulence->divDevRhoReff(U)==fvOptions(rho, U))
surfaceScalarField phid("phid", fvc::interpolate(psi)*( (mesh.Sf()&fvc::interpolate(HbyA)) +rhorAUf *fvc::ddtCorr(rho, U, phi)/fvc::interpolate(rho) ))
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
volScalarField & p
Definition: createFields.H:51
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
const volScalarField & psi
Definition: createFields.H:24
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcSnGrad.C:45
surfaceScalarField & phi
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
dimensionedScalar initialMass
Definition: createFields.H:83
PtrList< dimensionedScalar > rhoMin(fluidRegions.size())
dimensioned< Type > domainIntegrate(const GeometricField< Type, fvPatchField, volMesh > &vf)
const label pRefCell
const dictionary & simple
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
psiReactionThermo & thermo
Definition: createFields.H:32
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
volScalarField rhorAtU("rhorAtU", rho *rAtU)
PtrList< dimensionedScalar > rhoMax(fluidRegions.size())
Calculates and prints the continuity errors.
U
Definition: pcEqn.H:105