rhoSimpleFoam: added support for compressible liquid flows

rhoSimpleFoam now instantiates the lower-level fluidThermo which instantiates
either a psiThermo or rhoThermo according to the 'type' specification in
thermophysicalProperties, e.g.

thermoType
{
    type            hePsiThermo;
    mixture         pureMixture;
    transport       sutherland;
    thermo          janaf;
    equationOfState perfectGas;
    specie          specie;
    energy          sensibleInternalEnergy;
}

instantiates a psiThermo for a perfect gas with JANAF thermodynamics, whereas

thermoType
{
    type            heRhoThermo;
    mixture         pureMixture;
    properties      liquid;
    energy          sensibleInternalEnergy;
}

mixture
{
    H2O;
}

instantiates a rhoThermo for water, see new tutorial
compressible/rhoSimpleFoam/squareBendLiq.

In order to support complex equations of state the pressure can no longer be
unlimited and rhoSimpleFoam now limits the pressure rather than the density to
handle start-up more robustly.

For backward compatibility 'rhoMin' and 'rhoMax' can still be used in the SIMPLE
sub-dictionary of fvSolution which are converted into 'pMax' and 'pMin' but it
is better to set either 'pMax' and 'pMin' directly or use the more convenient
'pMinFactor' and 'pMinFactor' from which 'pMax' and 'pMin' are calculated using
the fixed boundary pressure or reference pressure e.g.

SIMPLE
{
    nNonOrthogonalCorrectors 0;

    pMinFactor      0.1;
    pMaxFactor      1.5;

    transonic       yes;
    consistent      yes;

    residualControl
    {
        p               1e-3;
        U               1e-4;
        e               1e-3;
        "(k|epsilon|omega)" 1e-3;
    }
}

applications/solvers/compressible/rhoSimpleFoam/createFields.H

@@ -1,10 +1,10 @@
 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 rho
@@ -38,35 +38,10 @@ volVectorField U
 
 #include "compressibleCreatePhi.H"
 
-
+pressureControl pressureControl(p, rho, simple.dict());
-label pRefCell = 0;
-scalar pRefValue = 0.0;
-setRefCell(p, simple.dict(), pRefCell, pRefValue);
 
 mesh.setFluxRequired(p.name());
 
-dimensionedScalar rhoMax
-(
-    dimensionedScalar::lookupOrDefault
-    (
-        "rhoMax",
-        simple.dict(),
-        dimDensity,
-        GREAT
-    )
-);
-
-dimensionedScalar rhoMin
-(
-    dimensionedScalar::lookupOrDefault
-    (
-        "rhoMin",
-        simple.dict(),
-        dimDensity,
-        0
-    )
-);
-
 Info<< "Creating turbulence model\n" << endl;
 autoPtr<compressible::turbulenceModel> turbulence
 (

applications/solvers/compressible/rhoSimpleFoam/pEqn.H

@@ -30,7 +30,11 @@
             // Relax the pressure equation to ensure diagonal-dominance
             pEqn.relax();
 
-            pEqn.setReference(pRefCell, pRefValue);
+            pEqn.setReference
+            (
+                pressureControl.refCell(),
+                pressureControl.refValue()
+            );
 
             pEqn.solve();
 
@@ -60,7 +64,11 @@
                 fvOptions(psi, p, rho.name())
             );
 
-            pEqn.setReference(pRefCell, pRefValue);
+            pEqn.setReference
+            (
+                pressureControl.refCell(),
+                pressureControl.refValue()
+            );
 
             pEqn.solve();
 
@@ -81,6 +89,8 @@
     U.correctBoundaryConditions();
     fvOptions.correct(U);
 
+    pressureControl.limit(p);
+
     // For closed-volume cases adjust the pressure and density levels
     // to obey overall mass continuity
     if (closedVolume)
@@ -89,16 +99,12 @@
             /fvc::domainIntegrate(psi);
     }
 
+    p.correctBoundaryConditions();
+
     rho = thermo.rho();
