From 50516486a4c01de0303c866da3961a8de3cb06fa Mon Sep 17 00:00:00 2001
From: Henry Weller <http://cfd.direct>
Date: Tue, 28 Feb 2017 11:14:59 +0000
Subject: [PATCH] rhoPimpleFoam: Added support for transonic flow of liquids
and real gases
Both stardard SIMPLE and the SIMPLEC (using the 'consistent' option in
fvSolution) are now supported for both subsonic and transonic flow of all
fluid types.
rhoPimpleFoam now instantiates the lower-level fluidThermo which instantiates
either a psiThermo or rhoThermo according to the 'type' specification in
thermophysicalProperties, see also commit a1c8cde310f48925143781a940811a5fca8a87c2
---
.../compressible/rhoPimpleFoam/createFields.H | 28 +--
.../solvers/compressible/rhoPimpleFoam/pEqn.H | 65 +++---
.../compressible/rhoPimpleFoam/pcEqn.H | 76 +++----
.../rhoPimpleFoam/rhoPimpleDyMFoam/pEqn.H | 64 +++---
.../rhoPimpleDyMFoam/rhoPimpleDyMFoam.C | 10 +-
.../rhoPimpleFoam/rhoPimpleFoam.C | 5 +-
.../compressible/rhoSimpleFoam/createFields.H | 4 +-
.../solvers/compressible/rhoSimpleFoam/pEqn.H | 18 +-
.../compressible/rhoSimpleFoam/pcEqn.H | 28 +--
.../rhoPorousSimpleFoam/createFields.H | 59 -----
.../general/pressureControl/pressureControl.C | 212 ++++++++++++------
.../general/pressureControl/pressureControl.H | 3 +-
.../annularThermalMixer/system/fvSolution | 5 +-
.../LES/pitzDaily/system/fvSolution | 5 +-
.../constant/thermophysicalProperties | 2 +-
.../RAS/angledDuct/system/fvSchemes | 3 +-
.../RAS/angledDuct/system/fvSolution | 14 +-
.../RAS/angledDuctLTS/system/fvSolution | 4 +-
.../RAS/cavity/system/fvSolution | 5 +-
.../RAS/mixerVessel2D/system/fvSolution | 15 +-
.../helmholtzResonance/system/fvSolution | 5 +-
21 files changed, 291 insertions(+), 339 deletions(-)
delete mode 100644 applications/solvers/compressible/rhoSimpleFoam/rhoPorousSimpleFoam/createFields.H
diff --git a/applications/solvers/compressible/rhoPimpleFoam/createFields.H b/applications/solvers/compressible/rhoPimpleFoam/createFields.H
index a7ee3eca457..84f70f2cd02 100644
--- a/applications/solvers/compressible/rhoPimpleFoam/createFields.H
+++ b/applications/solvers/compressible/rhoPimpleFoam/createFields.H
@@ -2,11 +2,11 @@
Info<< "Reading thermophysical properties\n" << endl;
-autoPtr<psiThermo> pThermo
+autoPtr<fluidThermo> pThermo
(
- psiThermo::New(mesh)
+ fluidThermo::New(mesh)
);
-psiThermo& thermo = pThermo();
+fluidThermo& thermo = pThermo();
thermo.validate(args.executable(), "h", "e");
volScalarField& p = thermo.p();
@@ -40,27 +40,7 @@ volVectorField U
#include "compressibleCreatePhi.H"
-dimensionedScalar rhoMax
-(
- dimensionedScalar::lookupOrDefault
- (
- "rhoMax",
- pimple.dict(),
- dimDensity,
- GREAT
- )
-);
-
-dimensionedScalar rhoMin
-(
- dimensionedScalar::lookupOrDefault
- (
- "rhoMin",
- pimple.dict(),
- dimDensity,
- 0
- )
-);
+pressureControl pressureControl(p, rho, pimple.dict(), false);
Info<< "Creating turbulence model\n" << endl;
autoPtr<compressible::turbulenceModel> turbulence
diff --git a/applications/solvers/compressible/rhoPimpleFoam/pEqn.H b/applications/solvers/compressible/rhoPimpleFoam/pEqn.H
index ac7107acf06..73ecafd89f8 100644
--- a/applications/solvers/compressible/rhoPimpleFoam/pEqn.H
+++ b/applications/solvers/compressible/rhoPimpleFoam/pEqn.H
@@ -1,8 +1,3 @@
-rho = thermo.rho();
-rho = max(rho, rhoMin);
-rho = min(rho, rhoMax);
-rho.relax();
-
volScalarField rAU(1.0/UEqn.A());
surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho*rAU));
volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p));
@@ -12,55 +7,54 @@ if (pimple.nCorrPISO() <= 1)
tUEqn.clear();
}
+surfaceScalarField phiHbyA
+(
+ "phiHbyA",
+ (
+ fvc::flux(rho*HbyA)
+ + rhorAUf*fvc::ddtCorr(rho, U, phi)
+ )
+);
+
+MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
+
+// Update the pressure BCs to ensure flux consistency
+constrainPressure(p, rho, U, phiHbyA, rhorAUf, MRF);
+
if (pimple.transonic())
{
surfaceScalarField phid
(
"phid",
- fvc::interpolate(psi)
- *(
- fvc::flux(HbyA)
- + rhorAUf*fvc::ddtCorr(rho, U, phi)/fvc::interpolate(rho)
- )
+ (fvc::interpolate(psi)/fvc::interpolate(rho))*phiHbyA
);
-
- MRF.makeRelative(fvc::interpolate(psi), phid);
+ phiHbyA -= fvc::interpolate(p)*phid;
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::ddt(psi, p)
+ + fvc::div(phiHbyA)
+ fvm::div(phid, p)
- fvm::laplacian(rhorAUf, p)
==
fvOptions(psi, p, rho.name())
);
+ // Relax the pressure equation to ensure diagonal-dominance
+ pEqn.relax();
+
pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
if (pimple.finalNonOrthogonalIter())
{
- phi == pEqn.flux();
+ phi = phiHbyA + pEqn.flux();
}
}
}
else
{
- surfaceScalarField phiHbyA
- (
- "phiHbyA",
- (
- fvc::flux(rho*HbyA)
- + rhorAUf*fvc::ddtCorr(rho, U, phi)
- )
- );
-
- MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
-
- // Update the pressure BCs to ensure flux consistency
- constrainPressure(p, rho, U, phiHbyA, rhorAUf, MRF);
-
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
@@ -87,19 +81,20 @@ else
// Explicitly relax pressure for momentum corrector
p.relax();
-// Recalculate density from the relaxed pressure
-rho = thermo.rho();
-rho = max(rho, rhoMin);
-rho = min(rho, rhoMax);
-rho.relax();
-Info<< "rho max/min : " << max(rho).value()
- << " " << min(rho).value() << endl;
-
U = HbyA - rAU*fvc::grad(p);
U.correctBoundaryConditions();
fvOptions.correct(U);
K = 0.5*magSqr(U);
+pressureControl.limit(p);
+p.correctBoundaryConditions();
+rho = thermo.rho();
+
+if (!pimple.transonic())
+{
+ rho.relax();
+}
+
if (thermo.dpdt())
{
dpdt = fvc::ddt(p);
diff --git a/applications/solvers/compressible/rhoPimpleFoam/pcEqn.H b/applications/solvers/compressible/rhoPimpleFoam/pcEqn.H
index 713f443fc5d..afbc2851e4a 100644
--- a/applications/solvers/compressible/rhoPimpleFoam/pcEqn.H
+++ b/applications/solvers/compressible/rhoPimpleFoam/pcEqn.H
@@ -1,8 +1,3 @@
-rho = thermo.rho();
-rho = max(rho, rhoMin);
-rho = min(rho, rhoMax);
-rho.relax();
-
volScalarField rAU(1.0/UEqn.A());
volScalarField rAtU(1.0/(1.0/rAU - UEqn.H1()));
volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p));
@@ -12,72 +7,64 @@ if (pimple.nCorrPISO() <= 1)
tUEqn.clear();
}
+surfaceScalarField phiHbyA
+(
+ "phiHbyA",
+ (
+ fvc::flux(rho*HbyA)
+ + fvc::interpolate(rho*rAU)*fvc::ddtCorr(rho, U, phi)
+ )
+);
+
+MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
+
+volScalarField rhorAtU("rhorAtU", rho*rAtU);
+
+// Update the pressure BCs to ensure flux consistency
+constrainPressure(p, rho, U, phiHbyA, rhorAtU, MRF);
+
if (pimple.transonic())
{
surfaceScalarField phid
(
"phid",
- fvc::interpolate(psi)
- *(
- fvc::flux(HbyA)
- + fvc::interpolate(rho*rAU)*fvc::ddtCorr(rho, U, phi)
- /fvc::interpolate(rho)
- )
+ (fvc::interpolate(psi)/fvc::interpolate(rho))*phiHbyA
);
- MRF.makeRelative(fvc::interpolate(psi), phid);
-
- surfaceScalarField phic
- (
- "phic",
+ phiHbyA +=
fvc::interpolate(rho*(rAtU - rAU))*fvc::snGrad(p)*mesh.magSf()
- );
+ - fvc::interpolate(p)*phid;
HbyA -= (rAU - rAtU)*fvc::grad(p);
- volScalarField rhorAtU("rhorAtU", rho*rAtU);
-
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::ddt(psi, p)
+ + fvc::div(phiHbyA)
+ fvm::div(phid, p)
- + fvc::div(phic)
- fvm::laplacian(rhorAtU, p)
==
fvOptions(psi, p, rho.name())
);
+ // Relax the pressure equation to ensure diagonal-dominance
+ pEqn.relax();
+
pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
if (pimple.finalNonOrthogonalIter())
{
- phi == phic + pEqn.flux();
+ phi = phiHbyA + pEqn.flux();
}
}
}
else
{
- surfaceScalarField phiHbyA
- (
- "phiHbyA",
- (
- fvc::flux(rho*HbyA)
- + fvc::interpolate(rho*rAU)*fvc::ddtCorr(rho, U, phi)
- )
- );
-
- MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
-
phiHbyA += fvc::interpolate(rho*(rAtU - rAU))*fvc::snGrad(p)*mesh.magSf();
HbyA -= (rAU - rAtU)*fvc::grad(p);
- volScalarField rhorAtU("rhorAtU", rho*rAtU);
-
- // Update the pressure BCs to ensure flux consistency
- constrainPressure(p, rho, U, phiHbyA, rhorAtU, MRF);
-
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
@@ -109,19 +96,16 @@ U.correctBoundaryConditions();
fvOptions.correct(U);
K = 0.5*magSqr(U);
-if (thermo.dpdt())
-{
- dpdt = fvc::ddt(p);
-}
-
-// Recalculate density from the relaxed pressure
+pressureControl.limit(p);
+p.correctBoundaryConditions();
rho = thermo.rho();
-rho = max(rho, rhoMin);
-rho = min(rho, rhoMax);
if (!pimple.transonic())
{
rho.relax();
}
-Info<< "rho max/min : " << max(rho).value() << " " << min(rho).value() << endl;
+if (thermo.dpdt())
+{
+ dpdt = fvc::ddt(p);
+}
diff --git a/applications/solvers/compressible/rhoPimpleFoam/rhoPimpleDyMFoam/pEqn.H b/applications/solvers/compressible/rhoPimpleFoam/rhoPimpleDyMFoam/pEqn.H
index 7bd540df40d..0f855626189 100644
--- a/applications/solvers/compressible/rhoPimpleFoam/rhoPimpleDyMFoam/pEqn.H
+++ b/applications/solvers/compressible/rhoPimpleFoam/rhoPimpleDyMFoam/pEqn.H
@@ -1,8 +1,3 @@
-rho = thermo.rho();
-rho = max(rho, rhoMin);
-rho = min(rho, rhoMax);
-rho.relax();
-
volScalarField rAU(1.0/UEqn.A());
surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho*rAU));
volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p));
@@ -12,55 +7,53 @@ if (pimple.nCorrPISO() <= 1)
tUEqn.clear();
}
+surfaceScalarField phiHbyA
+(
+ "phiHbyA",
+ fvc::flux(rho*HbyA)
+ + rhorAUf*fvc::ddtCorr(rho, U, rhoUf)
+);
+
+fvc::makeRelative(phiHbyA, rho, U);
+MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
+
+// Update the pressure BCs to ensure flux consistency
+constrainPressure(p, rho, U, phiHbyA, rhorAUf, MRF);
+
if (pimple.transonic())
{
surfaceScalarField phid
(
"phid",
- fvc::interpolate(psi)
- *(
- fvc::flux(HbyA)
- + rhorAUf*fvc::ddtCorr(rho, U, rhoUf)/fvc::interpolate(rho)
- )
+ (fvc::interpolate(psi)/fvc::interpolate(rho))*phiHbyA
);
-
- fvc::makeRelative(phid, psi, U);
- MRF.makeRelative(fvc::interpolate(psi), phid);
+ phiHbyA -= fvc::interpolate(p)*phid;
while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::ddt(psi, p)
+ + fvc::div(phiHbyA)
+ fvm::div(phid, p)
- fvm::laplacian(rhorAUf, p)
==
fvOptions(psi, p, rho.name())
);
+ // Relax the pressure equation to ensure diagonal-dominance
+ pEqn.relax();
+
pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
if (pimple.finalNonOrthogonalIter())
{
- phi == pEqn.flux();
+ phi = phiHbyA + pEqn.flux();
}
}
}
else
{
- surfaceScalarField phiHbyA
- (
- "phiHbyA",
- fvc::flux(rho*HbyA)
- + rhorAUf*fvc::ddtCorr(rho, U, rhoUf)
- );
-
- fvc::makeRelative(phiHbyA, rho, U);
- MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
-
- // Update the pressure BCs to ensure flux consistency
- constrainPressure(p, rho, U, phiHbyA, rhorAUf, MRF);
-
while (pimple.correctNonOrthogonal())
{
// Pressure corrector
@@ -88,19 +81,20 @@ else
// Explicitly relax pressure for momentum corrector
p.relax();
-// Recalculate density from the relaxed pressure
-rho = thermo.rho();
-rho = max(rho, rhoMin);
-rho = min(rho, rhoMax);
-rho.relax();
-Info<< "rho max/min : " << max(rho).value()
- << " " << min(rho).value() << endl;
-
U = HbyA - rAU*fvc::grad(p);
U.correctBoundaryConditions();
fvOptions.correct(U);
K = 0.5*magSqr(U);
+pressureControl.limit(p);
+p.correctBoundaryConditions();
+rho = thermo.rho();
+
+if (!pimple.transonic())
+{
+ rho.relax();
+}
+
{
rhoUf = fvc::interpolate(rho*U);
surfaceVectorField n(mesh.Sf()/mesh.magSf());
diff --git a/applications/solvers/compressible/rhoPimpleFoam/rhoPimpleDyMFoam/rhoPimpleDyMFoam.C b/applications/solvers/compressible/rhoPimpleFoam/rhoPimpleDyMFoam/rhoPimpleDyMFoam.C
index c7046bb9b40..1d20c3f10f7 100644
--- a/applications/solvers/compressible/rhoPimpleFoam/rhoPimpleDyMFoam/rhoPimpleDyMFoam.C
+++ b/applications/solvers/compressible/rhoPimpleFoam/rhoPimpleDyMFoam/rhoPimpleDyMFoam.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -22,10 +22,7 @@ License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Application
- rhoPimpleFoam
-
-Group
- grpCompressibleSolvers grpMovingMeshSolvers
+ rhoPimpleDyMFoam
Description
Transient solver for turbulent flow of compressible fluids for HVAC and
@@ -38,10 +35,11 @@ Description
#include "fvCFD.H"
#include "dynamicFvMesh.H"
-#include "psiThermo.H"
+#include "fluidThermo.H"
#include "turbulentFluidThermoModel.H"
#include "bound.H"
#include "pimpleControl.H"
+#include "pressureControl.H"
#include "CorrectPhi.H"
#include "fvOptions.H"
#include "localEulerDdtScheme.H"
diff --git a/applications/solvers/compressible/rhoPimpleFoam/rhoPimpleFoam.C b/applications/solvers/compressible/rhoPimpleFoam/rhoPimpleFoam.C
index 108ad9aa0e4..f30dae76a3b 100644
--- a/applications/solvers/compressible/rhoPimpleFoam/rhoPimpleFoam.C
+++ b/applications/solvers/compressible/rhoPimpleFoam/rhoPimpleFoam.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2017 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -34,10 +34,11 @@ Description
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
-#include "psiThermo.H"
+#include "fluidThermo.H"
#include "turbulentFluidThermoModel.H"
#include "bound.H"
#include "pimpleControl.H"
+#include "pressureControl.H"
#include "fvOptions.H"
#include "localEulerDdtScheme.H"
#include "fvcSmooth.H"
diff --git a/applications/solvers/compressible/rhoSimpleFoam/createFields.H b/applications/solvers/compressible/rhoSimpleFoam/createFields.H
index 4671347b66a..c7914f89d33 100644
--- a/applications/solvers/compressible/rhoSimpleFoam/createFields.H
+++ b/applications/solvers/compressible/rhoSimpleFoam/createFields.H
@@ -7,6 +7,8 @@ autoPtr<fluidThermo> pThermo
fluidThermo& thermo = pThermo();
thermo.validate(args.executable(), "h", "e");
+volScalarField& p = thermo.p();
+
volScalarField rho
(
IOobject
@@ -20,8 +22,6 @@ volScalarField rho
thermo.rho()
);
-volScalarField& p = thermo.p();
-
Info<< "Reading field U\n" << endl;
volVectorField U
(
diff --git a/applications/solvers/compressible/rhoSimpleFoam/pEqn.H b/applications/solvers/compressible/rhoSimpleFoam/pEqn.H
index 5f092ada2f6..c7edbb9e877 100644
--- a/applications/solvers/compressible/rhoSimpleFoam/pEqn.H
+++ b/applications/solvers/compressible/rhoSimpleFoam/pEqn.H
@@ -6,16 +6,18 @@
bool closedVolume = false;
+ surfaceScalarField phiHbyA("phiHbyA", fvc::flux(rho*HbyA));
+ MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
+
+ // Update the pressure BCs to ensure flux consistency
+ constrainPressure(p, rho, U, phiHbyA, rhorAUf, MRF);
+
if (simple.transonic())
{
- surfaceScalarField phiHbyA("phiHbyA", fvc::flux(rho*HbyA));
- surfaceScalarField rhof(fvc::interpolate(rho));
- MRF.makeRelative(rhof, phiHbyA);
-
surfaceScalarField phid
(
"phid",
- (fvc::interpolate(psi)/rhof)*phiHbyA
+ (fvc::interpolate(psi)/fvc::interpolate(rho))*phiHbyA
);
phiHbyA -= fvc::interpolate(p)*phid;
@@ -49,14 +51,8 @@
}
else
{
- surfaceScalarField phiHbyA("phiHbyA", fvc::flux(rho*HbyA));
- MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
-
closedVolume = adjustPhi(phiHbyA, U, p);
- // Update the pressure BCs to ensure flux consistency
- constrainPressure(p, rho, U, phiHbyA, rhorAUf, MRF);
-
while (simple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
diff --git a/applications/solvers/compressible/rhoSimpleFoam/pcEqn.H b/applications/solvers/compressible/rhoSimpleFoam/pcEqn.H
index af2005febd4..c7dc0f864d2 100644
--- a/applications/solvers/compressible/rhoSimpleFoam/pcEqn.H
+++ b/applications/solvers/compressible/rhoSimpleFoam/pcEqn.H
@@ -5,16 +5,20 @@ tUEqn.clear();
bool closedVolume = false;
+surfaceScalarField phiHbyA("phiHbyA", fvc::flux(rho*HbyA));
+MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
+
+volScalarField rhorAtU("rhorAtU", rho*rAtU);
+
+// Update the pressure BCs to ensure flux consistency
+constrainPressure(p, rho, U, phiHbyA, rhorAtU, MRF);
+
if (simple.transonic())
{
- surfaceScalarField phiHbyA("phiHbyA", fvc::flux(rho*HbyA));
- surfaceScalarField rhof(fvc::interpolate(rho));
- MRF.makeRelative(rhof, phiHbyA);
-
surfaceScalarField phid
(
"phid",
- (fvc::interpolate(psi)/rhof)*phiHbyA
+ (fvc::interpolate(psi)/fvc::interpolate(rho))*phiHbyA
);
phiHbyA +=
@@ -23,14 +27,12 @@ if (simple.transonic())
HbyA -= (rAU - rAtU)*fvc::grad(p);
- volScalarField rhorAtU("rhorAtU", rho*rAtU);
-
while (simple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
- fvm::div(phid, p)
- + fvc::div(phiHbyA)
+ fvc::div(phiHbyA)
+ + fvm::div(phid, p)
- fvm::laplacian(rhorAtU, p)
==
fvOptions(psi, p, rho.name())
@@ -55,19 +57,11 @@ if (simple.transonic())
}
else
{
- surfaceScalarField phiHbyA("phiHbyA", fvc::flux(rho*HbyA));
- MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
-
closedVolume = adjustPhi(phiHbyA, U, p);
phiHbyA += fvc::interpolate(rho*(rAtU - rAU))*fvc::snGrad(p)*mesh.magSf();
HbyA -= (rAU - rAtU)*fvc::grad(p);
- volScalarField rhorAtU("rhorAtU", rho*rAtU);
-
- // Update the pressure BCs to ensure flux consistency
- constrainPressure(p, rho, U, phiHbyA, rhorAtU, MRF);
-
while (simple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
diff --git a/applications/solvers/compressible/rhoSimpleFoam/rhoPorousSimpleFoam/createFields.H b/applications/solvers/compressible/rhoSimpleFoam/rhoPorousSimpleFoam/createFields.H
deleted file mode 100644
index 4671347b66a..00000000000
--- a/applications/solvers/compressible/rhoSimpleFoam/rhoPorousSimpleFoam/createFields.H
+++ /dev/null
@@ -1,59 +0,0 @@
-Info<< "Reading thermophysical properties\n" << endl;
-
-autoPtr<fluidThermo> pThermo
-(
- fluidThermo::New(mesh)
-);
-fluidThermo& thermo = pThermo();
-thermo.validate(args.executable(), "h", "e");
-
-volScalarField rho
-(
- IOobject
- (
- "rho",
- runTime.timeName(),
- mesh,
- IOobject::READ_IF_PRESENT,
- IOobject::AUTO_WRITE
- ),
- thermo.rho()
-);
-
-volScalarField& p = thermo.p();
-
-Info<< "Reading field U\n" << endl;
-volVectorField U
-(
- IOobject
- (
- "U",
- runTime.timeName(),
- mesh,
- IOobject::MUST_READ,
- IOobject::AUTO_WRITE
- ),
- mesh
-);
-
-#include "compressibleCreatePhi.H"
-
-pressureControl pressureControl(p, rho, simple.dict());
-
-mesh.setFluxRequired(p.name());
-
-Info<< "Creating turbulence model\n" << endl;
-autoPtr<compressible::turbulenceModel> turbulence
-(
- compressible::turbulenceModel::New
- (
- rho,
- U,
- phi,
- thermo
- )
-);
-
-dimensionedScalar initialMass = fvc::domainIntegrate(rho);
-
-#include "createMRF.H"
diff --git a/src/finiteVolume/cfdTools/general/pressureControl/pressureControl.C b/src/finiteVolume/cfdTools/general/pressureControl/pressureControl.C
index 101126a5e87..aca22e76990 100644
--- a/src/finiteVolume/cfdTools/general/pressureControl/pressureControl.C
+++ b/src/finiteVolume/cfdTools/general/pressureControl/pressureControl.C
@@ -32,105 +32,171 @@ Foam::pressureControl::pressureControl
(
const volScalarField& p,
const volScalarField& rho,
- const dictionary& dict
+ const dictionary& dict,
+ const bool pRefRequired
)
:
- refCell_(0),
+ refCell_(-1),
refValue_(0),
pMax_("pMax", dimPressure, 0),
pMin_("pMin", dimPressure, GREAT)
{
- if (setRefCell(p, dict, refCell_, refValue_))
- {
- pMax_.value() = refValue_;
- pMin_.value() = refValue_;
- }
-
- const volScalarField::Boundary& pbf = p.boundaryField();
- const volScalarField::Boundary& rhobf = rho.boundaryField();
+ bool pLimits = false;
- scalar rhoRefMax = -GREAT;
- scalar rhoRefMin = GREAT;
- bool rhoLimits = false;
-
- forAll(pbf, patchi)
+ // Set the reference cell and value for closed domain simulations
+ if (pRefRequired && setRefCell(p, dict, refCell_, refValue_))
{
- if (pbf[patchi].fixesValue())
- {
- rhoLimits = true;
+ pLimits = true;
- pMax_.value() = max(pMax_.value(), max(pbf[patchi]));
- pMin_.value() = min(pMin_.value(), min(pbf[patchi]));
-
- rhoRefMax = max(rhoRefMax, max(rhobf[patchi]));
- rhoRefMin = min(rhoRefMin, min(rhobf[patchi]));
- }
+ pMax_.value() = refValue_;
+ pMin_.value() = refValue_;
}
- if (dict.found("pMax"))
+ if (dict.found("pMax") && dict.found("pMin"))
{
pMax_.value() = readScalar(dict.lookup("pMax"));
+ pMin_.value() = readScalar(dict.lookup("pMin"));
}
- else if (dict.found("pMaxFactor"))
- {
- const scalar pMaxFactor(readScalar(dict.lookup("pMaxFactor")));
- pMax_ *= pMaxFactor;
- }
- else if (dict.found("rhoMax"))
+ else
{
- // For backward-compatibility infer the pMax from rhoMax
+ const volScalarField::Boundary& pbf = p.boundaryField();
+ const volScalarField::Boundary& rhobf = rho.boundaryField();
- IOWarningInFunction(dict)
- << "'rhoMax' specified rather than 'pMax' or 'pMaxFactor'" << nl
- << " This is supported for backward-compatibility but "
- "'pMax' or 'pMaxFactor' are more reliable." << endl;
+ scalar rhoRefMax = -GREAT;
+ scalar rhoRefMin = GREAT;
+ bool rhoLimits = false;
- if (!rhoLimits)
+ forAll(pbf, patchi)
{
- FatalIOErrorInFunction(dict)
- << "'rhoMax' specified rather than 'pMaxFactor'" << nl
- << " but the corresponding reference density cannot"
- " be evaluated from the boundary conditions." << nl
- << "Please specify 'pMaxFactor' rather than 'rhoMax'"
- << exit(FatalError);
- }
+ if (pbf[patchi].fixesValue())
+ {
+ pLimits = true;
+ rhoLimits = true;
- dimensionedScalar rhoMax("rhoMax", dimDensity, dict);
+ pMax_.value() = max(pMax_.value(), max(pbf[patchi]));
+ pMin_.value() = min(pMin_.value(), min(pbf[patchi]));
- pMax_ *= max(rhoMax.value()/rhoRefMax, 1);
- }
+ rhoRefMax = max(rhoRefMax, max(rhobf[patchi]));
+ rhoRefMin = min(rhoRefMin, min(rhobf[patchi]));
+ }
+ }
- if (dict.found("pMin"))
- {
- pMin_.value() = readScalar(dict.lookup("pMin"));
- }
- else if (dict.found("pMinFactor"))
- {
- const scalar pMinFactor(readScalar(dict.lookup("pMinFactor")));
- pMin_ *= pMinFactor;
- }
- else if (dict.found("rhoMin"))
- {
- // For backward-compatibility infer the pMin from rhoMin
+ reduce(rhoLimits, andOp<bool>());
+ if (rhoLimits)
+ {
+ reduce(pMax_.value(), maxOp<scalar>());
+ reduce(pMin_.value(), minOp<scalar>());
- IOWarningInFunction(dict)
- << "'rhoMin' specified rather than 'pMin' or 'pMinFactor'" << nl
- << " This is supported for backward-compatibility but"
- "'pMin' or 'pMinFactor' are more reliable." << endl;
+ reduce(rhoRefMax, maxOp<scalar>());
+ reduce(rhoRefMin, minOp<scalar>());
+ }
- if (!rhoLimits)
+ if (dict.found("pMax"))
{
- FatalIOErrorInFunction(dict)
- << "'rhoMin' specified rather than 'pMinFactor'" << nl
- << " but the corresponding reference density cannot"
- " be evaluated from the boundary conditions." << nl
- << "Please specify 'pMinFactor' rather than 'rhoMin'"
- << exit(FatalError);
+ pMax_.value() = readScalar(dict.lookup("pMax"));
+ }
+ else if (dict.found("pMaxFactor"))
+ {
+ if (!pLimits)
+ {
+ FatalIOErrorInFunction(dict)
+ << "'pMaxFactor' specified rather than 'pMax'" << nl
+ << " but the corresponding reference pressure cannot"
+ " be evaluated from the boundary conditions." << nl
+ << " Please specify 'pMax' rather than 'pMaxFactor'"
+ << exit(FatalIOError);
+ }
+
+ const scalar pMaxFactor(readScalar(dict.lookup("pMaxFactor")));
+ pMax_ *= pMaxFactor;
+ }
+ else if (dict.found("rhoMax"))
+ {
+ // For backward-compatibility infer the pMax from rhoMax
+
+ IOWarningInFunction(dict)
+ << "'rhoMax' specified rather than 'pMax' or 'pMaxFactor'"
+ << nl
+ << " This is supported for backward-compatibility but "
+ "'pMax' or 'pMaxFactor' are more reliable." << endl;
+
+ if (!pLimits)
+ {
+ FatalIOErrorInFunction(dict)
+ << "'rhoMax' specified rather than 'pMax'" << nl
+ << " but the corresponding reference pressure cannot"
+ " be evaluated from the boundary conditions." << nl
+ << " Please specify 'pMax' rather than 'rhoMax'"
+ << exit(FatalIOError);
+ }
+
+ if (!rhoLimits)
+ {
+ FatalIOErrorInFunction(dict)
+ << "'rhoMax' specified rather than 'pMaxFactor'" << nl
+ << " but the corresponding reference density cannot"
+ " be evaluated from the boundary conditions." << nl
+ << " Please specify 'pMaxFactor' rather than 'rhoMax'"
+ << exit(FatalIOError);
+ }
+
+ dimensionedScalar rhoMax("rhoMax", dimDensity, dict);
+
+ pMax_ *= max(rhoMax.value()/rhoRefMax, 1);
}
- dimensionedScalar rhoMin("rhoMin", dimDensity, dict);
-
- pMin_ *= min(rhoMin.value()/rhoRefMin, 1);
+ if (dict.found("pMin"))
+ {
+ pMin_.value() = readScalar(dict.lookup("pMin"));
+ }
+ else if (dict.found("pMinFactor"))
+ {
+ if (!pLimits)
+ {
+ FatalIOErrorInFunction(dict)
+ << "'pMinFactor' specified rather than 'pMin'" << nl
+ << " but the corresponding reference pressure cannot"
+ " be evaluated from the boundary conditions." << nl
+ << " Please specify 'pMin' rather than 'pMinFactor'"
+ << exit(FatalIOError);
+ }
+
+ const scalar pMinFactor(readScalar(dict.lookup("pMinFactor")));
+ pMin_ *= pMinFactor;
+ }
+ else if (dict.found("rhoMin"))
+ {
+ // For backward-compatibility infer the pMin from rhoMin
+
+ IOWarningInFunction(dict)
+ << "'rhoMin' specified rather than 'pMin' or 'pMinFactor'" << nl
+ << " This is supported for backward-compatibility but"
+ "'pMin' or 'pMinFactor' are more reliable." << endl;
+
+ if (!pLimits)
+ {
+ FatalIOErrorInFunction(dict)
+ << "'rhoMin' specified rather than 'pMin'" << nl
+ << " but the corresponding reference pressure cannot"
+ " be evaluated from the boundary conditions." << nl
+ << " Please specify 'pMin' rather than 'rhoMin'"
+ << exit(FatalIOError);
+ }
+
+ if (!rhoLimits)
+ {
+ FatalIOErrorInFunction(dict)
+ << "'rhoMin' specified rather than 'pMinFactor'" << nl
+ << " but the corresponding reference density cannot"
+ " be evaluated from the boundary conditions." << nl
+ << " Please specify 'pMinFactor' rather than 'rhoMin'"
+ << exit(FatalIOError);
+ }
+
+ dimensionedScalar rhoMin("rhoMin", dimDensity, dict);
+
+ pMin_ *= min(rhoMin.value()/rhoRefMin, 1);
+ }
}
Info<< "pressureControl" << nl
diff --git a/src/finiteVolume/cfdTools/general/pressureControl/pressureControl.H b/src/finiteVolume/cfdTools/general/pressureControl/pressureControl.H
index 19bd053d430..61f510ac603 100644
--- a/src/finiteVolume/cfdTools/general/pressureControl/pressureControl.H
+++ b/src/finiteVolume/cfdTools/general/pressureControl/pressureControl.H
@@ -76,7 +76,8 @@ public:
(
const volScalarField& p,
const volScalarField& rho,
- const dictionary& dict
+ const dictionary& dict,
+ const bool pRefRequired = true
);
diff --git a/tutorials/compressible/rhoPimpleDyMFoam/annularThermalMixer/system/fvSolution b/tutorials/compressible/rhoPimpleDyMFoam/annularThermalMixer/system/fvSolution
index 66786f37722..c771e4cfe36 100644
--- a/tutorials/compressible/rhoPimpleDyMFoam/annularThermalMixer/system/fvSolution
+++ b/tutorials/compressible/rhoPimpleDyMFoam/annularThermalMixer/system/fvSolution
@@ -63,8 +63,9 @@ PIMPLE
nOuterCorrectors 3;
nCorrectors 1;
nNonOrthogonalCorrectors 0;
- rhoMin 0.5;
- rhoMax 2.0;
+
+ pMaxFactor 1.2;
+ pMinFactor 0.8;
}
relaxationFactors
diff --git a/tutorials/compressible/rhoPimpleFoam/LES/pitzDaily/system/fvSolution b/tutorials/compressible/rhoPimpleFoam/LES/pitzDaily/system/fvSolution
index 4eddb387211..b2dc2e505ef 100644
--- a/tutorials/compressible/rhoPimpleFoam/LES/pitzDaily/system/fvSolution
+++ b/tutorials/compressible/rhoPimpleFoam/LES/pitzDaily/system/fvSolution
@@ -52,8 +52,9 @@ PIMPLE
nOuterCorrectors 3;
nCorrectors 1;
nNonOrthogonalCorrectors 0;
- rhoMin 0.5;
- rhoMax 2.0;
+
+ pMinFactor 0.5;
+ pMaxFactor 2.0;
}
relaxationFactors
diff --git a/tutorials/compressible/rhoPimpleFoam/RAS/angledDuct/constant/thermophysicalProperties b/tutorials/compressible/rhoPimpleFoam/RAS/angledDuct/constant/thermophysicalProperties
index 32542b344e4..7fc31fc96b5 100644
--- a/tutorials/compressible/rhoPimpleFoam/RAS/angledDuct/constant/thermophysicalProperties
+++ b/tutorials/compressible/rhoPimpleFoam/RAS/angledDuct/constant/thermophysicalProperties
@@ -23,7 +23,7 @@ thermoType
thermo hConst;
equationOfState perfectGas;
specie specie;
- energy sensibleEnthalpy;
+ energy sensibleInternalEnergy;
}
mixture
diff --git a/tutorials/compressible/rhoPimpleFoam/RAS/angledDuct/system/fvSchemes b/tutorials/compressible/rhoPimpleFoam/RAS/angledDuct/system/fvSchemes
index 359d6ad3e11..c0522b83d6a 100644
--- a/tutorials/compressible/rhoPimpleFoam/RAS/angledDuct/system/fvSchemes
+++ b/tutorials/compressible/rhoPimpleFoam/RAS/angledDuct/system/fvSchemes
@@ -30,8 +30,9 @@ divSchemes
default none;
div(phi,U) Gauss upwind;
div(phid,p) Gauss upwind;
+ div(phiv,p) Gauss linear;
div(phi,K) Gauss linear;
- div(phi,h) Gauss upwind;
+ div(phi,e) Gauss upwind;
div(phi,k) Gauss upwind;
div(phi,epsilon) Gauss upwind;
div(phi,R) Gauss upwind;
diff --git a/tutorials/compressible/rhoPimpleFoam/RAS/angledDuct/system/fvSolution b/tutorials/compressible/rhoPimpleFoam/RAS/angledDuct/system/fvSolution
index ff244067813..3c7d9238f77 100644
--- a/tutorials/compressible/rhoPimpleFoam/RAS/angledDuct/system/fvSolution
+++ b/tutorials/compressible/rhoPimpleFoam/RAS/angledDuct/system/fvSolution
@@ -28,11 +28,10 @@ solvers
pFinal
{
$p;
- tolerance 1e-07;
relTol 0;
}
- "(rho|U|h|k|epsilon|omega)"
+ "(rho|U|e|k|epsilon|omega)"
{
solver smoothSolver;
smoother symGaussSeidel;
@@ -40,10 +39,9 @@ solvers
relTol 0.1;
}
- "(rho|U|h|k|epsilon|omega)Final"
+ "(rho|U|e|k|epsilon|omega)Final"
{
$U;
- tolerance 1e-06;
relTol 0;
}
}
@@ -57,8 +55,8 @@ PIMPLE
nNonOrthogonalCorrectors 0;
consistent yes;
- rhoMin 0.4;
- rhoMax 2.0;
+ pMaxFactor 1.5;
+ pMinFactor 0.9;
residualControl
{
@@ -76,13 +74,13 @@ relaxationFactors
{
fields
{
- "p.*" 0.9;
+ "p.*" 1;
"rho.*" 1;
}
equations
{
"U.*" 0.9;
- "h.*" 0.7;
+ "e.*" 0.7;
"(k|epsilon|omega).*" 0.8;
}
}
diff --git a/tutorials/compressible/rhoPimpleFoam/RAS/angledDuctLTS/system/fvSolution b/tutorials/compressible/rhoPimpleFoam/RAS/angledDuctLTS/system/fvSolution
index 83e599258fd..8b4068e0d6a 100644
--- a/tutorials/compressible/rhoPimpleFoam/RAS/angledDuctLTS/system/fvSolution
+++ b/tutorials/compressible/rhoPimpleFoam/RAS/angledDuctLTS/system/fvSolution
@@ -56,8 +56,8 @@ PIMPLE
nCorrectors 1;
nNonOrthogonalCorrectors 0;
- rhoMin 0.5;
- rhoMax 2.0;
+ pMaxFactor 1.5;
+ pMinFactor 0.9;
maxCo 0.2;
rDeltaTSmoothingCoeff 0.1;
diff --git a/tutorials/compressible/rhoPimpleFoam/RAS/cavity/system/fvSolution b/tutorials/compressible/rhoPimpleFoam/RAS/cavity/system/fvSolution
index e7c64b955d2..17ed812b9dd 100644
--- a/tutorials/compressible/rhoPimpleFoam/RAS/cavity/system/fvSolution
+++ b/tutorials/compressible/rhoPimpleFoam/RAS/cavity/system/fvSolution
@@ -59,8 +59,9 @@ PIMPLE
nOuterCorrectors 1;
nCorrectors 2;
nNonOrthogonalCorrectors 0;
- rhoMin 0.5;
- rhoMax 2.0;
+
+ pMax 1.2e5;
+ pMin 0.8e5;
}
diff --git a/tutorials/compressible/rhoPimpleFoam/RAS/mixerVessel2D/system/fvSolution b/tutorials/compressible/rhoPimpleFoam/RAS/mixerVessel2D/system/fvSolution
index d4a55d13c47..02ea2958833 100644
--- a/tutorials/compressible/rhoPimpleFoam/RAS/mixerVessel2D/system/fvSolution
+++ b/tutorials/compressible/rhoPimpleFoam/RAS/mixerVessel2D/system/fvSolution
@@ -84,14 +84,13 @@ solvers
PIMPLE
{
- nOuterCorrectors 1;
- nCorrectors 2;
- nNonOrthogonalCorrectors 0;
- momentumPredictor yes;
- rhoMin 0.5;
- rhoMax 2.0;
- pRefCell 0;
- pRefValue 1e5;
+ nOuterCorrectors 1;
+ nCorrectors 2;
+ nNonOrthogonalCorrectors 0;
+ momentumPredictor yes;
+
+ pMax 1.2e5;
+ pMin 0.8e5;
}
relaxationFactors
diff --git a/tutorials/compressible/rhoPimpleFoam/laminar/helmholtzResonance/system/fvSolution b/tutorials/compressible/rhoPimpleFoam/laminar/helmholtzResonance/system/fvSolution
index 4eddb387211..b2dc2e505ef 100644
--- a/tutorials/compressible/rhoPimpleFoam/laminar/helmholtzResonance/system/fvSolution
+++ b/tutorials/compressible/rhoPimpleFoam/laminar/helmholtzResonance/system/fvSolution
@@ -52,8 +52,9 @@ PIMPLE
nOuterCorrectors 3;
nCorrectors 1;
nNonOrthogonalCorrectors 0;
- rhoMin 0.5;
- rhoMax 2.0;
+
+ pMinFactor 0.5;
+ pMaxFactor 2.0;
}
relaxationFactors
--
GitLab