diff --git a/src/functionObjects/field/Make/files b/src/functionObjects/field/Make/files
index ddcb030299..e11d10f36c 100644
--- a/src/functionObjects/field/Make/files
+++ b/src/functionObjects/field/Make/files
@@ -6,6 +6,8 @@ continuityError/continuityError.C
 
 derivedFields/derivedFields.C
 
+expressions/fvExpressionField.C
+
 fieldAverage/fieldAverage.C
 fieldAverage/fieldAverageItem/fieldAverageItem.C
 fieldAverage/fieldAverageItem/fieldAverageItemIO.C
diff --git a/src/functionObjects/field/expressions/fvExpressionField.C b/src/functionObjects/field/expressions/fvExpressionField.C
new file mode 100644
index 0000000000..fb395e7538
--- /dev/null
+++ b/src/functionObjects/field/expressions/fvExpressionField.C
@@ -0,0 +1,696 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2021 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+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 "fvExpressionField.H"
+#include "volFields.H"
+#include "surfaceFields.H"
+#include "pointMesh.H"
+#include "pointFields.H"
+#include "volumeExprDriver.H"
+#include "calculatedFvPatchField.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace functionObjects
+{
+    defineTypeNameAndDebug(fvExpressionField, 0);
+    addToRunTimeSelectionTable(functionObject, fvExpressionField, dictionary);
+}
+}
+
+
+const Foam::Enum
+<
+    Foam::functionObjects::fvExpressionField::actionType
+>
+Foam::functionObjects::fvExpressionField::actionNames_
+({
+    { actionType::opNone, "none" },
+    { actionType::opNew, "new" },
+    { actionType::opModify, "modify" },
+});
+
+
+// * * * * * * * * * * * * * * * Local Functions * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+word fieldGeoType(const expressions::FieldAssociation geoType)
+{
+    switch (geoType)
+    {
+        case expressions::FieldAssociation::POINT_DATA : return "points"; break;
+        case expressions::FieldAssociation::FACE_DATA : return "faces"; break;
+        case expressions::FieldAssociation::VOLUME_DATA : return "cells"; break;
+        default: break;
+    }
+    return "unknown";
+}
+
+
+template<class Type>
+static void doCorrectBoundaryConditions
+(
+    bool correctBCs,
+    GeometricField<Type, fvPatchField, volMesh>& field
+)
+{
+    if (correctBCs)
+    {
+        // Info<< "Correcting boundary conditions: " << field.name() << nl;
+        field.correctBoundaryConditions();
+
+        // Ensure that calculated patches are updated
+        for (auto& pf : field.boundaryFieldRef())
+        {
+            if (isA<calculatedFvPatchField<Type>>(pf))
+            {
+                pf = pf.patchInternalField();
+            }
+        }
+    }
+}
+
+
+template<class Type>
+void doCorrectBoundaryConditions
+(
+    bool correctBCs,
+    GeometricField<Type, pointPatchField, pointMesh>& field
+)
+{
+    if (correctBCs)
+    {
+        // Info<< "Correcting boundary conditions: " << field.name() << nl;
+        field.correctBoundaryConditions();
+    }
+}
+
+
+template<class Type>
+void doCorrectBoundaryConditions
+(
+    bool correctBCs,
+    GeometricField<Type, fvsPatchField, surfaceMesh>& field
+)
+{}
+
+} // End namespace Foam
+
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+template<class FieldType>
+bool Foam::functionObjects::fvExpressionField::loadAndStore(const IOobject& io)
+{
+    if (FieldType::typeName == io.headerClassName())
+    {
+        // Store field on mesh database
+        Log << "    Reading " << io.name()
+            << " (" << FieldType::typeName << ')' << endl;
+
+        mesh_.objectRegistry::store(new FieldType(io, mesh_));
+        return true;
+    }
+
+    return false;
+}
+
+
+template<class Type>
+bool Foam::functionObjects::fvExpressionField::loadField(const IOobject& io)
+{
+    typedef GeometricField<Type, fvPatchField, volMesh> VolFieldType;
+    // typedef typename VolFieldType::Internal IntVolFieldType;
+    typedef GeometricField<Type, fvsPatchField, surfaceMesh> SurfaceFieldType;
+
+    return
+    (
+        loadAndStore<VolFieldType>(io)
+        /// || loadAndStore<IntVolFieldType>(io)
+     || loadAndStore<SurfaceFieldType>(io)
+    );
+}
+
+
+Foam::label Foam::functionObjects::fvExpressionField::loadFields
+(
+    const UList<word>& fieldSet_
+)
+{
+    label nLoaded = 0;
+
+    for (const word& fieldName : fieldSet_)
+    {
+        // Already loaded?
+        const auto* ptr = mesh_.cfindObject<regIOobject>(fieldName);
+
+        if (ptr)
+        {
+            ++nLoaded;
+            DebugInfo
+                << "readFields : "
+                << ptr->name() << " (" << ptr->type()
+                << ") already in database" << endl;
+            continue;
+        }
+
+        // Load field as necessary
+        IOobject io
+        (
+            fieldName,
+            mesh_.time().timeName(),
+            mesh_,
+            IOobject::MUST_READ,
+            IOobject::NO_WRITE
+        );
+
+        const bool ok =
+        (
+            io.typeHeaderOk<regIOobject>(false) // Preload header info
+         && !io.headerClassName().empty()       // Extra safety
+         &&
+            (
+                loadField<scalar>(io)
+             || loadField<vector>(io)
+             || loadField<sphericalTensor>(io)
+             || loadField<symmTensor>(io)
+             || loadField<tensor>(io)
+            )
+        );
+
+        if (ok)
+        {
+            ++nLoaded;
+        }
+        else
+        {
+            DebugInfo
+                << "readFields : failed to load " << fieldName << endl;
+        }
+    }
+
+    return nLoaded;
+}
+
+
+template<class GeoField>
+bool Foam::functionObjects::fvExpressionField::setField
+(
+    GeoField& output,
+    const GeoField& evaluated,
+    const boolField& fieldMask
+)
+{
+    label numValuesChanged = 0;
+
+    // Internal field
+    if (fieldMask.empty())
+    {
+        // No field-mask - set all
+        numValuesChanged = output.size();
+
+        output.primitiveFieldRef() = evaluated;
+    }
+    else
+    {
+        auto& internal = output.primitiveFieldRef();
+
+        forAll(internal, idx)
+        {
+            if (fieldMask[idx])
+            {
+                internal[idx] = evaluated[idx];
+                ++numValuesChanged;
+            }
+        }
+    }
+
+    // Boundary fields
+    forAll(evaluated.boundaryField(), patchi)
+    {
+        auto& pf = output.boundaryFieldRef()[patchi];
+
+        if (pf.patch().coupled())
+        {
+            pf == evaluated.boundaryField()[patchi];
+        }
+    }
+
+    doCorrectBoundaryConditions(true, output);
+
+    if (action_ == actionType::opModify && log)
+    {
+        const label numTotal = returnReduce(output.size(), plusOp<label>());
+        reduce(numValuesChanged, plusOp<label>());
+
+        Info<< this->name() << ": set ";
+        if (numValuesChanged == numTotal)
+        {
+            Info<< "all ";
+        }
+        else
+        {
+            Info<< numValuesChanged << " of ";
+        }
+        Info<< numTotal << " values (field: "
+            << output.name() << ')' << nl << endl;
+    }
+
+    if (hasDimensions_)
+    {
+        // Log<< "Setting dimensions to " << dims << endl;
+        output.dimensions().reset(dimensions_);
+    }
+
+    return true;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::functionObjects::fvExpressionField::fvExpressionField
+(
+    const word& name,
+    const Time& runTime,
+    const dictionary& dict,
+    const bool loadFromFiles
+)
+:
+    fvMeshFunctionObject(name, runTime, dict),
+    dict_(dict),  // Deep copy
+    fieldName_(),
+    preloadFields_(),
+    maskExpr_(),
+    valueExpr_(),
+    dimensions_(),
+    action_(actionType::opNew),
+    autowrite_(false),
+    store_(true),
+    hasDimensions_(false),
+    loadFromFiles_(loadFromFiles)
+{
+    read(dict);
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::functionObjects::fvExpressionField::~fvExpressionField()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+Foam::word Foam::functionObjects::fvExpressionField::fieldName() const
+{
+    switch (action_)
+    {
+        case actionType::opNone:
+        {
+            break;  // No-op
+        }
+        case actionType::opNew:
+        {
+            return scopedName(fieldName_);
+        }
+        case actionType::opModify:
+        {
+            return fieldName_;
+        }
+    }
+
+    return word::null;
+}
+
+
+bool Foam::functionObjects::fvExpressionField::read(const dictionary& dict)
+{
+    fvMeshFunctionObject::read(dict);
+
+    action_ = actionNames_.getOrDefault("action", dict, actionType::opNew);
+
+    fieldName_ = dict.get<word>("field");
+    const word fldName = fieldName();
+
+    Log << type() << ' ' << this->name() << ':' << nl
+        << "    action  = " << actionNames_[action_] << nl
+        << "    field   = " << fldName << nl;
+
+    maskExpr_.clear();
+    valueExpr_.clear();
+
+    preloadFields_.clear();
+    dict.readIfPresent("readFields", preloadFields_);
+
+    switch (action_)
+    {
+        case actionType::opNone:
+        {
+            // No-op
+            break;
+        }
+        case actionType::opModify:
+        {
+            // Optional <fieldMask> for modify
+            maskExpr_.readEntry("fieldMask", dict, false);
+            [[fallthrough]];
+        }
+        case actionType::opNew:
+        {
+            // Mandatory <expression> for new and modify
+            valueExpr_.readEntry("expression", dict);
+            break;
+        }
+    }
+
+    autowrite_ = dict.getOrDefault("autowrite", false);
+    store_ = dict.getOrDefault("autowrite", true);
+
+    // "dimensions" is optional
+    dimensions_.clear();
+    hasDimensions_ = dimensions_.readEntry("dimensions", dict, false);
+
+    if (action_ == actionType::opNew)
+    {
+        if (!hasDimensions_)
+        {
+            Log << "    no 'dimensions' : treat '" << fldName
+                << "' as dimensionless" << endl;
+        }
+    }
+    else
+    {
+        // Ignore for none/modify
+        hasDimensions_  = false;
+    }
+
+
+    if (action_ == actionType::opNone)
+    {
+        driver_.reset(nullptr);
+        return true;  // Done
+    }
+
+    driver_.reset
+    (
+        new expressions::volumeExprDriver(mesh_, dict_)
+    );
+
+    driver_->setSearchBehaviour
+    (
+        expressions::exprDriver::searchControls
+        (
+            int(expressions::exprDriver::SEARCH_REGISTRY)
+          | (
+                loadFromFiles_
+              ? int(expressions::exprDriver::SEARCH_FILES)
+              : int(0)
+            )
+        ),
+        false  // No caching
+    );
+
+    driver_->readDict(dict_);
+
+    return true;
+}
+
+
+bool Foam::functionObjects::fvExpressionField::performAction(bool doWrite)
+{
+    using FieldAssociation = expressions::FieldAssociation;
+
+    if (!driver_ || action_ == actionType::opNone)
+    {
+        // No-op
+        return true;
+    }
+
+    const word fldName = fieldName();
+
+    if (loadFromFiles_)
+    {
+        loadFields(preloadFields_);
+    }
+
+    if (action_ == actionType::opModify && loadFromFiles_)
+    {
+        loadFields(wordList({fldName}));
+    }
+
+    auto& driver = *driver_;
+
+
+    // Current availability
+    auto* regIOobjectPtr = mesh_.getObjectPtr<regIOobject>(fldName);
+
+    if (action_ == actionType::opModify && !regIOobjectPtr)
+    {
+        // Cannot continue
+        FatalErrorInFunction
+            << type() << ' ' << this->name() << ':' << nl
+            << "    missing-field: " << fldName << nl
+            << exit(FatalError);
+
+        return false;
+    }
+
+
+    // Handle "field-mask" evaluation
+    bool evalFieldMask
+    (
+        (action_ == actionType::opModify)
+     && maskExpr_.size() && maskExpr_ != "true" && maskExpr_ != "1"
+    );
+
+    boolField fieldMask;
+    FieldAssociation maskFieldAssoc(FieldAssociation::NO_DATA);
+
+    if (evalFieldMask)
+    {
+        driver.parse(maskExpr_);
+
+        if (driver.isLogical())
+        {
+            auto& result = driver.result();
+            if (result.is_bool())
+            {
+                fieldMask = result.getResult<bool>();
+                maskFieldAssoc = driver.fieldAssociation();
+            }
+        }
+
+        // Slightly pedantic...
+        driver.clearField();
+        driver.clearResult();
+
+        evalFieldMask = (maskFieldAssoc != FieldAssociation::NO_DATA);
+
+        if (!evalFieldMask)
+        {
+            FatalErrorInFunction
+                << "field-mask: " << maskExpr_
+                << " does not evaluate to a logical expression: "
+                << driver.resultType() << nl
+                #ifdef FULLDEBUG
+                << "contents: " << fieldMask
+                #endif
+                << exit(FatalError);
+        }
+    }
+
+
+    // Start "expression" evaluation
+
+    bool applied = false;
+    autoPtr<regIOobject> toutputField;
+
+    {
+        driver.clearVariables();
+
+        driver.parse(valueExpr_);
+
+        if (evalFieldMask && maskFieldAssoc != driver.fieldAssociation())
+        {
+            FatalErrorInFunction
+                << "Mismatch between field-mask geometric type ("
+                << fieldGeoType(maskFieldAssoc) << ") and" << nl
+                << "expression geometric type ("
+                << fieldGeoType(driver.fieldAssociation()) << ')' << nl
+                << nl
+                << "Expression: " << valueExpr_ << nl
+                << "Field-mask: " << maskExpr_ << nl
+                << nl
+                << exit(FatalError);
+        }
+
+        // The output field does not appear to exist
+        // - create a new 'blank slate'
+        if (!regIOobjectPtr)
+        {
+            toutputField.reset(driver.dupZeroField());
+
+            if (toutputField)
+            {
+                toutputField->rename(fldName);
+
+                if (autowrite_)
+                {
+                    toutputField->writeOpt(IOobject::AUTO_WRITE);
+                }
+            }
+
+            if (!store_)
+            {
+                // Local (non-registered) field only
+                regIOobjectPtr = toutputField.get();
+            }
+            else
+            {
+                if (toutputField->checkIn() && toutputField->store())
+                {
+                    // Register and transfer ownership to registry
+                    toutputField.release();
+                }
+
+                regIOobjectPtr = mesh_.getObjectPtr<regIOobject>(fldName);
+            }
+        }
+
+
+        // Additional checks (TBD):
+
+        if (!regIOobjectPtr)
+        {
+            // Cannot continue
+            FatalErrorInFunction
+                << type() << ' ' << this->name() << ':' << nl
+                << "    missing-field: " << fldName << nl
+                << exit(FatalError);
+        }
+
+        // const word oldFieldType = regIOobjectPtr->type();
+
+        // if (driver.resultType() != oldFieldType)
+        // {
+        //     FatalErrorInFunction
+        //         << "Inconsistent types: " << fldName << " is  "
+        //         << oldFieldType
+        //         << " but the expression evaluates to "
+        //         << driver.resultType()
+        //         << exit(FatalError);
+        //     }
+
+        switch (driver.fieldAssociation())
+        {
+            #undef  doLocalCode
+            #define doLocalCode(GeoField)                                    \
+            {                                                                \
+                /* FieldType */                                              \
+                auto* outPtr = dynamic_cast<GeoField*>(regIOobjectPtr);      \
+                const auto* ptr = driver.isResultType<GeoField>();           \
+                                                                             \
+                if (outPtr && ptr)                                           \
+                {                                                            \
+                    applied = setField(*outPtr, *ptr, fieldMask);            \
+                    if (doWrite)                                             \
+                    {                                                        \
+                        outPtr->write();                                     \
+                    }                                                        \
+                    break;                                                   \
+                }                                                            \
+            }
+
+            case FieldAssociation::VOLUME_DATA:
+            {
+                doLocalCode(volScalarField);
+                doLocalCode(volVectorField);
+                doLocalCode(volTensorField);
+                doLocalCode(volSymmTensorField);
+                doLocalCode(volSphericalTensorField);
+                break;
+            }
+            case FieldAssociation::FACE_DATA:
+            {
+                doLocalCode(surfaceScalarField);
+                doLocalCode(surfaceVectorField);
+                doLocalCode(surfaceTensorField);
+                doLocalCode(surfaceSymmTensorField);
+                doLocalCode(surfaceSphericalTensorField);
+                break;
+            }
+            case FieldAssociation::POINT_DATA:
+            {
+                doLocalCode(pointScalarField);
+                doLocalCode(pointVectorField);
+                doLocalCode(pointTensorField);
+                doLocalCode(pointSymmTensorField);
+                doLocalCode(pointSphericalTensorField);
+                break;
+            }
+
+            default: break;
+            #undef doLocalCode
+        }
+    }
+
+
+    // Clear out heavier data
+    driver.clearResult();
+    driver.clearField();
+
+    if (!applied)
+    {
+        // Or error?
+        WarningInFunction
+            << type() << ' ' << this->name() << ": Failed to apply "
+            << actionNames_[action_] << " for " << fldName
+            << nl;
+    }
+
+    return true;
+}
+
+
+bool Foam::functionObjects::fvExpressionField::execute()
+{
+    return performAction(false);
+}
+
+
+bool Foam::functionObjects::fvExpressionField::write()
+{
+    return performAction(true);
+}
+
+
+// ************************************************************************* //
diff --git a/src/functionObjects/field/expressions/fvExpressionField.H b/src/functionObjects/field/expressions/fvExpressionField.H
new file mode 100644
index 0000000000..8b3f5b64fd
--- /dev/null
+++ b/src/functionObjects/field/expressions/fvExpressionField.H
@@ -0,0 +1,248 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2021 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+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::functionObjects::fvExpressionField
+
+Group
+    grpFieldFunctionObjects
+
+Description
+    Function object that generates or modifies a field based on expressions.
+
+Usage
+    A minimal example:
+    \verbatim
+    <name1>
+    {
+        type    exprField;
+        libs    (fieldFunctionObjects);
+        field   pTotal;
+
+        expression  "p + 0.5*(rho*magSqr(U))";
+        dimensions  [ Pa ];
+    }
+
+    // Modify an existing field
+    <name2>
+    {
+        type    exprField;
+        libs    (fieldFunctionObjects);
+        field   pTotal;
+        action  modify;
+
+        // Static pressure only in these regions
+        fieldMask   "(mag(pos()) < 0.05) && (pos().y() > 0)";
+        expression  "p";
+    }
+    \endverbatim
+
+    where the entries mean:
+    \table
+      Property     | Description                        | Type   | Reqd | Dflt
+      type         | Type name: exprField               | word   | yes  |
+      libs         | Libraries: fieldFunctionObjects    | wordList | yes |
+      field        | Name of input or output field      | word   | yes |
+      expression   | Field evaluation expression        | string | yes |
+      action       | Type of operation: see below       | word   | no  | new
+      autowrite    | Add AUTO_WRITE to created field    | bool   | no  | false
+      store        | Store calculated field             | bool   | no  | true
+      fieldMask    | Masking as logical expression      | string | no  | ""
+      dimensions   | Apply specified dimensions to created field | dim-spec | no |
+      readFields   | Preload named fields (post-process mode) | wordList | no |
+      useNamePrefix | Add prefix scoping to output name | bool   | no  | false
+    \endtable
+
+    Options for the \c action entry:
+    \plaintable
+      none          | No operation
+      new           | Define field based on expression (default)
+      modify        | Adjust field according to expression and fieldMask
+    \endplaintable
+
+Note
+    The \c useNamePrefix entry is always ignored for the \c modify action.
+
+SourceFiles
+    fvExpressionField.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef functionObjects_fvExpressionField_H
+#define functionObjects_fvExpressionField_H
+
+#include "fvMeshFunctionObject.H"
+#include "volumeExprDriver.H"
+#include "Enum.H"
+#include "dimensionSet.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace functionObjects
+{
+
+/*---------------------------------------------------------------------------*\
+                      Class fvExpressionField Declaration
+\*---------------------------------------------------------------------------*/
+
+class fvExpressionField
+:
+    public fvMeshFunctionObject
+{
+public:
+
+    // Public Data Types
+
+        //- Action type enumeration
+        enum actionType : unsigned char
+        {
+            opNone = 0,     //!< No operation
+            opNew,          //!< Create/overwrite field (default)
+            opModify        //!< Modify existing field
+        };
+
+        //- Action type names
+        static const Enum<actionType> actionNames_;
+
+
+protected:
+
+    // Private Data
+
+        //- The context dictionary
+        dictionary dict_;
+
+        //- Name of the field
+        word fieldName_;
+
+        //- Names fields to preload
+        wordList preloadFields_;
+
+        //- The field-mask expression (modify mode)
+        expressions::exprString maskExpr_;
+
+        //- Expression to evaluate
+        expressions::exprString valueExpr_;
+
+        //- Dimensions for new field
+        dimensionSet dimensions_;
+
+        //- Operation mode
+        actionType action_;
+
+        //- Set AUTO_WRITE for new field
+        bool autowrite_;
+
+        //- Store calculated field
+        bool store_;
+
+        //- True if dimensions_ should be used (creation)
+        bool hasDimensions_;
+
+        //- Load fields from files (not from objectRegistry)
+        bool loadFromFiles_;
+
+        autoPtr<expressions::volumeExprDriver> driver_;
+
+
+    // Private Member Functions
+
+        //- Attempt load from io, store on database if successful
+        template<class FieldType>
+        bool loadAndStore(const IOobject& io);
+
+        //- Forward to loadAndStore for supported types
+        template<class Type>
+        bool loadField(const IOobject& io);
+
+        //- Attempt to load specified fields
+        label loadFields(const UList<word>& fieldSet_);
+
+        template<class GeoField>
+        bool setField
+        (
+            GeoField& output,
+            const GeoField& evaluated,
+            const boolField& cond
+        );
+
+        bool performAction(bool doWrite);
+
+public:
+
+    //- Runtime type information
+    TypeName("exprField");
+
+
+    // Constructors
+
+        //- Construct from Time and dictionary
+        fvExpressionField
+        (
+            const word& name,
+            const Time& runTime,
+            const dictionary& dict,
+            const bool loadFromFiles = false
+        );
+
+        //- No copy construct
+        fvExpressionField(const fvExpressionField&) = delete;
+
+        //- No copy assignment
+        void operator=(const fvExpressionField&) = delete;
+
+
+    //- Destructor
+    virtual ~fvExpressionField();
+
+
+    // Member Functions
+
+        //- Qualified/unqualified field name (depends on action)
+        virtual word fieldName() const;
+
+        //- Read the data
+        virtual bool read(const dictionary& dict);
+
+        //- Execute
+        virtual bool execute();
+
+        //- Write
+        virtual bool write();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace functionObjects
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/tutorials/compressible/rhoSimpleFoam/squareBend/system/derivedFields b/tutorials/compressible/rhoSimpleFoam/squareBend/system/derivedFields
index e0ab565deb..829aff42e0 100644
--- a/tutorials/compressible/rhoSimpleFoam/squareBend/system/derivedFields
+++ b/tutorials/compressible/rhoSimpleFoam/squareBend/system/derivedFields
@@ -24,4 +24,42 @@ derivedFields
 }
 
 
+exprField1
+{
+    // Mandatory entries
+    type     exprField;
+    libs     (fieldFunctionObjects);
+    field    pTot;
+    /// action new;
+    ///autowrite       false;
+    ///store           true;
+    useNamePrefix   false;
+
+    // Use expression to define total pressure
+
+    executeControl  timeStep;
+    writeControl    none;
+
+    expression  "p + 0.5*(rho * magSqr(U))";
+    dimensions  [ Pa ];
+}
+
+
+exprField1.mod0
+{
+    // Mandatory entries
+    type    exprField;
+    libs    (fieldFunctionObjects);
+    field   pTot;
+    action  modify;
+
+    executeControl  timeStep;
+    writeControl    writeTime;
+
+    // Static pressure only in these regions
+    fieldMask   "(mag(pos()) < 0.05) && (pos().y() > 0) || cellZone(inlet)";
+    expression  "p";
+}
+
+
 // ************************************************************************* //
diff --git a/tutorials/compressible/rhoSimpleFoam/squareBend/system/sampling b/tutorials/compressible/rhoSimpleFoam/squareBend/system/sampling
index caaa2d8cca..85f696a6da 100644
--- a/tutorials/compressible/rhoSimpleFoam/squareBend/system/sampling
+++ b/tutorials/compressible/rhoSimpleFoam/squareBend/system/sampling
@@ -69,7 +69,7 @@ sampled
 
     surfaceFormat   none;
 
-    fields      (p rho U T rhoU pTotal);
+    fields      (p rho U T rhoU pTotal pTot);
 
     // Remove derived fields we created prior
     removeFieldsOnExecute   (rhoU pTotal);
@@ -182,7 +182,7 @@ areaAverage
     operation       weightedAreaAverage;
     weightField     rhoU;
     weightFields    ( rho U none );  // 2012 and later
-    fields          ( p pTotal );
+    fields          ( p pTotal pTot );
 }
 
 areaIntegrate