-    rho = max(rho, rhoMin);
-    rho = min(rho, rhoMax);
 
     if (!simple.transonic())
     {
         rho.relax();
     }
-
-    Info<< "rho max/min : "
-        << max(rho).value() << " "
-        << min(rho).value() << endl;
 }

applications/solvers/compressible/rhoSimpleFoam/pcEqn.H

@@ -40,7 +40,11 @@ if (simple.transonic())
         // Relax the pressure equation to maintain diagonal dominance
         pEqn.relax();
 
-        pEqn.setReference(pRefCell, pRefValue);
+        pEqn.setReference
+        (
+            pressureControl.refCell(),
+            pressureControl.refValue()
+        );
 
         pEqn.solve();
 
@@ -75,7 +79,11 @@ else
             fvOptions(psi, p, rho.name())
         );
 
-        pEqn.setReference(pRefCell, pRefValue);
+        pEqn.setReference
+        (
+            pressureControl.refCell(),
+            pressureControl.refValue()
+        );
 
         pEqn.solve();
 
@@ -97,6 +105,8 @@ U = HbyA - rAtU*fvc::grad(p);
 U.correctBoundaryConditions();
 fvOptions.correct(U);
 
+pressureControl.limit(p);
+
 // For closed-volume cases adjust the pressure and density levels
 // to obey overall mass continuity
 if (closedVolume)
@@ -105,14 +115,12 @@ if (closedVolume)
         /fvc::domainIntegrate(psi);
 }
 
+p.correctBoundaryConditions();
+
 // Recalculate density from the relaxed pressure
 rho = thermo.rho();
-rho = max(rho, rhoMin);
-rho = min(rho, rhoMax);
 
 if (!simple.transonic())
 {
     rho.relax();
 }
-
-Info<< "rho max/min : " << max(rho).value() << " " << min(rho).value() << endl;

applications/solvers/compressible/rhoSimpleFoam/rhoPorousSimpleFoam/createFields.H

@@ -1,10 +1,10 @@
 Info<< "Reading thermophysical properties\n" << endl;
 
-autoPtr<rhoThermo> pThermo
+autoPtr<fluidThermo> pThermo
 (
-    rhoThermo::New(mesh)
+    fluidThermo::New(mesh)
 );
-rhoThermo& thermo = pThermo();
+fluidThermo& thermo = pThermo();
 thermo.validate(args.executable(), "h", "e");
 
 volScalarField rho
@@ -38,46 +38,10 @@ volVectorField U
 
 #include "compressibleCreatePhi.H"
 
-
+pressureControl pressureControl(p, rho, simple.dict());
-label pRefCell = 0;
-scalar pRefValue = 0.0;
-setRefCell(p, simple.dict(), pRefCell, pRefValue);
 
 mesh.setFluxRequired(p.name());
 
-dimensionedScalar rhoMax
-(
-    dimensionedScalar::lookupOrDefault
-    (
-        "rhoMax",
-        simple.dict(),
-        dimDensity,
-        GREAT
-    )
-);
-
-dimensionedScalar rhoMin
-(
-    dimensionedScalar::lookupOrDefault
-    (
-        "rhoMin",
-        simple.dict(),
-        dimDensity,
-        0
-    )
-);
-
-dimensionedScalar pMin
-(
-    dimensionedScalar::lookupOrDefault
-    (
-        "pMin",
-        simple.dict(),
-        dimPressure,
-        10
-    )
-);
-
 Info<< "Creating turbulence model\n" << endl;
 autoPtr<compressible::turbulenceModel> turbulence
 (

applications/solvers/compressible/rhoSimpleFoam/rhoPorousSimpleFoam/pEqn.H

@@ -48,7 +48,11 @@
 
         fvScalarMatrix& pEqn = tpEqn.ref();
 
-        pEqn.setReference(pRefCell, pRefValue);
+        pEqn.setReference
+        (
+            pressureControl.refCell(),
+            pressureControl.refValue()
+        );
 
         pEqn.solve();
 
@@ -75,6 +79,8 @@
     U.correctBoundaryConditions();
     fvOptions.correct(U);
 
+    pressureControl.limit(p);
+
     // For closed-volume cases adjust the pressure and density levels
     // to obey overall mass continuity
     if (closedVolume)
@@ -83,13 +89,6 @@
             /fvc::domainIntegrate(psi);
     }
 
