ENH: Pressure function object - added optional inclusion of hydrostatic pressure

Example - create p-rgh from p:

    pressure-p-rgh
    {
        type            pressure;
        libs            (fieldFunctionObjects);
        writeControl    writeTime;
        mode            static;
        rho             rhoInf;
        rhoInf          1;
        hydrostaticMode subtract;
        g               (0 -9.81 0);
        hRef            0;
    }

- the hydrostaticMode entry is optional - setting is shown during construction
- g and/or hRef values are retrieved from the database if not specified

src/functionObjects/field/pressure/pressure.C

@@ -29,6 +29,7 @@ License
 #include "pressure.H"
 #include "volFields.H"
 #include "basicThermo.H"
+#include "uniformDimensionedFields.H"
 #include "addToRunTimeSelectionTable.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@@ -55,6 +56,17 @@ Foam::functionObjects::pressure::modeNames
     { TOTAL_COEFF, "totalCoeff" },
 });
 
+const Foam::Enum
+<
+    Foam::functionObjects::pressure::hydrostaticMode
+>
+Foam::functionObjects::pressure::hydrostaticModeNames
+({
+    { NONE, "none" },
+    { ADD, "add" },
+    { SUBTRACT, "subtract" },
+});
+
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
 Foam::word Foam::functionObjects::pressure::resultName() const
@@ -80,6 +92,26 @@ Foam::word Foam::functionObjects::pressure::resultName() const
             << abort(FatalError);
     }
 
+    switch (hydrostaticMode_)
+    {
+        case NONE:
+        {
+            break;
+        }
+        case ADD:
+        {
+            rName = rName + "+rgh";
+
+            break;
+        }
+        case SUBTRACT:
+        {
+            rName = rName + "-rgh";
+
+            break;
+        }
+    }
+
     if (mode_ & COEFF)
     {
         rName += "_coeff";
@@ -141,18 +173,90 @@ Foam::tmp<Foam::volScalarField> Foam::functionObjects::pressure::rhoScale
 }
 
 
+void Foam::functionObjects::pressure::addHydrostaticContribution
+(
+    volScalarField& p
+) const
+{
+    // Add/subtract hydrostatic contribution
+
+    if (hydrostaticMode_ == NONE)
+    {
+        return;
+    }
+
+    if (!gInitialised_)
+    {
+        g_ = mesh_.time().lookupObject<uniformDimensionedVectorField>("g");
+    }
+
+    if (!hRefInitialised_)
+    {
+        hRef_ = mesh_.lookupObject<uniformDimensionedScalarField>("hRef");
+    }
+
+    const dimensionedScalar ghRef
+    (
+        (g_ & (cmptMag(g_.value())/mag(g_.value())))*hRef_
+    );
+
+    tmp<volScalarField> rgh = rhoScale(p, (g_ & mesh_.C()) - ghRef);
+
+    switch (hydrostaticMode_)
+    {
+        case ADD:
+        {
+            p += rgh;
+            break;
+        }
+        case SUBTRACT:
+        {
+            p -= rgh;
+            break;
+        }
+        default:
+        {}
+    }
+}
+
+
 Foam::tmp<Foam::volScalarField> Foam::functionObjects::pressure::calcPressure
 (
     const volScalarField& p,
     const tmp<volScalarField>& tp
 ) const
 {
+    // Initialise to the pressure reference level
+    auto tresult =
+        tmp<volScalarField>::New
+        (
+            IOobject
+            (
+                name() + ":p",
+                mesh_.time().timeName(),
+                mesh_,
+                IOobject::NO_READ
+            ),
+            mesh_,
+            dimensionedScalar("p", dimPressure, pRef_)
+        );
+
+    volScalarField& result = tresult.ref();
+
+    addHydrostaticContribution(result);
+
+    if (mode_ & STATIC)
+    {
+        result += tp;
+        return tresult;
+    }
+
     if (mode_ & TOTAL)
     {
-        return
+        result +=
             tp
-          + dimensionedScalar("pRef", dimPressure, pRef_)
           + rhoScale(p, 0.5*magSqr(lookupObject<volVectorField>(UName_)));
+        return tresult;
     }
 
     if (mode_ & ISENTROPIC)
@@ -175,10 +279,12 @@ Foam::tmp<Foam::volScalarField> Foam::functionObjects::pressure::calcPressure
            /sqrt(gamma*tp.ref()/thermoPtr->rho())
         );
 
-        return tp()*(pow(1 + (gamma - 1)/2*sqr(Mb), gamma/(gamma - 1)));
+        result += tp*(pow(1 + (gamma - 1)/2*sqr(Mb), gamma/(gamma - 1)));
+        return tresult;
     }
 
-    return tp + dimensionedScalar("pRef", dimPressure, pRef_);
+
+    return tresult;
 }
 
 
@@ -246,13 +352,18 @@ Foam::functionObjects::pressure::pressure
 :
     fieldExpression(name, runTime, dict, "p"),
     mode_(STATIC),