-    p = max(p, pMin);
-
     rho = thermo.rho();
-    rho = max(rho, rhoMin);
-    rho = min(rho, rhoMax);
     rho.relax();
-    Info<< "rho max/min : "
-        << max(rho).value() << " "
-        << min(rho).value() << endl;
 }

applications/solvers/compressible/rhoSimpleFoam/rhoPorousSimpleFoam/rhoPorousSimpleFoam.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
@@ -31,11 +31,12 @@ Description
 \*---------------------------------------------------------------------------*/
 
 #include "fvCFD.H"
-#include "rhoThermo.H"
+#include "fluidThermo.H"
 #include "turbulentFluidThermoModel.H"
+#include "simpleControl.H"
+#include "pressureControl.H"
 #include "fvOptions.H"
 #include "IOporosityModelList.H"
-#include "simpleControl.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 

applications/solvers/compressible/rhoSimpleFoam/rhoSimpleFoam.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
@@ -30,9 +30,10 @@ Description
 \*---------------------------------------------------------------------------*/
 
 #include "fvCFD.H"
-#include "psiThermo.H"
+#include "fluidThermo.H"
 #include "turbulentFluidThermoModel.H"
 #include "simpleControl.H"
+#include "pressureControl.H"
 #include "fvOptions.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

etc/templates/compressibleInflowOutflow/system/fvSchemes

@@ -44,9 +44,9 @@ divSchemes
     div(phi,e)      $turbulence;
     div(phi,h)      $turbulence;
     div(phi,K)      $turbulence;
-    div(phi,Ekp)    $turbulence; 
+    div(phi,Ekp)    $turbulence;
 
-    div(phid,p)     bounded Gauss upwind;
+    div(phid,p)     Gauss upwind;
     div((phi|interpolate(rho)),p)  bounded Gauss upwind;
 
     div(((rho*nuEff)*dev2(T(grad(U)))))    Gauss linear;

etc/templates/compressibleInflowOutflow/system/fvSolution

@@ -43,8 +43,8 @@ SIMPLE
     }
 
     nNonOrthogonalCorrectors 0;
-    rhoMin          rhoMin [ 1 -3 0 0 0 ] 0.1;
+    pMinFactor      0.1;
-    rhoMax          rhoMax [ 1 -3 0 0 0 ] 1.5;
+    pMaxFactor      1.5;
 }
 
 relaxationFactors
@@ -56,7 +56,7 @@ relaxationFactors
     }
     equations
     {
-     	U               0.7;
+        U               0.7;
         "(e|h)"         0.7;
         "(k|epsilon|omega)" 0.7;
     }

src/TurbulenceModels/compressible/CompressibleTurbulenceModel/CompressibleTurbulenceModel.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2013-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2013-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -115,6 +115,19 @@ public:
             return this->transport_.mu(patchi);
         }
 
+        //- Return the laminar viscosity
+        virtual tmp<volScalarField> nu() const
+        {
+            return this->transport_.mu()/this->rho_;
+        }
+
+        //- Return the laminar viscosity on patchi
+        virtual tmp<scalarField> nu(const label patchi) const
+        {
+            return
+                this->transport_.mu(patchi)/this->rho_.boundaryField()[patchi];
+        }
+
         //- Return the turbulence dynamic viscosity
         virtual tmp<volScalarField> mut() const
         {

src/TurbulenceModels/turbulenceModels/TurbulenceModel/TurbulenceModel.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2013-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2013-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -158,13 +158,13 @@ public:
         }
 
         //- Return the laminar viscosity
-        tmp<volScalarField> nu() const
+        virtual tmp<volScalarField> nu() const
         {
             return transport_.nu();
         }
 
         //- Return the laminar viscosity on patchi