+    hydrostaticMode_(NONE),
     UName_("U"),
     rhoName_("rho"),
     pRef_(0),
     pInf_(0),
     UInf_(Zero),
     rhoInf_(1),
-    rhoInfInitialised_(false)
+    rhoInfInitialised_(false),
+    g_(dimAcceleration),
+    gInitialised_(false),
+    hRef_(dimLength),
+    hRefInitialised_(false)
 {
     read(dict);
 }
@@ -300,10 +411,34 @@ bool Foam::functionObjects::pressure::read(const dictionary& dict)
         }
     }
 
-    Info<< "    operating mode: " << modeNames[mode_] << nl;
+    Info<< "    Operating mode: " << modeNames[mode_] << nl;
 
     pRef_ = dict.lookupOrDefault<scalar>("pRef", 0);
 
+    if
+    (
+        hydrostaticModeNames.readIfPresent
+        (
+            "hydrostaticMode",
+            dict,
+            hydrostaticMode_
+        )
+     && hydrostaticMode_
+    )
+    {
+        Info<< "    Hydrostatic mode: "
+            << hydrostaticModeNames[hydrostaticMode_]
+            << nl;
+        gInitialised_ = dict.readIfPresent("g", g_);
+        hRefInitialised_ = dict.readIfPresent("hRef", hRef_);
+    }
+    else
+    {
+        Info<< "    Not including hydrostatic effects" << nl;
+    }
+
+
+
     if (mode_ & COEFF)
     {
         dict.readEntry("pInf", pInf_);

src/functionObjects/field/pressure/pressure.H

@@ -100,6 +100,7 @@ Usage
         pInf         | Freestream pressure for coefficient calculation | no |
         UInf         | Freestream velocity for coefficient calculation | no |
         rhoInf       | Freestream density for coefficient calculation  | no |
+        hydrostaticMode | Hydrostatic contributions (see below) | no | none
     \endtable
 
     The \c mode entry is used to select the type of pressure that is calculated.
@@ -110,6 +111,13 @@ Usage
     - staticCoeff
     - totalCoeff
 
+    The optional \c hydrostaticMode entry provides handling for the term
+    \f$ \rho (\vec{g} \dot \vec{h})\f$ where options include
+    - \c none : not included
+    - \c add : add the term, e.g. to convert from p_rgh to p
+    - \c subtract : subtract the term, e.g. to convert from p to p_rgh
+
+
 See also
     Foam::functionObjects::fieldExpression
     Foam::functionObjects::fvMeshFunctionObject
@@ -124,6 +132,7 @@ SourceFiles
 
 #include "fieldExpression.H"
 #include "volFieldsFwd.H"
+#include "dimensionedVector.H"
 #include "dimensionedScalar.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -148,16 +157,26 @@ public:
         //- Enumeration for pressure calculation mode
         enum mode : unsigned
         {
-            STATIC = 0x1,           //!< Static pressure
+            STATIC       = (1 << 0), //!< Static pressure
-            TOTAL = 0x2,            //!< Total pressure
+            TOTAL        = (1 << 1), //!< Total pressure
-            ISENTROPIC = 0x4,       //!< Isentropic pressure
+            ISENTROPIC   = (1 << 2), //!< Isentropic pressure
-            COEFF = 0x8,            //!< Coefficient manipulator
+            COEFF        = (1 << 3), //!< Coefficient manipulator
             STATIC_COEFF = (STATIC | COEFF),
-            TOTAL_COEFF = (TOTAL | COEFF)
+            TOTAL_COEFF  = (TOTAL | COEFF),
         };
 
         static const Enum<mode> modeNames;
 
+        //- Enumeration for hydrostatic contributions
+        enum hydrostaticMode : unsigned
+        {
+            NONE = 0,
+            ADD,
+            SUBTRACT
+        };
+
+        static const Enum<hydrostaticMode> hydrostaticModeNames;
+
 
 private:
 
@@ -166,6 +185,9 @@ private:
         //- Calculation mode
         mode mode_;
 
+        //- Hydrostatic constribution mode
+        hydrostaticMode hydrostaticMode_;
+
         //- Name of velocity field, default is "U"
         word UName_;
 
@@ -194,6 +216,21 @@ private:
             bool rhoInfInitialised_;
 
 
+        //- p +/- rgh calculation
+
+            //- Gravity vector
+            mutable dimensionedVector g_;
+
+            //- Flag to show whether g has been initialised
+            bool gInitialised_;
+
+            //- Reference height
+            mutable dimensionedScalar hRef_;
+
+            //- Flag to show whether hRef has been initialised
+            bool hRefInitialised_;
+
+
     // Private Member Functions
 
         //- Return the name of the derived pressure field
@@ -209,8 +246,8 @@ private:
             const tmp<volScalarField>& tsf
         ) const;
 
-        //- Return the reference pressure
+        //- Add the hydrostatic contribution
-        tmp<volScalarField> pRef(const tmp<volScalarField>& tp) const;
+        void addHydrostaticContribution(volScalarField& p) const;
 
         //- Calculate and return the pressure
         tmp<volScalarField> calcPressure