-        tmp<scalarField> nu(const label patchi) const
+        virtual tmp<scalarField> nu(const label patchi) const
         {
             return transport_.nu(patchi);
         }

src/finiteVolume/Make/files

@@ -402,6 +402,7 @@ $(general)/constrainHbyA/constrainHbyA.C
 $(general)/adjustPhi/adjustPhi.C
 $(general)/bound/bound.C
 $(general)/CorrectPhi/correctUphiBCs.C
+$(general)/pressureControl/pressureControl.C
 
 solutionControl = $(general)/solutionControl
 $(solutionControl)/solutionControl/solutionControl.C

src/finiteVolume/cfdTools/general/findRefCell/findRefCell.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -27,7 +27,7 @@ License
 
 // * * * * * * * * * * * * * * * Global Functions  * * * * * * * * * * * * * //
 
-void Foam::setRefCell
+bool Foam::setRefCell
 (
     const volScalarField& field,
     const volScalarField& fieldRef,
@@ -109,11 +109,17 @@ void Foam::setRefCell
         }
 
         refValue = readScalar(dict.lookup(refValueName));
+
+        return true;
+    }
+    else
+    {
+        return false;
     }
 }
 
 
-void Foam::setRefCell
+bool Foam::setRefCell
 (
     const volScalarField& field,
     const dictionary& dict,
@@ -122,7 +128,7 @@ void Foam::setRefCell
     const bool forceReference
 )
 {
-    setRefCell(field, field, dict, refCelli, refValue, forceReference);
+    return setRefCell(field, field, dict, refCelli, refValue, forceReference);
 }
 
 

src/finiteVolume/cfdTools/general/findRefCell/findRefCell.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2017 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -46,8 +46,8 @@ namespace Foam
 
 //- If the field fieldRef needs referencing find the reference cell nearest
 //  (in index) to the given cell looked-up for field, but which is not on a
-//  cyclic, symmetry or processor patch.
+//  cyclic, symmetry or processor patch and return true, otherwise return false.
-void setRefCell
+bool setRefCell
 (
     const volScalarField& field,
     const volScalarField& fieldRef,
@@ -59,8 +59,8 @@ void setRefCell
 
 //- If the field needs referencing find the reference cell nearest
 //  (in index) to the given cell looked-up for field, but which is not on a
-//  cyclic, symmetry or processor patch.
+//  cyclic, symmetry or processor patch and return true, otherwise return false.
-void setRefCell
+bool setRefCell
 (
     const volScalarField& field,
     const dictionary& dict,

src/finiteVolume/cfdTools/general/pressureControl/pressureControl.C

@@ -0,0 +1,155 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2017 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "pressureControl.H"
+#include "findRefCell.H"
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::pressureControl::pressureControl
+(
+    const volScalarField& p,
+    const volScalarField& rho,
+    const dictionary& dict
+)
+:
+    refCell_(0),
+    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();
+
+    scalar rhoRefMax = -GREAT;
+    scalar rhoRefMin = GREAT;
+    bool rhoLimits = false;
+
+    forAll(pbf, patchi)
+    {
+        if (pbf[patchi].fixesValue())
+        {
+            rhoLimits = 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]));
+        }
+    }
+
+    if (dict.found("pMax"))
+    {
+        pMax_.value() = readScalar(dict.lookup("pMax"));
+    }
+    else if (dict.found("pMaxFactor"))
+    {
+        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 (!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(FatalError);
+        }
+
+        dimensionedScalar rhoMax("rhoMax", dimDensity, dict);
+
+        pMax_ *= max(rhoMax.value()/rhoRefMax, 1);
+    }
+
+    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
+
+        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 (!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(FatalError);
+        }
+
+        dimensionedScalar rhoMin("rhoMin", dimDensity, dict);
+
+        pMin_ *= min(rhoMin.value()/rhoRefMin, 1);
+    }
+
+    Info<< "pressureControl" << nl
+        << "    pMax/pMin " << pMax_.value() << " " << pMin_.value()
+        << nl << endl;
+}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+void Foam::pressureControl::limit(volScalarField& p) const
+{
+    Info<< "pressureControl: p max/min "
+        << max(p).value() << " "
+        << min(p).value() << endl;
+
+    p = max(p, pMin_);
+    p = min(p, pMax_);
+}
+
+
+// ************************************************************************* //

src/finiteVolume/cfdTools/general/pressureControl/pressureControl.H

@@ -0,0 +1,108 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2017 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::pressureControl
+
+Description
+    Provides controls for the pressure reference is closed-volume simulations
+    and a general method for limiting the pressure during the startup of
+    steady-state simulations.
+
+SourceFiles
+    pressureControlI.H
+    pressureControl.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef pressureControl_H
+#define pressureControl_H
+
+#include "dimensionedScalar.H"
+#include "volFieldsFwd.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                         Class pressureControl Declaration
+\*---------------------------------------------------------------------------*/
+
+class pressureControl
+{
+    // Private data
+
+        //- Optional cell in which the reference pressure is set
+        label refCell_;
+
+        //- Optional pressure reference level
+        scalar refValue_;
+
+        //- Pressure lower-limit
+        dimensionedScalar pMax_;
+
+        //- Pressure upper-limit
+        dimensionedScalar pMin_;
+
+
+public:
+
+    // Constructors
+
+        //- Construct from the SIMPLE/PIMPLE sub-dictionary
+        pressureControl
+        (
+            const volScalarField& p,
+            const volScalarField& rho,
+            const dictionary& dict
+        );
+
+
+    // Member Functions
+
+        //- Return the cell in which the reference pressure is set
+        inline label refCell() const;
+
+        //- Return the pressure reference level
+        inline scalar refValue() const;
+
+        //- Limit the pressure
+        void limit(volScalarField& p) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#include "pressureControlI.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/finiteVolume/cfdTools/general/pressureControl/pressureControlI.H

@@ -0,0 +1,40 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2017 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+inline Foam::label Foam::pressureControl::refCell() const
+{
+    return refCell_;
+}
+
+
+inline Foam::scalar Foam::pressureControl::refValue() const
+{
+    return refValue_;
+}
+
+
+// ************************************************************************* //

tutorials/compressible/rhoPorousSimpleFoam/angledDuctExplicit/system/fvSolution

@@ -54,8 +54,8 @@ solvers
 SIMPLE
 {
     nNonOrthogonalCorrectors 0;
-    rhoMin          0.5;
+    pMinFactor      0.5;
-    rhoMax          1.5;
+    pMaxFactor      1.5;
 
     residualControl
     {

tutorials/compressible/rhoPorousSimpleFoam/angledDuctImplicit/system/fvSolution

@@ -48,8 +48,8 @@ SIMPLE
 {
     nUCorrectors        2;
     nNonOrthogonalCorrectors 0;
-    rhoMin              0.5;
+    pMinFactor          0.5;
-    rhoMax              2.0;
+    pMaxFactor          2.0;
 
     residualControl
     {

tutorials/compressible/rhoSimpleFoam/angledDuctExplicitFixedCoeff/system/fvSolution

@@ -54,8 +54,8 @@ solvers
 SIMPLE
 {
     nNonOrthogonalCorrectors 0;
-    rhoMin          0.4;
+    pMinFactor      0.4;
-    rhoMax          1.5;
+    pMaxFactor      1.5;
 
     residualControl
     {

tutorials/compressible/rhoSimpleFoam/squareBend/system/fvSolution

@@ -39,8 +39,8 @@ solvers
 SIMPLE
 {
     nNonOrthogonalCorrectors 0;
-    rhoMin          0.1;
+    pMinFactor      0.1;
-    rhoMax          1.0;
+    pMaxFactor      1.5;
     transonic       yes;
     consistent      yes;
 
@@ -60,7 +60,6 @@ relaxationFactors
     fields
     {
         p               1;
-        rho             1;
     }
     equations
     {

tutorials/compressible/rhoSimpleFoam/squareBendLiq/0/T

@@ -0,0 +1,42 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    object      T;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 0 0 1 0 0 0];
+
+internalField   uniform 300;
+
+boundaryField
+{
+    Default_Boundary_Region
+    {
+        type            fixedValue;
+        value           uniform 350;
+    }
+
+    inlet
+    {
+        type            fixedValue;
+        value           $internalField;
+    }
+
+    outlet
+    {
+        type            inletOutlet;
+        inletValue      $internalField;
+    }
+}
+
+// ************************************************************************* //

tutorials/compressible/rhoSimpleFoam/squareBendLiq/0/U

@@ -0,0 +1,42 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volVectorField;
+    object      U;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 1 -1 0 0 0 0];
+
+internalField   uniform (0 0 0);
+
+boundaryField
+{
+    Default_Boundary_Region
+    {
+        type            noSlip;
+    }
+    inlet
+    {
+        type            flowRateInletVelocity;
+        massFlowRate    constant 5;
+        rhoInlet        1000;    // Guess for rho
+    }
+    outlet
+    {
+        type            inletOutlet;
+        value           uniform (0 0 0);
+        inletValue      uniform (0 0 0);
+    }
+}
+
+
+// ************************************************************************* //

tutorials/compressible/rhoSimpleFoam/squareBendLiq/0/alphat

@@ -0,0 +1,43 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      alphat;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [1 -1 -1 0 0 0 0];
+
+internalField   uniform 0;
+
+boundaryField
+{
+    Default_Boundary_Region
+    {
+        type            compressible::alphatWallFunction;
+        Prt             0.85;
+        value           uniform 0;
+    }
+    inlet
+    {
+        type            calculated;
+        value           uniform 0;
+    }
+    outlet
+    {
+        type            calculated;
+        value           uniform 0;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/compressible/rhoSimpleFoam/squareBendLiq/0/epsilon

@@ -0,0 +1,47 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      epsilon;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 2 -3 0 0 0 0];
+
+internalField   uniform 200;
+
+boundaryField
+{
+    Default_Boundary_Region
+    {
+        type            epsilonWallFunction;
+        Cmu             0.09;
+        kappa           0.41;
+        E               9.8;
+        value           uniform 200;
+    }
+    inlet
+    {
+        type            turbulentMixingLengthDissipationRateInlet;
+        mixingLength    0.005;
+        value           uniform 200;
+    }
+    outlet
+    {
+        type            inletOutlet;
+        inletValue      uniform 200;
+        value           uniform 200;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/compressible/rhoSimpleFoam/squareBendLiq/0/k

@@ -0,0 +1,44 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      k;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 2 -2 0 0 0 0];
+
+internalField   uniform 1;
+
+boundaryField
+{
+    Default_Boundary_Region
+    {
+        type            kqRWallFunction;
+        value           uniform 1;
+    }
+    inlet
+    {
+        type            turbulentIntensityKineticEnergyInlet;
+        intensity       0.05;
+        value           uniform 1;
+    }
+    outlet
+    {
+        type            inletOutlet;
+        inletValue      uniform 1;
+        value           uniform 1;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/compressible/rhoSimpleFoam/squareBendLiq/0/nut

@@ -0,0 +1,45 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      nut;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 2 -1 0 0 0 0];
+
+internalField   uniform 0;
+
+boundaryField
+{
+    Default_Boundary_Region
+    {
+        type            nutkWallFunction;
+        Cmu             0.09;
+        kappa           0.41;
+        E               9.8;
+        value           uniform 0;
+    }
+    inlet
+    {
+        type            calculated;
+        value           uniform 0;
+    }
+    outlet
+    {
+        type            calculated;
+        value           uniform 0;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/compressible/rhoSimpleFoam/squareBendLiq/0/p

@@ -0,0 +1,39 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    object      p;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [1 -1 -2 0 0 0 0];
+
+internalField   uniform 1e5;
+
+boundaryField
+{
+    Default_Boundary_Region
+    {
+        type            zeroGradient;
+    }
+    inlet
+    {
+        type            zeroGradient;
+    }
+    outlet
+    {
+        type            fixedValue;
+        value           $internalField;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/compressible/rhoSimpleFoam/squareBendLiq/constant/thermophysicalProperties

@@ -0,0 +1,31 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "constant";
+    object      thermophysicalProperties;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+thermoType
+{
+    type            heRhoThermo;
+    mixture         pureMixture;
+    properties      liquid;
+    energy          sensibleInternalEnergy;
+}
+
+mixture
+{
+    H2O;
+}
+
+// ************************************************************************* //

tutorials/compressible/rhoSimpleFoam/squareBendLiq/constant/turbulenceProperties

@@ -0,0 +1,30 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "constant";
+    object      turbulenceProperties;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+simulationType  RAS;
+
+RAS
+{
+    RASModel        kEpsilon;
+
+    turbulence      on;
+
+    printCoeffs     on;
+}
+
+
+// ************************************************************************* //

tutorials/compressible/rhoSimpleFoam/squareBendLiq/system/blockMeshDict

@@ -0,0 +1,127 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    object      blockMeshDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+convertToMeters 0.001;
+
+vertices
+(
+    // front-plane: z = +25mm
+    // inlet region
+    (   -50  25   25)           // pt 0
+    (     0  25   25)           // pt 1
+    (   -50  75   25)           // pt 2
+    (     0  75   25)           // pt 3
+    // outlet region
+    (  -500 -75   25)           // pt 4
+    (     0 -75   25)           // pt 5
+    (  -500 -25   25)           // pt 6
+    (     0 -25   25)           // pt 7
+    // bend mid-points
+    (    25   0   25)           // pt 8
+    (    75   0   25)           // pt 9
+    // back-plane: z = -25mm
+    // inlet region
+    (   -50  25   -25)          // pt 0 + 10
+    (     0  25   -25)          // pt 1 + 10
+    (   -50  75   -25)          // pt 2 + 10
+    (     0  75   -25)          // pt 3 + 10
+    // outlet region
+    (  -500 -75   -25)          // pt 4 + 10
+    (     0 -75   -25)          // pt 5 + 10
+    (  -500 -25   -25)          // pt 7 + 10
+    (     0 -25   -25)          // pt 8 + 10
+    // bend mid-points
+    (    25   0   -25)          // pt 8 + 10
+    (    75   0   -25)          // pt 9 + 10
+);
+
+blocks
+(
+    hex (0 1 11 10  2 3 13 12) inlet  ( 20 20 20)  simpleGrading (1 1 1)
+    hex (4 5 15 14  6 7 17 16) outlet (200 20 20)  simpleGrading (1 1 1)
+
+    hex (1 8 18 11  3 9 19 13) bend1  ( 30 20 20)  simpleGrading (1 1 1)
+    hex (5 9 19 15  7 8 18 17) bend2  ( 30 20 20)  simpleGrading (1 1 1)
+);
+
+edges
+(
+   // block 2
+   arc  1  8  ( 17.678  17.678  25)
+   arc 11 18  ( 17.678  17.678 -25)
+   arc  3  9  ( 53.033  53.033  25)
+   arc 13 19  ( 53.033  53.033 -25)
+   // block 3
+   arc  7  8  ( 17.678  -17.678  25)
+   arc 17 18  ( 17.678  -17.678 -25)
+   arc  5  9  ( 53.033  -53.033  25)
+   arc 15 19  ( 53.033  -53.033 -25)
+);
+
+boundary
+(
+        // is there no way of defining all my 'defaultFaces' to be 'wall'?
+    Default_Boundary_Region
+    {
+        type wall;
+        faces
+        (
+            // block0
+            ( 0 1 3 2 )
+            ( 11 10 12 13 )
+            ( 0 10 11 1 )
+            ( 2 3 13 12 )
+            // block1
+            ( 4 5 7 6 )
+            ( 15 14 16 17 )
+            ( 4 14 15 5 )
+            ( 6 7 17 16 )
+            // block2
+            ( 1 8 9 3 )
+            ( 18 11 13 19 )
+            ( 3 9 19 13 )
+            ( 1 11 18 8 )
+            // block3
+            ( 5 9 8 7 )
+            ( 19 15 17 18 )
+            ( 5 15 19 9 )
+            ( 7 8 18 17 )
+        );
+    }
+    inlet
+    {
+        type patch;
+        faces
+        (
+            (0 2 12 10)
+        );
+    }
+
+    outlet
+    {
+        type patch;
+        faces
+        (
+            (4 6 16 14)
+        );
+    }
+);
+
+mergePatchPairs
+(
+);
+
+// ************************************************************************* //

tutorials/compressible/rhoSimpleFoam/squareBendLiq/system/controlDict

@@ -0,0 +1,51 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "system";
+    object      controlDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+application     rhoSimpleFoam;
+
+startFrom       startTime;
+
+startTime       0;
+
+stopAt          endTime;
+
+endTime         500;
+
+deltaT          1;
+
+writeControl    timeStep;
+
+writeInterval   100;
+
+purgeWrite      0;
+
+writeFormat     ascii;
+
+writePrecision  6;
+
+writeCompression off;
+
+timeFormat      general;
+
+timePrecision   6;
+
+graphFormat     raw;
+
+runTimeModifiable true;
+
+
+// ************************************************************************* //

tutorials/compressible/rhoSimpleFoam/squareBendLiq/system/decomposeParDict

@@ -0,0 +1,45 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "system";
+    object      decomposeParDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+numberOfSubdomains 8;
+
+method          hierarchical;
+
+simpleCoeffs
+{
+    n               (8 1 1);
+    delta           0.001;
+}
+
+hierarchicalCoeffs
+{
+    n               (4 2 1);
+    delta           0.001;
+    order           xyz;
+}
+
+manualCoeffs
+{
+    dataFile        "";
+}
+
+distributed     no;
+
+roots           ( );
+
+
+// ************************************************************************* //

tutorials/compressible/rhoSimpleFoam/squareBendLiq/system/fvSchemes

@@ -0,0 +1,59 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "system";
+    object      fvSchemes;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+ddtSchemes
+{
+    default             steadyState;
+}
+
+gradSchemes
+{
+    default             Gauss linear;
+
+    limited             cellLimited Gauss linear 1;
+}
+
+divSchemes
+{
+    default             none;
+
+    div(phi,U)          bounded Gauss linearUpwind limited;
+    div(phi,e)          bounded Gauss linearUpwind limited;
+    div(phi,epsilon)    bounded Gauss linearUpwind limited;
+    div(phi,k)          bounded Gauss linearUpwind limited;
+    div(phi,Ekp)        bounded Gauss linearUpwind limited;
+
+    div(((rho*nuEff)*dev2(T(grad(U)))))      Gauss linear;
+}
+
+laplacianSchemes
+{
+    default         Gauss linear corrected;
+}
+
+interpolationSchemes
+{
+    default         linear;
+}
+
+snGradSchemes
+{
+    default         corrected;
+}
+
+
+// ************************************************************************* //

tutorials/compressible/rhoSimpleFoam/squareBendLiq/system/fvSolution

@@ -0,0 +1,71 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "system";
+    object      fvSolution;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+solvers
+{
+    p
+    {
+        solver          GAMG;
+        smoother        GaussSeidel;
+
+        tolerance       1e-7;
+        relTol          0.1;
+    }
+
+    "(U|e|k|epsilon)"
+    {
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
+
+        tolerance       1e-7;
+        relTol          0.1;
+    }
+}
+
+SIMPLE
+{
+    nNonOrthogonalCorrectors 0;
+    pMinFactor      0.1;
+    pMaxFactor      1.5;
+
+    transonic       no;
+    consistent      no;
+
+    residualControl
+    {
+        p               1e-4;
+        U               1e-4;
+        "(k|omega|epsilon|e|h)" 1e-3;
+    }
+}
+
+relaxationFactors
+{
+    fields
+    {
+        p               0.3;
+    }
+    equations
+    {
+        U               0.7;
+        e               0.7;
+        k               0.7;
+        epsilon         0.7;
+    }
+}
+
+// ************************************************************************* //