diff --git a/applications/utilities/postProcessing/graphics/PVReaders/PVFoamReader/CMakeLists.txt b/applications/utilities/postProcessing/graphics/PVReaders/PVFoamReader/CMakeLists.txt
index 4815a13d8d87539bd06373031cc139bb43120a68..18302332ae81e2acbcdb04b4cf74eaeb5ade3a62 100644
--- a/applications/utilities/postProcessing/graphics/PVReaders/PVFoamReader/CMakeLists.txt
+++ b/applications/utilities/postProcessing/graphics/PVReaders/PVFoamReader/CMakeLists.txt
@@ -13,6 +13,7 @@ include_directories(
     $ENV{WM_PROJECT_DIR}/src/OpenFOAM/lnInclude
     $ENV{WM_PROJECT_DIR}/src/OSspecific/$ENV{WM_OSTYPE}/lnInclude
     $ENV{WM_PROJECT_DIR}/src/conversion/lnInclude
+    $ENV{WM_PROJECT_DIR}/src/finiteArea/lnInclude
     $ENV{WM_PROJECT_DIR}/src/finiteVolume/lnInclude
     ${PROJECT_SOURCE_DIR}/../foamPv
     ${PROJECT_SOURCE_DIR}/../vtkPVFoam
@@ -66,6 +67,7 @@ target_link_libraries(
     vtkPVFoam-pv${PARAVIEW_VERSION_MAJOR}.${PARAVIEW_VERSION_MINOR}
     foamPv-pv${PARAVIEW_VERSION_MAJOR}.${PARAVIEW_VERSION_MINOR}
     conversion
+    finiteArea
     finiteVolume
     OpenFOAM
 )
diff --git a/applications/utilities/postProcessing/graphics/PVReaders/PVFoamReader/vtkPVFoamReader.cxx b/applications/utilities/postProcessing/graphics/PVReaders/PVFoamReader/vtkPVFoamReader.cxx
index c6ec1e676cd5fa567bee9614c8237c870ec17a6c..ac4df9384aba36969cc1f926a9244d5d99978910 100644
--- a/applications/utilities/postProcessing/graphics/PVReaders/PVFoamReader/vtkPVFoamReader.cxx
+++ b/applications/utilities/postProcessing/graphics/PVReaders/PVFoamReader/vtkPVFoamReader.cxx
@@ -309,16 +309,16 @@ int vtkPVFoamReader::RequestData
     // Get the requested time step.
     // We only support requests for a single time step
 
-    int nRequestTime = 0;
-    double requestTime[nInfo];
+    std::vector<double> requestTime;
+    requestTime.reserve(nInfo);
 
     // taking port0 as the lead for other outputs would be nice, but fails when
     // a filter is added - we need to check everything
     // but since PREVIOUS_UPDATE_TIME_STEPS() is protected, relay the logic
     // to the vtkPVFoam::setTime() method
-    for (int infoI = 0; infoI < nInfo; ++infoI)
+    for (int infoi = 0; infoi < nInfo; ++infoi)
     {
-        vtkInformation *outInfo = outputVector->GetInformationObject(infoI);
+        vtkInformation *outInfo = outputVector->GetInformationObject(infoi);
 
         const int nsteps =
             outInfo->Length(vtkStreamingDemandDrivenPipeline::TIME_STEPS());
@@ -329,7 +329,7 @@ int vtkPVFoamReader::RequestData
          && nsteps > 0
         )
         {
-            requestTime[nRequestTime] =
+            const double timeValue =
             (
                 1 == nsteps
                 // Only one time-step available, UPDATE_TIME_STEP is unreliable
@@ -337,16 +337,13 @@ int vtkPVFoamReader::RequestData
               : outInfo->Get(vtkStreamingDemandDrivenPipeline::UPDATE_TIME_STEP())
             );
 
-            // outInfo->Set(vtkDataObject::DATA_TIME_STEP(), requestTime[nRequestTime]);
-            // this->SetTimeValue(requestedTimeValue);
-            ++nRequestTime;
+            // outInfo->Set(vtkDataObject::DATA_TIME_STEP(), timeValue);
+            // this->SetTimeValue(timeValue);
+            requestTime.push_back(timeValue);
         }
     }
 
-    if (nRequestTime)
-    {
-        backend_->setTime(nRequestTime, requestTime);
-    }
+    backend_->setTime(requestTime);
 
     vtkMultiBlockDataSet* output = vtkMultiBlockDataSet::SafeDownCast
     (
diff --git a/applications/utilities/postProcessing/graphics/PVReaders/PVFoamReader/vtkPVFoamReader.h b/applications/utilities/postProcessing/graphics/PVReaders/PVFoamReader/vtkPVFoamReader.h
index 550c1167bbb887c1bfc8c81b607e3cf72e634b8a..1b1f4f3962a5ff0c1d09b59b5993b31153786779 100644
--- a/applications/utilities/postProcessing/graphics/PVReaders/PVFoamReader/vtkPVFoamReader.h
+++ b/applications/utilities/postProcessing/graphics/PVReaders/PVFoamReader/vtkPVFoamReader.h
@@ -25,9 +25,9 @@ Class
     vtkPVFoamReader
 
 Description
-    reads a dataset in OpenFOAM format
+    Reads a dataset in OpenFOAM format
 
-    vtkPVblockMeshReader creates an multiblock dataset.
+    vtkPVFoamReader creates an multiblock dataset.
     It uses the OpenFOAM infrastructure (fvMesh, etc) to handle mesh and
     field data.
 
diff --git a/applications/utilities/postProcessing/graphics/PVReaders/foamPv/foamPvCore.C b/applications/utilities/postProcessing/graphics/PVReaders/foamPv/foamPvCore.C
index 57e52af03dbed1396cf74dc95d89fa9a85c2e041..347ebf4383f7e4c203a7df00b1507722a6c2cdf4 100644
--- a/applications/utilities/postProcessing/graphics/PVReaders/foamPv/foamPvCore.C
+++ b/applications/utilities/postProcessing/graphics/PVReaders/foamPv/foamPvCore.C
@@ -43,6 +43,35 @@ namespace Foam
     defineTypeNameAndDebug(foamPvCore, 0);
 }
 
+
+// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
+
+Foam::Ostream& Foam::foamPvCore::printDataArraySelection
+(
+    Ostream& os,
+    vtkDataArraySelection* select
+)
+{
+    if (!select)
+    {
+        return os;
+    }
+
+    const int n = select->GetNumberOfArrays();
+
+    os << n << '(';
+    for (int i=0; i < n; ++i)
+    {
+        if (i) os << ' ';
+        os  << select->GetArrayName(i) << '='
+            << (select->GetArraySetting(i) ? 1 : 0);
+    }
+    os << ')';
+
+    return os;
+}
+
+
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 void Foam::foamPvCore::addToBlock
@@ -229,10 +258,8 @@ Foam::word Foam::foamPvCore::getFoamName(const std::string& str)
         // Already checked for valid/invalid chars
         return word(str.substr(beg, beg+end), false);
     }
-    else
-    {
-        return word::null;
-    }
+
+    return word::null;
 }
 
 
diff --git a/applications/utilities/postProcessing/graphics/PVReaders/foamPv/foamPvCore.H b/applications/utilities/postProcessing/graphics/PVReaders/foamPv/foamPvCore.H
index e7feca1f9009f5e25c6c6b35ed752f37a462c0b4..59ca049ea58fca4e9c826e0f37530ac071a1a7ce 100644
--- a/applications/utilities/postProcessing/graphics/PVReaders/foamPv/foamPvCore.H
+++ b/applications/utilities/postProcessing/graphics/PVReaders/foamPv/foamPvCore.H
@@ -119,13 +119,51 @@ public:
             setStart(startAt);
         }
 
-
         //- Increment the size
         void operator+=(label n)
         {
             setSize(size() + n);
         }
 
+        //- True if the labelRange intersects any key in the Map
+        template<class T>
+        bool intersects(const Map<T>& map) const
+        {
+            for (const label idx : *this)
+            {
+                if (map.found(idx))
+                {
+                    return true;
+                }
+            }
+
+            return false;
+        }
+
+        //- The intersection ids with keys in the Map
+        template<class T>
+        List<label> intersection(const Map<T>& map) const
+        {
+            List<label> indices(Foam::min(map.size(), this->size()));
+
+            if (indices.size())
+            {
+                label n = 0;
+
+                for (const label idx : *this)
+                {
+                    if (map.found(idx))
+                    {
+                        indices[n++] = idx;
+                    }
+                }
+
+                indices.setSize(n);
+            }
+
+            return indices;
+        }
+
     }; // End class arrayRange
 
 
@@ -139,6 +177,7 @@ private:
         //- Disallow default bitwise assignment
         void operator=(const foamPvCore&) = delete;
 
+
 public:
 
     //- Static data members
@@ -146,8 +185,15 @@ public:
 
 
     //- Construct null
-    foamPvCore()
-    {}
+    foamPvCore() = default;
+
+
+    //- Print information about vtkDataArraySelection
+    static Ostream& printDataArraySelection
+    (
+        Ostream& os,
+        vtkDataArraySelection* select
+    );
 
 
     //- Convenience method for the VTK multiblock API
diff --git a/applications/utilities/postProcessing/graphics/PVReaders/foamPv/foamPvCoreTemplates.C b/applications/utilities/postProcessing/graphics/PVReaders/foamPv/foamPvCoreTemplates.C
index cd9e2302b23d80aa907e3b197081ae550463ff36..8d40eee752adc6da7df56b6a8f3ec0ee0287831a 100644
--- a/applications/utilities/postProcessing/graphics/PVReaders/foamPv/foamPvCoreTemplates.C
+++ b/applications/utilities/postProcessing/graphics/PVReaders/foamPv/foamPvCoreTemplates.C
@@ -152,7 +152,8 @@ void Foam::foamPvCore::setSelectedArrayEntries
 )
 {
     const int n = select->GetNumberOfArrays();
-    // disable everything not explicitly enabled
+
+    // Disable everything not explicitly enabled
     select->DisableAllArrays();
 
     // Loop through entries, enabling as required
diff --git a/applications/utilities/postProcessing/graphics/PVReaders/foamPv/foamPvFields.H b/applications/utilities/postProcessing/graphics/PVReaders/foamPv/foamPvFields.H
index 69bb841fcf01a985ea170140403cc844ff7f32aa..725c8cf29d0804b50d9c1057fb3c918b4a884c95 100644
--- a/applications/utilities/postProcessing/graphics/PVReaders/foamPv/foamPvFields.H
+++ b/applications/utilities/postProcessing/graphics/PVReaders/foamPv/foamPvFields.H
@@ -29,7 +29,7 @@ Description
 #ifndef foamPvFields_H
 #define foamPvFields_H
 
-#include "volFields.H"
+#include "symmTensor.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -40,27 +40,17 @@ namespace Foam
                         Class foamPvFields Declaration
 \*---------------------------------------------------------------------------*/
 
-class foamPvFields
+struct foamPvFields
 {
-    // Private Member Functions
-
-        //- Disallow default bitwise copy construct
-        foamPvFields(const foamPvFields&) = delete;
-
-        //- Disallow default bitwise assignment
-        void operator=(const foamPvFields&) = delete;
-
-public:
-
     //- Remapping for some data types
     template<class Type>
     inline static void remapTuple(float vec[])
     {}
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+};
 
-}; // End class foamPvFields
 
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 // Template specialization for symmTensor
 template<>
diff --git a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/Make/options b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/Make/options
index 8c7fdce47261a7acbac1c83b7dbe1854937b2a83..ec74dc17d9a37d8d6158b4f3d54733ebed776fcc 100644
--- a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/Make/options
+++ b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/Make/options
@@ -4,6 +4,7 @@ EXE_INC = \
     ${c++LESSWARN} \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
+    -I$(LIB_SRC)/finiteArea/lnInclude \
     -I$(LIB_SRC)/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/lagrangian/basic/lnInclude \
     -I$(LIB_SRC)/fileFormats/lnInclude \
@@ -14,6 +15,7 @@ EXE_INC = \
     -I../PVFoamReader
 
 LIB_LIBS = \
+    -lfiniteArea \
     -ldynamicMesh \
     -lconversion \
     -lgenericPatchFields \
diff --git a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoam.C b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoam.C
index 67857bc0bf1a87d049bc137026b9ab263ad2f498..5da86f74521634301b5a315c7f5dce399070031a 100644
--- a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoam.C
+++ b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoam.C
@@ -27,9 +27,12 @@ License
 #include "vtkPVFoamReader.h"
 
 // OpenFOAM includes
+#include "areaFaMesh.H"
+#include "faMesh.H"
 #include "fvMesh.H"
 #include "Time.H"
 #include "patchZones.H"
+#include "IOobjectList.H"
 
 // VTK includes
 #include "vtkDataArraySelection.h"
@@ -108,7 +111,7 @@ void Foam::vtkPVFoam::resetCounters()
     // Reset array range information (ids and sizes)
     rangeVolume_.reset();
     rangePatches_.reset();
-    rangeLagrangian_.reset();
+    rangeClouds_.reset();
     rangeCellZones_.reset();
     rangeFaceZones_.reset();
     rangePointZones_.reset();
@@ -131,50 +134,49 @@ bool Foam::vtkPVFoam::addOutputBlock
     vtkSmartPointer<vtkMultiBlockDataSet> block;
     int datasetNo = 0;
 
-    for (auto partId : selector)
+    const List<label> partIds = selector.intersection(selectedPartIds_);
+
+    for (const auto partId : partIds)
     {
-        if (selectedPartIds_.found(partId))
+        const auto& longName = selectedPartIds_[partId];
+        const word shortName = getFoamName(longName);
+
+        auto iter = cache.find(longName);
+        if (iter.found() && iter.object().dataset)
         {
-            const auto& longName = selectedPartIds_[partId];
-            const word shortName = getFoamName(longName);
+            auto dataset = iter.object().dataset;
 
-            auto iter = cache.find(longName);
-            if (iter.found() && iter.object().dataset)
+            if (singleDataset)
             {
-                auto dataset = iter.object().dataset;
-
-                if (singleDataset)
-                {
-                    output->SetBlock(blockNo, dataset);
-                    output->GetMetaData(blockNo)->Set
-                    (
-                        vtkCompositeDataSet::NAME(),
-                        shortName.c_str()
-                    );
-
-                    ++datasetNo;
-                    break;
-                }
-                else if (datasetNo == 0)
-                {
-                    block = vtkSmartPointer<vtkMultiBlockDataSet>::New();
-                    output->SetBlock(blockNo, block);
-                    output->GetMetaData(blockNo)->Set
-                    (
-                        vtkCompositeDataSet::NAME(),
-                        selector.name()
-                    );
-                }
-
-                block->SetBlock(datasetNo, dataset);
-                block->GetMetaData(datasetNo)->Set
+                output->SetBlock(blockNo, dataset);
+                output->GetMetaData(blockNo)->Set
                 (
                     vtkCompositeDataSet::NAME(),
                     shortName.c_str()
                 );
 
                 ++datasetNo;
+                break;
             }
+            else if (datasetNo == 0)
+            {
+                block = vtkSmartPointer<vtkMultiBlockDataSet>::New();
+                output->SetBlock(blockNo, block);
+                output->GetMetaData(blockNo)->Set
+                (
+                    vtkCompositeDataSet::NAME(),
+                    selector.name()
+                );
+            }
+
+            block->SetBlock(datasetNo, dataset);
+            block->GetMetaData(datasetNo)->Set
+            (
+                vtkCompositeDataSet::NAME(),
+                shortName.c_str()
+            );
+
+            ++datasetNo;
         }
     }
 
@@ -182,17 +184,22 @@ bool Foam::vtkPVFoam::addOutputBlock
 }
 
 
-int Foam::vtkPVFoam::setTime(int nRequest, const double requestTimes[])
+int Foam::vtkPVFoam::setTime(const std::vector<double>& requestTimes)
 {
+    if (requestTimes.empty())
+    {
+        return 0;
+    }
+
     Time& runTime = dbPtr_();
 
     // Get times list
     instantList Times = runTime.times();
 
     int nearestIndex = timeIndex_;
-    for (int requestI = 0; requestI < nRequest; ++requestI)
+    for (const double& timeValue : requestTimes)
     {
-        int index = Time::findClosestTimeIndex(Times, requestTimes[requestI]);
+        const int index = Time::findClosestTimeIndex(Times, timeValue);
         if (index >= 0 && index != timeIndex_)
         {
             nearestIndex = index;
@@ -208,13 +215,15 @@ int Foam::vtkPVFoam::setTime(int nRequest, const double requestTimes[])
     if (debug)
     {
         Info<< "<beg> setTime(";
-        for (int requestI = 0; requestI < nRequest; ++requestI)
+        unsigned reqi = 0;
+        for (const double& timeValue : requestTimes)
         {
-            if (requestI) Info<< ", ";
-            Info<< requestTimes[requestI];
+            if (reqi) Info<< ", ";
+            Info<< timeValue;
+            ++reqi;
         }
         Info<< ") - previousIndex = " << timeIndex_
-            << ", nearestIndex = " << nearestIndex << endl;
+            << ", nearestIndex = " << nearestIndex << nl;
     }
 
     // See what has changed
@@ -223,8 +232,13 @@ int Foam::vtkPVFoam::setTime(int nRequest, const double requestTimes[])
         timeIndex_ = nearestIndex;
         runTime.setTime(Times[nearestIndex], nearestIndex);
 
-        // When the changes, so do the fields
-        meshState_ = meshPtr_ ? meshPtr_->readUpdate() : polyMesh::TOPO_CHANGE;
+        // When mesh changes, so do fields
+        meshState_ =
+        (
+            volMeshPtr_
+          ? volMeshPtr_->readUpdate()
+          : polyMesh::TOPO_CHANGE
+        );
 
         reader_->UpdateProgress(0.05);
 
@@ -237,7 +251,7 @@ int Foam::vtkPVFoam::setTime(int nRequest, const double requestTimes[])
         Info<< "<end> setTime() - selectedTime="
             << Times[nearestIndex].name() << " index=" << timeIndex_
             << "/" << Times.size()
-            << " meshUpdateState=" << updateStateName(meshState_) << endl;
+            << " meshUpdateState=" << updateStateName(meshState_) << nl;
     }
 
     return nearestIndex;
@@ -260,15 +274,17 @@ Foam::vtkPVFoam::vtkPVFoam
 :
     reader_(reader),
     dbPtr_(nullptr),
-    meshPtr_(nullptr),
+    volMeshPtr_(nullptr),
+    areaMeshPtr_(nullptr),
     meshRegion_(polyMesh::defaultRegion),
     meshDir_(polyMesh::meshSubDir),
     timeIndex_(-1),
     decomposePoly_(false),
     meshState_(polyMesh::TOPO_CHANGE),
-    rangeVolume_("unzoned"),
+    rangeVolume_("volMesh"),
+    rangeArea_("areaMesh"),
     rangePatches_("patch"),
-    rangeLagrangian_("lagrangian"),
+    rangeClouds_("lagrangian"),
     rangeCellZones_("cellZone"),
     rangeFaceZones_("faceZone"),
     rangePointZones_("pointZone"),
@@ -278,7 +294,7 @@ Foam::vtkPVFoam::vtkPVFoam
 {
     if (debug)
     {
-        Info<< "vtkPVFoam - " << FileName << endl;
+        Info<< "vtkPVFoam - " << FileName << nl;
         printMemory();
     }
 
@@ -315,13 +331,13 @@ Foam::vtkPVFoam::vtkPVFoam
     // could be stringent and insist the prefix match the directory name...
     // Note: cannot use fileName::name() due to the embedded '{}'
     string caseName(fileName(FileName).lessExt());
-    string::size_type beg = caseName.find_last_of("/{");
-    string::size_type end = caseName.find('}', beg);
+    const auto beg = caseName.find_last_of("/{");
+    const auto end = caseName.find('}', beg);
 
     if
     (
-        beg != string::npos && caseName[beg] == '{'
-     && end != string::npos && end == caseName.size()-1
+        beg != std::string::npos && caseName[beg] == '{'
+     && end != std::string::npos && end == caseName.size()-1
     )
     {
         meshRegion_ = caseName.substr(beg+1, end-beg-1);
@@ -343,7 +359,7 @@ Foam::vtkPVFoam::vtkPVFoam
         Info<< "fullCasePath=" << fullCasePath << nl
             << "FOAM_CASE=" << getEnv("FOAM_CASE") << nl
             << "FOAM_CASENAME=" << getEnv("FOAM_CASENAME") << nl
-            << "region=" << meshRegion_ << endl;
+            << "region=" << meshRegion_ << nl;
     }
 
     // Create time object
@@ -369,10 +385,11 @@ Foam::vtkPVFoam::~vtkPVFoam()
 {
     if (debug)
     {
-        Info<< "~vtkPVFoam" << endl;
+        Info<< "~vtkPVFoam" << nl;
     }
 
-    delete meshPtr_;
+    delete volMeshPtr_;
+    delete areaMeshPtr_;
 }
 
 
@@ -383,64 +400,62 @@ void Foam::vtkPVFoam::updateInfo()
     if (debug)
     {
         Info<< "<beg> updateInfo"
-            << " [meshPtr=" << (meshPtr_ ? "set" : "nullptr") << "] timeIndex="
-            << timeIndex_ << endl;
+            << " [volMeshPtr=" << (volMeshPtr_ ? "set" : "nullptr")
+            << "] timeIndex="
+            << timeIndex_ << nl;
     }
 
     resetCounters();
 
-    vtkDataArraySelection* partSelection = reader_->GetPartSelection();
+    // Part selection
+    {
+        vtkDataArraySelection* select = reader_->GetPartSelection();
 
-    // there are two ways to ensure we have the correct list of parts:
-    // 1. remove everything and then set particular entries 'on'
-    // 2. build a 'char **' list and call SetArraysWithDefault()
-    //
-    // Nr. 2 has the potential advantage of not touching the modification
-    // time of the vtkDataArraySelection, but the qt/paraview proxy
-    // layer doesn't care about that anyhow.
+        // There are two ways to ensure we have the correct list of parts:
+        // 1. remove everything and then set particular entries 'on'
+        // 2. build a 'char **' list and call SetArraysWithDefault()
+        //
+        // Nr. 2 has the potential advantage of not touching the modification
+        // time of the vtkDataArraySelection, but the qt/paraview proxy
+        // layer doesn't care about that anyhow.
 
-    HashSet<string> enabled;
-    if (!partSelection->GetNumberOfArrays() && !meshPtr_)
-    {
-        // Fake enable 'internalMesh' on the first call
-        enabled = { "internalMesh" };
-    }
-    else
-    {
-        // preserve the enabled selections
-        enabled = getSelectedArraySet(partSelection);
-    }
+        HashSet<string> enabled;
+        if (!select->GetNumberOfArrays() && !volMeshPtr_)
+        {
+            // Fake enable 'internalMesh' on the first call
+            enabled = { "internalMesh" };
+        }
+        else
+        {
+            // Preserve the enabled selections
+            enabled = getSelectedArraySet(select);
+        }
 
-    // Clear current mesh parts list
-    partSelection->RemoveAllArrays();
+        select->RemoveAllArrays();   // Clear existing list
 
-    // Update mesh parts list - add Lagrangian at the bottom
-    updateInfoInternalMesh(partSelection);
-    updateInfoPatches(partSelection, enabled);
-    updateInfoSets(partSelection);
-    updateInfoZones(partSelection);
-    updateInfoLagrangian(partSelection);
+        // Update mesh parts list - add Lagrangian at the bottom
+        updateInfoInternalMesh(select);
+        updateInfoPatches(select, enabled);
+        updateInfoSets(select);
+        updateInfoZones(select);
+        updateInfoAreaMesh(select);
+        updateInfoLagrangian(select);
 
-    // Adjust/restore the enabled selections
-    setSelectedArrayEntries(partSelection, enabled);
+        setSelectedArrayEntries(select, enabled);  // Adjust/restore selected
+    }
 
-    // Update volume, point and lagrangian fields
-    updateInfoFields<fvPatchField, volMesh>
-    (
-        reader_->GetVolFieldSelection()
-    );
-    updateInfoFields<pointPatchField, pointMesh>
-    (
-        reader_->GetPointFieldSelection()
-    );
-    updateInfoLagrangianFields
-    (
-        reader_->GetLagrangianFieldSelection()
-    );
+    // Volume and area fields - includes save/restore of selected
+    updateInfoContinuumFields(reader_->GetVolFieldSelection());
+
+    // Point fields - includes save/restore of selected
+    updateInfoPointFields(reader_->GetPointFieldSelection());
+
+    // Lagrangian fields - includes save/restore of selected
+    updateInfoLagrangianFields(reader_->GetLagrangianFieldSelection());
 
     if (debug)
     {
-        Info<< "<end> updateInfo" << endl;
+        Info<< "<end> updateInfo" << nl;
     }
 }
 
@@ -531,26 +546,29 @@ void Foam::vtkPVFoam::Update
     {
         if (debug)
         {
-            Info<< "<beg> updateFoamMesh" << endl;
+            Info<< "<beg> updateFoamMesh" << nl;
             printMemory();
         }
 
         if (!caching)
         {
-            delete meshPtr_;
-            meshPtr_ = nullptr;
+            delete volMeshPtr_;
+            delete areaMeshPtr_;
+
+            volMeshPtr_ = nullptr;
+            areaMeshPtr_ = nullptr;
         }
 
         // Check to see if the OpenFOAM mesh has been created
-        if (!meshPtr_)
+        if (!volMeshPtr_)
         {
             if (debug)
             {
                 Info<< "Creating OpenFOAM mesh for region " << meshRegion_
-                    << " at time=" << dbPtr_().timeName() << endl;
+                    << " at time=" << dbPtr_().timeName() << nl;
             }
 
-            meshPtr_ = new fvMesh
+            volMeshPtr_ = new fvMesh
             (
                 IOobject
                 (
@@ -567,13 +585,27 @@ void Foam::vtkPVFoam::Update
         {
             if (debug)
             {
-                Info<< "Using existing OpenFOAM mesh" << endl;
+                Info<< "Using existing OpenFOAM mesh" << nl;
             }
         }
 
+        if (rangeArea_.intersects(selectedPartIds_))
+        {
+            if (!areaMeshPtr_)
+            {
+                areaMeshPtr_ = new faMesh(*volMeshPtr_);
+            }
+        }
+        else
+        {
+            delete areaMeshPtr_;
+
+            areaMeshPtr_ = nullptr;
+        }
+
         if (debug)
         {
-            Info<< "<end> updateFoamMesh" << endl;
+            Info<< "<end> updateFoamMesh" << nl;
             printMemory();
         }
     }
@@ -600,6 +632,8 @@ void Foam::vtkPVFoam::Update
         reader_->UpdateProgress(0.7);
     }
 
+    convertMeshArea();
+
     convertMeshLagrangian();
 
     reader_->UpdateProgress(0.8);
@@ -607,8 +641,9 @@ void Foam::vtkPVFoam::Update
     // Update fields
     convertVolFields();
     convertPointFields();
-    convertLagrangianFields();
+    convertAreaFields();
 
+    convertLagrangianFields();
 
     // Assemble multiblock output
     addOutputBlock(output, cachedVtu_, rangeVolume_, true); // One dataset
@@ -619,16 +654,17 @@ void Foam::vtkPVFoam::Update
     addOutputBlock(output, cachedVtu_, rangeCellSets_);
     addOutputBlock(output, cachedVtp_, rangeFaceSets_);
     addOutputBlock(output, cachedVtp_, rangePointSets_);
+    addOutputBlock(output, cachedVtp_, rangeArea_);
     addOutputBlock
     (
         (outputLagrangian ? outputLagrangian : output),
         cachedVtp_,
-        rangeLagrangian_
+        rangeClouds_
     );
 
     if (debug)
     {
-        Info<< "done reader part" << nl << endl;
+        Info<< "done reader part" << nl << nl;
     }
     reader_->UpdateProgress(0.95);
 
@@ -655,8 +691,11 @@ void Foam::vtkPVFoam::UpdateFinalize()
 {
     if (!reader_->GetMeshCaching())
     {
-        delete meshPtr_;
-        meshPtr_ = nullptr;
+        delete volMeshPtr_;
+        delete areaMeshPtr_;
+
+        volMeshPtr_ = nullptr;
+        areaMeshPtr_ = nullptr;
     }
 
     reader_->UpdateProgress(1.0);
@@ -674,20 +713,20 @@ std::vector<double> Foam::vtkPVFoam::findTimes(const bool skipZero) const
 
         // find the first time for which this mesh appears to exist
         label begIndex = timeLst.size();
-        forAll(timeLst, timeI)
+        forAll(timeLst, timei)
         {
             if
             (
                 IOobject
                 (
                     "points",
-                    timeLst[timeI].name(),
+                    timeLst[timei].name(),
                     meshDir_,
                     runTime
                 ).typeHeaderOk<pointIOField>(false, false)
             )
             {
-                begIndex = timeI;
+                begIndex = timei;
                 break;
             }
         }
@@ -735,15 +774,15 @@ void Foam::vtkPVFoam::renderPatchNames
     const bool show
 )
 {
-    // always remove old actors first
+    // Always remove old actors first
 
-    forAll(patchTextActors_, patchi)
+    for (auto& actor : patchTextActors_)
     {
-        renderer->RemoveViewProp(patchTextActors_[patchi]);
+        renderer->RemoveViewProp(actor);
     }
     patchTextActors_.clear();
 
-    if (show && meshPtr_)
+    if (show && volMeshPtr_)
     {
         // get the display patches, strip off any prefix/suffix
         hashedWordList selectedPatches = getSelected
@@ -757,7 +796,7 @@ void Foam::vtkPVFoam::renderPatchNames
             return;
         }
 
-        const polyBoundaryMesh& pbMesh = meshPtr_->boundaryMesh();
+        const polyBoundaryMesh& pbMesh = volMeshPtr_->boundaryMesh();
 
         // Find the total number of zones
         // Each zone will take the patch name
@@ -784,19 +823,19 @@ void Foam::vtkPVFoam::renderPatchNames
 
             boolList featEdge(pp.nEdges(), false);
 
-            forAll(edgeFaces, edgeI)
+            forAll(edgeFaces, edgei)
             {
-                const labelList& eFaces = edgeFaces[edgeI];
+                const labelList& eFaces = edgeFaces[edgei];
 
                 if (eFaces.size() == 1)
                 {
                     // Note: could also do ones with > 2 faces but this gives
                     // too many zones for baffles
-                    featEdge[edgeI] = true;
+                    featEdge[edgei] = true;
                 }
                 else if (mag(n[eFaces[0]] & n[eFaces[1]]) < 0.5)
                 {
-                    featEdge[edgeI] = true;
+                    featEdge[edgei] = true;
                 }
             }
 
@@ -808,7 +847,7 @@ void Foam::vtkPVFoam::renderPatchNames
             labelList zoneNFaces(pZones.nZones(), 0);
 
             // Create storage for additional zone centres
-            forAll(zoneNFaces, zoneI)
+            forAll(zoneNFaces, zonei)
             {
                 zoneCentre[patchi].append(Zero);
             }
@@ -816,14 +855,14 @@ void Foam::vtkPVFoam::renderPatchNames
             // Do averaging per individual zone
             forAll(pp, facei)
             {
-                label zoneI = pZones[facei];
-                zoneCentre[patchi][zoneI] += pp[facei].centre(pp.points());
-                zoneNFaces[zoneI]++;
+                const label zonei = pZones[facei];
+                zoneCentre[patchi][zonei] += pp[facei].centre(pp.points());
+                zoneNFaces[zonei]++;
             }
 
-            forAll(zoneCentre[patchi], zoneI)
+            forAll(zoneCentre[patchi], zonei)
             {
-                zoneCentre[patchi][zoneI] /= zoneNFaces[zoneI];
+                zoneCentre[patchi][zonei] /= zoneNFaces[zonei];
             }
         }
 
@@ -842,7 +881,7 @@ void Foam::vtkPVFoam::renderPatchNames
         if (debug)
         {
             Info<< "displayed zone centres = " << displayZoneI << nl
-                << "zones per patch = " << nZones << endl;
+                << "zones per patch = " << nZones << nl;
         }
 
         // Set the size of the patch labels to max number of zones
@@ -850,7 +889,7 @@ void Foam::vtkPVFoam::renderPatchNames
 
         if (debug)
         {
-            Info<< "constructing patch labels" << endl;
+            Info<< "constructing patch labels" << nl;
         }
 
         // Actor index
@@ -875,7 +914,7 @@ void Foam::vtkPVFoam::renderPatchNames
                 {
                     Info<< "patch name = " << pp.name() << nl
                         << "anchor = " << zoneCentre[patchi][globalZoneI] << nl
-                        << "globalZoneI = " << globalZoneI << endl;
+                        << "globalZoneI = " << globalZoneI << nl;
                 }
 
                 // Into a list for later removal
@@ -892,9 +931,9 @@ void Foam::vtkPVFoam::renderPatchNames
     }
 
     // Add text to each renderer
-    forAll(patchTextActors_, actori)
+    for (auto& actor : patchTextActors_)
     {
-        renderer->AddViewProp(patchTextActors_[actori]);
+        renderer->AddViewProp(actor);
     }
 }
 
@@ -902,10 +941,10 @@ void Foam::vtkPVFoam::renderPatchNames
 void Foam::vtkPVFoam::PrintSelf(ostream& os, vtkIndent indent) const
 {
     os  << indent << "Number of nodes: "
-        << (meshPtr_ ? meshPtr_->nPoints() : 0) << "\n";
+        << (volMeshPtr_ ? volMeshPtr_->nPoints() : 0) << "\n";
 
     os  << indent << "Number of cells: "
-        << (meshPtr_ ? meshPtr_->nCells() : 0) << "\n";
+        << (volMeshPtr_ ? volMeshPtr_->nCells() : 0) << "\n";
 
     os  << indent << "Number of available time steps: "
         << (dbPtr_.valid() ? dbPtr_().times().size() : 0) << "\n";
@@ -918,8 +957,8 @@ void Foam::vtkPVFoam::printInfo() const
 {
     std::cout
         << "Region: " << meshRegion_ << "\n"
-        << "nPoints: " << (meshPtr_ ? meshPtr_->nPoints() : 0) << "\n"
-        << "nCells: "  << (meshPtr_ ? meshPtr_->nCells() : 0) << "\n"
+        << "nPoints: " << (volMeshPtr_ ? volMeshPtr_->nPoints() : 0) << "\n"
+        << "nCells: "  << (volMeshPtr_ ? volMeshPtr_->nCells() : 0) << "\n"
         << "nTimes: "
         << (dbPtr_.valid() ? dbPtr_().times().size() : 0) << "\n";
 
diff --git a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoam.H b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoam.H
index 226ff6b9bdb9c92a736c7fe15c0c7d8d9003be8b..3eb719499362248f37eb96086649e58de56bffd5 100644
--- a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoam.H
+++ b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoam.H
@@ -99,6 +99,7 @@ namespace Foam
 // OpenFOAM class forward declarations
 class argList;
 class Time;
+class faMesh;
 class fvMesh;
 class IOobjectList;
 class polyPatch;
@@ -127,9 +128,8 @@ class vtkPVFoam
         //  with the output fields.
         //  The original copy is reused for different timestep
         template<class DataType>
-        class foamVtkCaching
+        struct foamVtkCaching
         {
-        public:
             typedef DataType dataType;
 
             //- The geometry, without any cell/point data
@@ -208,7 +208,7 @@ class vtkPVFoam
 
 
         //- Bookkeeping for vtkPolyData
-        class foamVtpData
+        struct foamVtpData
         :
             public foamVtkCaching<vtkPolyData>,
             public foamVtkMeshMaps
@@ -216,7 +216,7 @@ class vtkPVFoam
 
 
         //- Bookkeeping for vtkUnstructuredGrid
-        class foamVtuData
+        struct foamVtuData
         :
             public foamVtkCaching<vtkUnstructuredGrid>,
             public foamVtkMeshMaps
@@ -231,8 +231,11 @@ class vtkPVFoam
         //- OpenFOAM time control
         autoPtr<Time> dbPtr_;
 
-        //- OpenFOAM mesh
-        fvMesh* meshPtr_;
+        //- OpenFOAM finite volume mesh
+        fvMesh* volMeshPtr_;
+
+        //- OpenFOAM finite area mesh
+        faMesh* areaMeshPtr_;
 
         //- The mesh region
         word meshRegion_;
@@ -262,8 +265,9 @@ class vtkPVFoam
         //  used to index into selectedPartIds and thus indirectly into
         //  cachedVtp, cachedVtu
         arrayRange rangeVolume_;
+        arrayRange rangeArea_;
         arrayRange rangePatches_;
-        arrayRange rangeLagrangian_;
+        arrayRange rangeClouds_;
         arrayRange rangeCellZones_;
         arrayRange rangeFaceZones_;
         arrayRange rangePointZones_;
@@ -295,6 +299,9 @@ class vtkPVFoam
         //- Internal mesh info
         void updateInfoInternalMesh(vtkDataArraySelection* select);
 
+        //- Finite area mesh info
+        void updateInfoAreaMesh(vtkDataArraySelection* select);
+
         //- Lagrangian info
         void updateInfoLagrangian(vtkDataArraySelection* select);
 
@@ -326,18 +333,28 @@ class vtkPVFoam
         template<template<class> class patchType, class meshType>
         void updateInfoFields
         (
-            vtkDataArraySelection* select
+            vtkDataArraySelection* select,
+            const IOobjectList& objects
         );
 
+        //- Volume/Area field info
+        void updateInfoContinuumFields(vtkDataArraySelection* select);
+
+        //- Point field info
+        void updateInfoPointFields(vtkDataArraySelection* select);
+
         //- Lagrangian field info
         void updateInfoLagrangianFields(vtkDataArraySelection* select);
 
 
       // Mesh conversion functions
 
-        //- Convert InternalMesh
+        //- Convert internalMesh
         void convertMeshVolume();
 
+        //- Convert areaMesh
+        void convertMeshArea();
+
         //- Convert Lagrangian points
         void convertMeshLagrangian();
 
@@ -490,9 +507,12 @@ class vtkPVFoam
         ) const;
 
 
-        //- Convert volume fields
+        //- Convert finite volume fields
         void convertVolFields();
 
+        //- Convert finite area fields
+        void convertAreaFields();
+
         //- Convert point fields
         void convertPointFields();
 
@@ -528,6 +548,14 @@ class vtkPVFoam
             const IOobjectList& objects
         );
 
+        //- Area fields - all types
+        template<class Type>
+        void convertAreaFields
+        (
+            const faMesh& mesh,
+            const IOobjectList& objects
+        );
+
         //- Volume field - all selected parts
         template<class Type>
         void convertVolFieldBlock
@@ -631,7 +659,7 @@ public:
         //- Set the runTime to the first plausible request time,
         //  returns the timeIndex
         //  sets to "constant" on error
-        int setTime(int count, const double requestTimes[]);
+        int setTime(const std::vector<double>& requestTimes);
 
 
         //- The current time index
diff --git a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamFieldTemplates.C b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamFieldTemplates.C
index fdab4abe23ba496bce96e8885cb57cc2e3f22ba4..5a2ddfcb50d90deff91a9386b85eb1b194512ef0 100644
--- a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamFieldTemplates.C
+++ b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamFieldTemplates.C
@@ -38,6 +38,8 @@ InClass
 #include "zeroGradientFvPatchField.H"
 #include "interpolatePointToCell.H"
 #include "foamPvFields.H"
+#include "areaFaMesh.H"
+#include "areaFields.H"
 
 // vtk includes
 #include "vtkFloatArray.h"
@@ -69,7 +71,7 @@ void Foam::vtkPVFoam::convertVolField
     {
         if (debug)
         {
-            Info<< "convertVolField interpolating:" << fld.name() << endl;
+            Info<< "convertVolField interpolating:" << fld.name() << nl;
         }
 
         ptfPtr.reset
@@ -83,12 +85,12 @@ void Foam::vtkPVFoam::convertVolField
     convertVolFieldBlock(fld, ptfPtr, rangeCellSets_);  // cellSets
 
     // Patches - currently skip field conversion for groups
-    for (auto partId : rangePatches_)
+    for
+    (
+        const auto partId
+      : rangePatches_.intersection(selectedPartIds_)
+    )
     {
-        if (!selectedPartIds_.found(partId))
-        {
-            continue;
-        }
         const auto& longName = selectedPartIds_[partId];
 
         auto iter = cachedVtp_.find(longName);
@@ -97,10 +99,12 @@ void Foam::vtkPVFoam::convertVolField
             // Should not happen, but for safety require a vtk geometry
             continue;
         }
+
         foamVtpData& vtpData = iter.object();
         auto dataset = vtpData.dataset;
 
         const labelList& patchIds = vtpData.additionalIds();
+
         if (patchIds.empty())
         {
             continue;
@@ -184,12 +188,12 @@ void Foam::vtkPVFoam::convertVolField
 
 
     // Face Zones
-    for (auto partId : rangeFaceZones_)
+    for
+    (
+        const auto partId
+      : rangeFaceZones_.intersection(selectedPartIds_)
+    )
     {
-        if (!selectedPartIds_.found(partId))
-        {
-            continue;
-        }
         const auto& longName = selectedPartIds_[partId];
         const word zoneName = getFoamName(longName);
 
@@ -199,6 +203,7 @@ void Foam::vtkPVFoam::convertVolField
             // Should not happen, but for safety require a vtk geometry
             continue;
         }
+
         foamVtpData& vtpData = iter.object();
         auto dataset = vtpData.dataset;
 
@@ -224,12 +229,12 @@ void Foam::vtkPVFoam::convertVolField
 
 
     // Face Sets
-    for (auto partId : rangeFaceSets_)
+    for
+    (
+        const auto partId
+      : rangeFaceSets_.intersection(selectedPartIds_)
+    )
     {
-        if (!selectedPartIds_.found(partId))
-        {
-            continue;
-        }
         const auto& longName = selectedPartIds_[partId];
         const word selectName = getFoamName(longName);
 
@@ -349,12 +354,12 @@ void Foam::vtkPVFoam::convertVolFieldBlock
     const arrayRange& range
 )
 {
-    for (auto partId : range)
+    for
+    (
+        const auto partId
+      : range.intersection(selectedPartIds_)
+    )
     {
-        if (!selectedPartIds_.found(partId))
-        {
-            continue;
-        }
         const auto& longName = selectedPartIds_[partId];
 
         auto iter = cachedVtu_.find(longName);
@@ -363,6 +368,7 @@ void Foam::vtkPVFoam::convertVolFieldBlock
             // Should not happen, but for safety require a vtk geometry
             continue;
         }
+
         foamVtuData& vtuData = iter.object();
         auto dataset = vtuData.dataset;
 
@@ -387,6 +393,65 @@ void Foam::vtkPVFoam::convertVolFieldBlock
 }
 
 
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+//
+// area-fields
+//
+
+template<class Type>
+void Foam::vtkPVFoam::convertAreaFields
+(
+    const faMesh& mesh,
+    const IOobjectList& objects
+)
+{
+    typedef GeometricField<Type, faPatchField, areaMesh> FieldType;
+
+    const List<label> partIds =
+        rangeArea_.intersection(selectedPartIds_);
+
+    if (partIds.empty())
+    {
+        return;
+    }
+
+    forAllConstIters(objects, iter)
+    {
+        // Restrict to GeometricField<Type, ...>
+        const auto& ioobj = *(iter.object());
+
+        if (ioobj.headerClassName() == FieldType::typeName)
+        {
+            // Load field
+            FieldType fld(ioobj, mesh);
+
+            // Convert
+            for (const auto partId : partIds)
+            {
+                const auto& longName = selectedPartIds_[partId];
+
+                auto iter = cachedVtp_.find(longName);
+                if (!iter.found() || !iter.object().dataset)
+                {
+                    // Should not happen, but for safety require a vtk geometry
+                    continue;
+                }
+
+                foamVtpData& vtpData = iter.object();
+                auto dataset = vtpData.dataset;
+
+                vtkSmartPointer<vtkFloatArray> cdata = convertFieldToVTK
+                (
+                    fld.name(),
+                    fld
+                );
+                dataset->GetCellData()->AddArray(cdata);
+            }
+        }
+    }
+}
+
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 //
@@ -416,7 +481,7 @@ void Foam::vtkPVFoam::convertPointFields
 
         if (debug)
         {
-            Info<< "convertPointFields : " << fieldName << endl;
+            Info<< "convertPointFields : " << fieldName << nl;
         }
 
         FieldType pfld(ioobj, pMesh);
@@ -426,12 +491,12 @@ void Foam::vtkPVFoam::convertPointFields
         convertPointFieldBlock(pfld, rangeCellSets_);   // cellSets
 
         // Patches - currently skip field conversion for groups
-        for (auto partId : rangePatches_)
+        for
+        (
+            const auto partId
+          : rangePatches_.intersection(selectedPartIds_)
+        )
         {
-            if (!selectedPartIds_.found(partId))
-            {
-                continue;
-            }
             const auto& longName = selectedPartIds_[partId];
 
             auto iter = cachedVtp_.find(longName);
@@ -440,6 +505,7 @@ void Foam::vtkPVFoam::convertPointFields
                 // Should not happen, but for safety require a vtk geometry
                 continue;
             }
+
             foamVtpData& vtpData = iter.object();
             auto dataset = vtpData.dataset;
 
@@ -461,12 +527,12 @@ void Foam::vtkPVFoam::convertPointFields
         }
 
         // Face Zones
-        for (auto partId : rangeFaceZones_)
+        for
+        (
+            const auto partId
+          : rangeFaceZones_.intersection(selectedPartIds_)
+        )
         {
-            if (!selectedPartIds_.found(partId))
-            {
-                continue;
-            }
             const auto& longName = selectedPartIds_[partId];
             const word zoneName = getFoamName(longName);
 
@@ -476,6 +542,7 @@ void Foam::vtkPVFoam::convertPointFields
                 // Should not happen, but for safety require a vtk geometry
                 continue;
             }
+
             foamVtpData& vtpData = iter.object();
             auto dataset = vtpData.dataset;
 
@@ -513,12 +580,12 @@ void Foam::vtkPVFoam::convertPointFieldBlock
     const arrayRange& range
 )
 {
-    for (auto partId : range)
+    for
+    (
+        const auto partId
+      : range.intersection(selectedPartIds_)
+    )
     {
-        if (!selectedPartIds_.found(partId))
-        {
-            continue;
-        }
         const auto& longName = selectedPartIds_[partId];
 
         auto iter = cachedVtu_.find(longName);
@@ -527,6 +594,7 @@ void Foam::vtkPVFoam::convertPointFieldBlock
             // Should not happen, but for safety require a vtk geometry
             continue;
         }
+
         foamVtuData& vtuData = iter.object();
         auto dataset = vtuData.dataset;
 
@@ -579,7 +647,7 @@ vtkSmartPointer<vtkFloatArray> Foam::vtkPVFoam::convertPointField
             << pfld.name()
             << " size="  << (nPoints + addPointCellLabels.size())
             << " (" << nPoints << " + " << addPointCellLabels.size()
-            << ") nComp=" << nComp << endl;
+            << ") nComp=" << nComp << nl;
     }
 
     float vec[nComp];
@@ -728,7 +796,7 @@ Foam::label Foam::vtkPVFoam::transcribeFloatData
             << pTraits<Type>::typeName
             << "' : target array has " << array->GetNumberOfComponents()
             << " components instead of " << nComp
-            << endl;
+            << nl;
     }
 
     const vtkIdType maxSize = array->GetNumberOfTuples();
@@ -740,7 +808,7 @@ Foam::label Foam::vtkPVFoam::transcribeFloatData
             << "vtk array '" << array->GetName()
             << "' copy with out-of-range (0-" << long(maxSize-1) << ")"
             << " starting location"
-            << endl;
+            << nl;
 
         return 0;
     }
@@ -751,7 +819,7 @@ Foam::label Foam::vtkPVFoam::transcribeFloatData
             << "' copy ends out-of-range (" << long(maxSize) << ")"
             << " using sizing (start,size) = ("
             << start << "," << input.size() << ")"
-            << endl;
+            << nl;
 
         return 0;
     }
@@ -788,7 +856,7 @@ Foam::vtkPVFoam::convertFieldToVTK
         Info<< "convert UList<" << pTraits<Type>::typeName << "> "
             << name
             << " size="  << fld.size()
-            << " nComp=" << nComp << endl;
+            << " nComp=" << nComp << nl;
     }
 
     auto data = vtkSmartPointer<vtkFloatArray>::New();
@@ -816,7 +884,7 @@ Foam::vtkPVFoam::convertFaceFieldToVTK
         Info<< "convert face field: "
             << fld.name()
             << " size="  << faceLabels.size()
-            << " nComp=" << int(pTraits<Type>::nComponents) << endl;
+            << " nComp=" << int(pTraits<Type>::nComponents) << nl;
     }
 
     const fvMesh& mesh = fld.mesh();
@@ -887,7 +955,7 @@ Foam::vtkPVFoam::convertVolFieldToVTK
             << " size=" << cellMap.size()
             << " (" << fld.size() << " + "
             << (cellMap.size() - fld.size())
-            << ") nComp=" << nComp << endl;
+            << ") nComp=" << nComp << nl;
     }
 
     float scratch[nComp];
@@ -906,6 +974,7 @@ Foam::vtkPVFoam::convertVolFieldToVTK
     return data;
 }
 
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 #endif
diff --git a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamFields.C b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamFields.C
index 9da212916f2e3ab1d1b0be919c5916f277155b16..a51295b3bf34da446039b83b995d71d740812c56 100644
--- a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamFields.C
+++ b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamFields.C
@@ -26,6 +26,7 @@ License
 #include "vtkPVFoam.H"
 
 // OpenFOAM includes
+#include "areaFaMesh.H"
 #include "Cloud.H"
 #include "IOobjectList.H"
 #include "vtkPVFoamReader.h"
@@ -42,7 +43,7 @@ License
 
 void Foam::vtkPVFoam::convertVolFields()
 {
-    const fvMesh& mesh = *meshPtr_;
+    const fvMesh& mesh = *volMeshPtr_;
 
     const bool interpFields = reader_->GetInterpolateVolFields();
     hashedWordList selectedFields = getSelected
@@ -67,11 +68,11 @@ void Foam::vtkPVFoam::convertVolFields()
 
     if (debug)
     {
-        Info<< "<beg> " << FUNCTION_NAME << endl;
+        Info<< "<beg> " << FUNCTION_NAME << nl;
         forAllConstIters(objects, iter)
         {
             Info<< "  " << iter()->name()
-                << " == " << iter()->objectPath() << endl;
+                << " == " << iter()->objectPath() << nl;
         }
         printMemory();
     }
@@ -110,7 +111,7 @@ void Foam::vtkPVFoam::convertVolFields()
 
     if (debug)
     {
-        Info<< "<end> " << FUNCTION_NAME << endl;
+        Info<< "<end> " << FUNCTION_NAME << nl;
         printMemory();
     }
 }
@@ -118,7 +119,7 @@ void Foam::vtkPVFoam::convertVolFields()
 
 void Foam::vtkPVFoam::convertPointFields()
 {
-    const fvMesh& mesh = *meshPtr_;
+    const fvMesh& mesh = *volMeshPtr_;
 
     hashedWordList selectedFields = getSelected
     (
@@ -129,7 +130,7 @@ void Foam::vtkPVFoam::convertPointFields()
     {
         if (debug)
         {
-            Info<< "no point fields selected" << endl;
+            Info<< "no point fields selected" << nl;
         }
         return;
     }
@@ -146,11 +147,11 @@ void Foam::vtkPVFoam::convertPointFields()
 
     if (debug)
     {
-        Info<< "<beg> convert volume -> point fields" << endl;
+        Info<< "<beg> convert volume -> point fields" << nl;
         forAllConstIters(objects, iter)
         {
             Info<< "  " << iter()->name()
-                << " == " << iter()->objectPath() << endl;
+                << " == " << iter()->objectPath() << nl;
         }
         printMemory();
     }
@@ -166,7 +167,61 @@ void Foam::vtkPVFoam::convertPointFields()
 
     if (debug)
     {
-        Info<< "<end> convert volume -> point fields" << endl;
+        Info<< "<end> convert volume -> point fields" << nl;
+        printMemory();
+    }
+}
+
+
+void Foam::vtkPVFoam::convertAreaFields()
+{
+    if (!areaMeshPtr_)
+    {
+        return;
+    }
+
+    const faMesh& mesh = *areaMeshPtr_;
+
+    vtkDataArraySelection* select = reader_->GetVolFieldSelection();
+
+    hashedWordList selectedFields = getSelected(select);
+
+    if (selectedFields.empty())
+    {
+        return;
+    }
+
+    // Get objects (fields) for this time - only keep selected fields
+    // the region name is already in the mesh db
+    IOobjectList objects(mesh.mesh(), dbPtr_().timeName());
+
+    objects.filterKeys(selectedFields);
+
+    if (objects.empty())
+    {
+        return;
+    }
+
+    if (debug)
+    {
+        Info<< "<beg> " << FUNCTION_NAME << nl;
+        forAllConstIters(objects, iter)
+        {
+            Info<< "  " << iter()->name()
+                << " == " << iter()->objectPath() << nl;
+        }
+        printMemory();
+    }
+
+    convertAreaFields<scalar>(mesh, objects);
+    convertAreaFields<vector>(mesh, objects);
+    convertAreaFields<sphericalTensor>(mesh, objects);
+    convertAreaFields<symmTensor>(mesh, objects);
+    convertAreaFields<tensor>(mesh, objects);
+
+    if (debug)
+    {
+        Info<< "<end> " << FUNCTION_NAME << nl;
         printMemory();
     }
 }
@@ -174,8 +229,10 @@ void Foam::vtkPVFoam::convertPointFields()
 
 void Foam::vtkPVFoam::convertLagrangianFields()
 {
-    const arrayRange& range = rangeLagrangian_;
-    const fvMesh& mesh = *meshPtr_;
+    const List<label> partIds =
+        rangeClouds_.intersection(selectedPartIds_);
+
+    const fvMesh& mesh = *volMeshPtr_;
 
     hashedWordList selectedFields = getSelected
     (
@@ -189,17 +246,12 @@ void Foam::vtkPVFoam::convertLagrangianFields()
 
     if (debug)
     {
-        Info<< "<beg> " << FUNCTION_NAME << endl;
+        Info<< "<beg> " << FUNCTION_NAME << nl;
         printMemory();
     }
 
-    for (auto partId : range)
+    for (const auto partId : partIds)
     {
-        if (!selectedPartIds_.found(partId))
-        {
-            continue;
-        }
-
         const auto& longName = selectedPartIds_[partId];
         const word cloudName = getFoamName(longName);
 
@@ -229,11 +281,11 @@ void Foam::vtkPVFoam::convertLagrangianFields()
 
         if (debug)
         {
-            Info<< "converting OpenFOAM lagrangian fields" << endl;
+            Info<< "converting OpenFOAM lagrangian fields" << nl;
             forAllConstIters(objects, iter)
             {
                 Info<< "  " << iter()->name()
-                    << " == " << iter()->objectPath() << endl;
+                    << " == " << iter()->objectPath() << nl;
             }
         }
 
@@ -247,7 +299,7 @@ void Foam::vtkPVFoam::convertLagrangianFields()
 
     if (debug)
     {
-        Info<< "<end> " << FUNCTION_NAME << endl;
+        Info<< "<end> " << FUNCTION_NAME << nl;
         printMemory();
     }
 }
diff --git a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamMesh.C b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamMesh.C
index 3ddd5d09569a9e533155a441273a2c48999b59b8..1b873f397209f9c3c97a8032f20278aae9e59bc7 100644
--- a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamMesh.C
+++ b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamMesh.C
@@ -29,6 +29,7 @@ License
 #include "cellSet.H"
 #include "faceSet.H"
 #include "pointSet.H"
+#include "faMesh.H"
 #include "fvMeshSubset.H"
 #include "vtkPVFoamReader.h"
 #include "uindirectPrimitivePatch.H"
@@ -44,60 +45,58 @@ License
 
 void Foam::vtkPVFoam::convertMeshVolume()
 {
-    const arrayRange& range = rangeVolume_;
-    const fvMesh& mesh = *meshPtr_;
+    // Convert internalMesh - actually only one part
+    const List<label> partIds =
+        rangeVolume_.intersection(selectedPartIds_);
+
+    const fvMesh& mesh = *volMeshPtr_;
 
     if (debug)
     {
-        Info<< "<beg> " << FUNCTION_NAME << endl;
+        Info<< "<beg> " << FUNCTION_NAME << nl;
         printMemory();
     }
 
-    // Convert the internalMesh
-    // this looks like more than one part, but it isn't
-    for (auto partId : range)
+    for (const auto partId : partIds)
     {
-        if (selectedPartIds_.found(partId))
-        {
-            const auto& longName = selectedPartIds_[partId];
-            foamVtuData& vtuData = cachedVtu_(longName);
+        const auto& longName = selectedPartIds_[partId];
+        foamVtuData& vtuData = cachedVtu_(longName);
 
-            vtkSmartPointer<vtkUnstructuredGrid> vtkgeom;
-            if (vtuData.nPoints())
+        vtkSmartPointer<vtkUnstructuredGrid> vtkgeom;
+        if (vtuData.nPoints())
+        {
+            if (meshState_ == polyMesh::UNCHANGED)
             {
-                if (meshState_ == polyMesh::UNCHANGED)
-                {
-                    if (debug)
-                    {
-                        Info<< "reuse " << longName << nl;
-                    }
-                    vtuData.reuse(); // reuse
-                    continue;
-                }
-                else if (meshState_ == polyMesh::POINTS_MOVED)
+                if (debug)
                 {
-                    if (debug)
-                    {
-                        Info<< "move points " << longName << nl;
-                    }
-                    vtkgeom = vtuData.getCopy();
-                    vtkgeom->SetPoints(movePoints(mesh, vtuData));
+                    Info<< "reuse " << longName << nl;
                 }
+                vtuData.reuse(); // reuse
+                continue;
             }
-
-            if (!vtkgeom)
+            else if (meshState_ == polyMesh::POINTS_MOVED)
             {
-                // Nothing usable from cache - create new geometry
-                vtkgeom = volumeVTKMesh(mesh, vtuData);
+                if (debug)
+                {
+                    Info<< "move points " << longName << nl;
+                }
+                vtkgeom = vtuData.getCopy();
+                vtkgeom->SetPoints(movePoints(mesh, vtuData));
             }
+        }
 
-            vtuData.set(vtkgeom);
+        if (!vtkgeom)
+        {
+            // Nothing usable from cache - create new geometry
+            vtkgeom = volumeVTKMesh(mesh, vtuData);
         }
+
+        vtuData.set(vtkgeom);
     }
 
     if (debug)
     {
-        Info<< "<end> " << FUNCTION_NAME << endl;
+        Info<< "<end> " << FUNCTION_NAME << nl;
         printMemory();
     }
 }
@@ -105,33 +104,32 @@ void Foam::vtkPVFoam::convertMeshVolume()
 
 void Foam::vtkPVFoam::convertMeshLagrangian()
 {
-    const arrayRange& range = rangeLagrangian_;
-    const fvMesh& mesh = *meshPtr_;
+    const List<label> partIds =
+        rangeClouds_.intersection(selectedPartIds_);
+
+    const fvMesh& mesh = *volMeshPtr_;
 
     if (debug)
     {
-        Info<< "<beg> " << FUNCTION_NAME << endl;
+        Info<< "<beg> " << FUNCTION_NAME << nl;
         printMemory();
     }
 
-    for (auto partId : range)
+    for (const auto partId : partIds)
     {
-        if (selectedPartIds_.found(partId))
-        {
-            const auto& longName = selectedPartIds_[partId];
-            const word cloudName = getFoamName(longName);
-
-            // Direct conversion, no caching for Lagrangian
-            cachedVtp_(longName).set
-            (
-                lagrangianVTKMesh(mesh, cloudName)
-            );
-        }
+        const auto& longName = selectedPartIds_[partId];
+        const word cloudName = getFoamName(longName);
+
+        // Direct conversion, no caching for Lagrangian
+        cachedVtp_(longName).set
+        (
+            lagrangianVTKMesh(mesh, cloudName)
+        );
     }
 
     if (debug)
     {
-        Info<< "<end> " << FUNCTION_NAME << endl;
+        Info<< "<end> " << FUNCTION_NAME << nl;
         printMemory();
     }
 }
@@ -139,22 +137,20 @@ void Foam::vtkPVFoam::convertMeshLagrangian()
 
 void Foam::vtkPVFoam::convertMeshPatches()
 {
-    const arrayRange& range = rangePatches_;
-    const fvMesh& mesh = *meshPtr_;
+    const List<label> partIds =
+        rangePatches_.intersection(selectedPartIds_);
+
+    const fvMesh& mesh = *volMeshPtr_;
     const polyBoundaryMesh& patches = mesh.boundaryMesh();
 
     if (debug)
     {
-        Info<< "<beg> " << FUNCTION_NAME << endl;
+        Info<< "<beg> " << FUNCTION_NAME << nl;
         printMemory();
     }
 
-    for (auto partId : range)
+    for (const auto partId : partIds)
     {
-        if (!selectedPartIds_.found(partId))
-        {
-            continue;
-        }
         const auto& longName = selectedPartIds_[partId];
         const word partName = getFoamName(longName);
         foamVtpData& vtpData = cachedVtp_(longName);
@@ -198,14 +194,14 @@ void Foam::vtkPVFoam::convertMeshPatches()
             if (debug)
             {
                 Info<< "Creating VTK mesh for patches [" << patchIds <<"] "
-                    << longName << endl;
+                    << longName << nl;
             }
 
             // Store good patch ids as additionalIds
             vtpData.additionalIds().reserve(patchIds.size());
 
             label sz = 0;
-            for (auto id : patchIds)
+            for (const auto id : patchIds)
             {
                 const label n = patches[id].size();
 
@@ -221,7 +217,7 @@ void Foam::vtkPVFoam::convertMeshPatches()
             DynamicList<label>& faceLabels = vtpData.cellMap();
             faceLabels.reserve(sz);
 
-            for (auto id : vtpData.additionalIds())
+            for (const auto id : vtpData.additionalIds())
             {
                 const auto& pp = patches[id];
                 forAll(pp, i)
@@ -252,7 +248,7 @@ void Foam::vtkPVFoam::convertMeshPatches()
             if (debug)
             {
                 Info<< "Creating VTK mesh for patch [" << patchId <<"] "
-                    << partName << endl;
+                    << partName << nl;
             }
 
             if (patchId >= 0)
@@ -276,36 +272,112 @@ void Foam::vtkPVFoam::convertMeshPatches()
 
     if (debug)
     {
-        Info<< "<end> " << FUNCTION_NAME << endl;
+        Info<< "<end> " << FUNCTION_NAME << nl;
         printMemory();
     }
 }
 
 
-void Foam::vtkPVFoam::convertMeshCellZones()
+void Foam::vtkPVFoam::convertMeshArea()
 {
-    const arrayRange& range = rangeCellZones_;
-    const fvMesh& mesh = *meshPtr_;
+    // Convert areaMesh - actually only one part
+    const List<label> partIds =
+        rangeArea_.intersection(selectedPartIds_);
 
-    if (range.empty())
+    if (partIds.empty())
     {
         return;
     }
+    else if (!areaMeshPtr_ && volMeshPtr_)
+    {
+        areaMeshPtr_ = new faMesh(*volMeshPtr_);
+    }
+
+    if (!areaMeshPtr_)
+    {
+        return;
+    }
+
+    const faMesh& mesh = *areaMeshPtr_;
 
     if (debug)
     {
-        Info<< "<beg> " << FUNCTION_NAME << endl;
+        Info<< "<beg> " << FUNCTION_NAME << nl;
         printMemory();
     }
 
-    const cellZoneMesh& zMesh = mesh.cellZones();
-    for (auto partId : range)
+    for (const auto partId : partIds)
     {
-        if (!selectedPartIds_.found(partId))
+        const auto& longName = selectedPartIds_[partId];
+        foamVtpData& vtpData = cachedVtp_(longName);
+
+        vtkSmartPointer<vtkPolyData> vtkgeom;
+
+// This needs checking
+// #if 0
+//         if (vtpData.nPoints())
+//         {
+//             if (meshState_ == polyMesh::UNCHANGED)
+//             {
+//                 if (debug)
+//                 {
+//                     Info<< "reuse " << longName << nl;
+//                 }
+//                 vtpData.reuse(); // reuse
+//                 continue;
+//             }
+//             else if (meshState_ == polyMesh::POINTS_MOVED)
+//             {
+//                 if (debug)
+//                 {
+//                     Info<< "move points " << longName << nl;
+//                 }
+//                 vtkgeom = vtpData.getCopy();
+//                 vtkgeom->SetPoints(movePoints(mesh, vtpData));
+//             }
+//         }
+// #endif
+
+        if (!vtkgeom)
         {
-            continue;
+            vtpData.clear(); // Remove any old mappings
+
+            // Nothing usable from cache - create new geometry
+            vtkgeom = patchVTKMesh(mesh.patch());
         }
 
+        vtpData.set(vtkgeom);
+    }
+
+    if (debug)
+    {
+        Info<< "<end> " << FUNCTION_NAME << nl;
+        printMemory();
+    }
+}
+
+
+void Foam::vtkPVFoam::convertMeshCellZones()
+{
+    const List<label> partIds =
+        rangeCellZones_.intersection(selectedPartIds_);
+
+    const fvMesh& mesh = *volMeshPtr_;
+
+    if (partIds.empty())
+    {
+        return;
+    }
+
+    if (debug)
+    {
+        Info<< "<beg> " << FUNCTION_NAME << nl;
+        printMemory();
+    }
+
+    const cellZoneMesh& zMesh = mesh.cellZones();
+    for (const auto partId : partIds)
+    {
         const auto& longName = selectedPartIds_[partId];
         const word zoneName = getFoamName(longName);
         const label  zoneId = zMesh.findZoneID(zoneName);
@@ -318,7 +390,7 @@ void Foam::vtkPVFoam::convertMeshCellZones()
         if (debug)
         {
             Info<< "Creating VTK mesh for cellZone[" << zoneId << "] "
-                << zoneName << endl;
+                << zoneName << nl;
         }
 
         foamVtuData& vtuData = cachedVtu_(longName);
@@ -362,7 +434,7 @@ void Foam::vtkPVFoam::convertMeshCellZones()
 
     if (debug)
     {
-        Info<< "<end> " << FUNCTION_NAME << endl;
+        Info<< "<end> " << FUNCTION_NAME << nl;
         printMemory();
     }
 }
@@ -370,28 +442,25 @@ void Foam::vtkPVFoam::convertMeshCellZones()
 
 void Foam::vtkPVFoam::convertMeshCellSets()
 {
-    const arrayRange& range = rangeCellSets_;
-    const fvMesh& mesh = *meshPtr_;
+    const List<label> partIds =
+        rangeCellSets_.intersection(selectedPartIds_);
+
+    const fvMesh& mesh = *volMeshPtr_;
 
     if (debug)
     {
-        Info<< "<beg> " << FUNCTION_NAME << endl;
+        Info<< "<beg> " << FUNCTION_NAME << nl;
         printMemory();
     }
 
-    for (auto partId : range)
+    for (const auto partId : partIds)
     {
-        if (!selectedPartIds_.found(partId))
-        {
-            continue;
-        }
-
         const auto& longName = selectedPartIds_[partId];
         const word partName = getFoamName(longName);
 
         if (debug)
         {
-            Info<< "Creating VTK mesh for cellSet=" << partName << endl;
+            Info<< "Creating VTK mesh for cellSet=" << partName << nl;
         }
 
         foamVtuData& vtuData = cachedVtu_(longName);
@@ -435,7 +504,7 @@ void Foam::vtkPVFoam::convertMeshCellSets()
 
     if (debug)
     {
-        Info<< "<end> " << FUNCTION_NAME << endl;
+        Info<< "<end> " << FUNCTION_NAME << nl;
         printMemory();
     }
 }
@@ -443,28 +512,25 @@ void Foam::vtkPVFoam::convertMeshCellSets()
 
 void Foam::vtkPVFoam::convertMeshFaceZones()
 {
-    const arrayRange& range = rangeFaceZones_;
-    const fvMesh& mesh = *meshPtr_;
+    const List<label> partIds =
+        rangeFaceZones_.intersection(selectedPartIds_);
+
+    const fvMesh& mesh = *volMeshPtr_;
 
-    if (range.empty())
+    if (partIds.empty())
     {
         return;
     }
 
     if (debug)
     {
-        Info<< "<beg> " << FUNCTION_NAME << endl;
+        Info<< "<beg> " << FUNCTION_NAME << nl;
         printMemory();
     }
 
     const faceZoneMesh& zMesh = mesh.faceZones();
-    for (auto partId : range)
+    for (const auto partId : partIds)
     {
-        if (!selectedPartIds_.found(partId))
-        {
-            continue;
-        }
-
         const auto& longName = selectedPartIds_[partId];
         const word zoneName = getFoamName(longName);
         const label  zoneId = zMesh.findZoneID(zoneName);
@@ -476,7 +542,7 @@ void Foam::vtkPVFoam::convertMeshFaceZones()
         if (debug)
         {
             Info<< "Creating VTKmesh for faceZone[" << zoneId << "] "
-                << zoneName << endl;
+                << zoneName << nl;
         }
 
         foamVtpData& vtpData = cachedVtp_(longName);
@@ -513,7 +579,7 @@ void Foam::vtkPVFoam::convertMeshFaceZones()
 
     if (debug)
     {
-        Info<< "<end> " << FUNCTION_NAME << endl;
+        Info<< "<end> " << FUNCTION_NAME << nl;
         printMemory();
     }
 }
@@ -521,28 +587,25 @@ void Foam::vtkPVFoam::convertMeshFaceZones()
 
 void Foam::vtkPVFoam::convertMeshFaceSets()
 {
-    const arrayRange& range = rangeFaceSets_;
-    const fvMesh& mesh = *meshPtr_;
+    const List<label> partIds =
+        rangeFaceSets_.intersection(selectedPartIds_);
+
+    const fvMesh& mesh = *volMeshPtr_;
 
     if (debug)
     {
-        Info<< "<beg> " << FUNCTION_NAME << endl;
+        Info<< "<beg> " << FUNCTION_NAME << nl;
         printMemory();
     }
 
-    for (auto partId : range)
+    for (const auto partId : partIds)
     {
-        if (!selectedPartIds_.found(partId))
-        {
-            continue;
-        }
-
         const auto& longName = selectedPartIds_[partId];
         const word partName = getFoamName(longName);
 
         if (debug)
         {
-            Info<< "Creating VTK mesh for faceSet=" << partName << endl;
+            Info<< "Creating VTK mesh for faceSet=" << partName << nl;
         }
 
         foamVtpData& vtpData = cachedVtp_(longName);
@@ -597,7 +660,7 @@ void Foam::vtkPVFoam::convertMeshFaceSets()
 
     if (debug)
     {
-        Info<< "<end> " << FUNCTION_NAME << endl;
+        Info<< "<end> " << FUNCTION_NAME << nl;
         printMemory();
     }
 }
@@ -605,84 +668,83 @@ void Foam::vtkPVFoam::convertMeshFaceSets()
 
 void Foam::vtkPVFoam::convertMeshPointZones()
 {
-    const arrayRange& range = rangePointZones_;
-    const fvMesh& mesh = *meshPtr_;
+    const List<label> partIds =
+        rangePointZones_.intersection(selectedPartIds_);
+
+    const fvMesh& mesh = *volMeshPtr_;
+
+    if (partIds.empty())
+    {
+        return;
+    }
 
     if (debug)
     {
-        Info<< "<beg> " << FUNCTION_NAME << endl;
+        Info<< "<beg> " << FUNCTION_NAME << nl;
         printMemory();
     }
 
-    if (range.size())
+    const pointZoneMesh& zMesh = mesh.pointZones();
+    for (const auto partId : partIds)
     {
-        const pointZoneMesh& zMesh = mesh.pointZones();
-        for (auto partId : range)
+        const auto& longName = selectedPartIds_[partId];
+        const word zoneName = getFoamName(longName);
+        const label zoneId = zMesh.findZoneID(zoneName);
+
+        if (zoneId < 0)
         {
-            if (!selectedPartIds_.found(partId))
-            {
-                continue;
-            }
+            continue;
+        }
 
-            const auto& longName = selectedPartIds_[partId];
-            const word zoneName = getFoamName(longName);
-            const label zoneId = zMesh.findZoneID(zoneName);
+        foamVtpData& vtpData = cachedVtp_(longName);
 
-            if (zoneId < 0)
+        vtkSmartPointer<vtkPolyData> vtkgeom;
+        if (vtpData.nPoints() && vtpData.pointMap().size())
+        {
+            if (meshState_ == polyMesh::UNCHANGED)
             {
+                if (debug)
+                {
+                    Info<< "reusing " << longName << nl;
+                }
+                vtpData.reuse();
                 continue;
             }
-
-            foamVtpData& vtpData = cachedVtp_(longName);
-
-            vtkSmartPointer<vtkPolyData> vtkgeom;
-            if (vtpData.nPoints() && vtpData.pointMap().size())
+            else if (meshState_ == polyMesh::POINTS_MOVED)
             {
-                if (meshState_ == polyMesh::UNCHANGED)
-                {
-                    if (debug)
-                    {
-                        Info<< "reusing " << longName << nl;
-                    }
-                    vtpData.reuse();
-                    continue;
-                }
-                else if (meshState_ == polyMesh::POINTS_MOVED)
+                if (debug)
                 {
-                    if (debug)
-                    {
-                        Info<< "move points " << longName << nl;
-                    }
-                    vtkgeom = vtpData.getCopy();
+                    Info<< "move points " << longName << nl;
                 }
+                vtkgeom = vtpData.getCopy();
             }
+        }
 
-            if (!vtkgeom)
-            {
-                // First time, or topo change
-                vtkgeom = vtkSmartPointer<vtkPolyData>::New();
-                vtpData.pointMap() = zMesh[zoneId];
-            }
-
-            const pointField& points = mesh.points();
-            const labelUList& pointMap = vtpData.pointMap();
+        if (!vtkgeom)
+        {
+            // First time, or topo change
+            vtkgeom = vtkSmartPointer<vtkPolyData>::New();
+            vtpData.pointMap() = zMesh[zoneId];
+        }
 
-            auto vtkpoints = vtkSmartPointer<vtkPoints>::New();
+        const pointField& points = mesh.points();
+        const labelUList& pointMap = vtpData.pointMap();
 
-            vtkpoints->SetNumberOfPoints(pointMap.size());
-            forAll(pointMap, pointi)
-            {
-                vtkpoints->SetPoint(pointi, points[pointMap[pointi]].v_);
-            }
+        auto vtkpoints = vtkSmartPointer<vtkPoints>::New();
 
-            vtkgeom->SetPoints(vtkpoints);
-            vtpData.set(vtkgeom);
+        vtkpoints->SetNumberOfPoints(pointMap.size());
+        forAll(pointMap, pointi)
+        {
+            vtkpoints->SetPoint(pointi, points[pointMap[pointi]].v_);
         }
+
+        vtkgeom->SetPoints(vtkpoints);
+        vtpData.set(vtkgeom);
     }
 
     if (debug)
     {
-        Info<< "<end> " << FUNCTION_NAME << endl;
+        Info<< "<end> " << FUNCTION_NAME << nl;
         printMemory();
     }
 }
@@ -690,22 +752,19 @@ void Foam::vtkPVFoam::convertMeshPointZones()
 
 void Foam::vtkPVFoam::convertMeshPointSets()
 {
-    const arrayRange& range = rangePointSets_;
-    const fvMesh& mesh = *meshPtr_;
+    const List<label> partIds =
+        rangePointSets_.intersection(selectedPartIds_);
+
+    const fvMesh& mesh = *volMeshPtr_;
 
     if (debug)
     {
-        Info<< "<beg> " << FUNCTION_NAME << endl;
+        Info<< "<beg> " << FUNCTION_NAME << nl;
         printMemory();
     }
 
-    for (auto partId : range)
+    for (const auto partId : partIds)
     {
-        if (!selectedPartIds_.found(partId))
-        {
-            continue;
-        }
-
         const auto& longName = selectedPartIds_[partId];
         const word partName = getFoamName(longName);
 
@@ -757,7 +816,7 @@ void Foam::vtkPVFoam::convertMeshPointSets()
 
     if (debug)
     {
-        Info<< "<end> " << FUNCTION_NAME << endl;
+        Info<< "<end> " << FUNCTION_NAME << nl;
         printMemory();
     }
 }
diff --git a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamMeshLagrangian.C b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamMeshLagrangian.C
index f7be26d4be3f389e481facc0f58ed8d4759a804d..d7706d66565ebad567329167703198e10c9248c6 100644
--- a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamMeshLagrangian.C
+++ b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamMeshLagrangian.C
@@ -50,7 +50,7 @@ vtkSmartPointer<vtkPolyData> Foam::vtkPVFoam::lagrangianVTKMesh
     if (debug)
     {
         Info<< "<beg> lagrangianVTKMesh - timePath "
-            << mesh.time().timePath()/cloud::prefix/cloudName << endl;
+            << mesh.time().timePath()/cloud::prefix/cloudName << nl;
         printMemory();
     }
 
@@ -71,7 +71,7 @@ vtkSmartPointer<vtkPolyData> Foam::vtkPVFoam::lagrangianVTKMesh
 
         if (debug)
         {
-            Info<< "cloud with " << parcels.size() << " parcels" << endl;
+            Info<< "cloud with " << parcels.size() << " parcels" << nl;
         }
 
         auto vtkpoints = vtkSmartPointer<vtkPoints>::New();
@@ -91,7 +91,7 @@ vtkSmartPointer<vtkPolyData> Foam::vtkPVFoam::lagrangianVTKMesh
 
     if (debug)
     {
-        Info<< "<end> lagrangianVTKMesh" << endl;
+        Info<< "<end> lagrangianVTKMesh" << nl;
         printMemory();
     }
 
diff --git a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamMeshVolume.C b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamMeshVolume.C
index 18beae4583e23308ffccc604ece0a8a66191f0d5..f1d1ebd36b60e7b676000c3fb90ebb64dc951c72 100644
--- a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamMeshVolume.C
+++ b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamMeshVolume.C
@@ -115,7 +115,7 @@ vtkSmartPointer<vtkUnstructuredGrid> Foam::vtkPVFoam::volumeVTKMesh
 {
     if (debug)
     {
-        Info<< "<beg> " << FUNCTION_NAME << endl;
+        Info<< "<beg> " << FUNCTION_NAME << nl;
         printMemory();
     }
 
@@ -209,7 +209,7 @@ vtkSmartPointer<vtkUnstructuredGrid> Foam::vtkPVFoam::volumeVTKMesh
 
     if (debug)
     {
-        Info<< "<end> " << FUNCTION_NAME << endl;
+        Info<< "<end> " << FUNCTION_NAME << nl;
         printMemory();
     }
 
diff --git a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamTemplates.C b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamTemplates.C
index ec64378d2f37b454211a0f009b086e4329a92741..986477320d86f0b1e04b8e303c66a656cc338968 100644
--- a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamTemplates.C
+++ b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamTemplates.C
@@ -69,9 +69,9 @@ vtkSmartPointer<vtkCellArray> Foam::vtkPVFoam::patchFacesVTKCells
     const faceList& faces = p.localFaces();
 
     label nAlloc = faces.size();
-    forAll(faces, facei)
+    for (const face& f : faces)
     {
-        nAlloc += faces[facei].size();
+        nAlloc += f.size();
     }
 
     auto cells = vtkSmartPointer<vtkCellArray>::New();
@@ -87,15 +87,13 @@ vtkSmartPointer<vtkCellArray> Foam::vtkPVFoam::patchFacesVTKCells
     // Cell connectivity for polygons
     // [size, verts..., size, verts... ]
     label idx = 0;
-    forAll(faces, facei)
+    for (const face& f : faces)
     {
-        const face& f = faces[facei];
-
         cellsUL[idx++] = f.size();
 
-        forAll(f, fp)
+        for (const label verti : f)
         {
-            cellsUL[idx++] = f[fp];
+            cellsUL[idx++] = verti;
         }
     }
 
diff --git a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamUpdateInfo.C b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamUpdateInfo.C
index bae3336c120f5620798235d8f419f550712f7749..290f4ad41be85d59e820f8aa9e24f21dd438da52 100644
--- a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamUpdateInfo.C
+++ b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamUpdateInfo.C
@@ -34,11 +34,17 @@ License
 #include "polyBoundaryMeshEntries.H"
 #include "entry.H"
 #include "Cloud.H"
-#include "vtkPVFoamReader.h"
+#include "areaFaMesh.H"
 
 // VTK includes
 #include "vtkDataArraySelection.h"
 
+// OpenFOAM/VTK interface
+#include "vtkPVFoamReader.h"
+
+// Templates (only needed here)
+#include "vtkPVFoamUpdateTemplates.C"
+
 // * * * * * * * * * * * * * * * Private Classes * * * * * * * * * * * * * * //
 
 namespace Foam
@@ -142,17 +148,57 @@ void Foam::vtkPVFoam::updateInfoInternalMesh
 {
     if (debug)
     {
-        Info<< "<beg> updateInfoInternalMesh" << endl;
+        Info<< "<beg> " << FUNCTION_NAME << nl;
     }
 
     // Determine mesh parts (internalMesh, patches...)
-    //- Add internal mesh as first entry
+    // Add internal mesh as first entry
     rangeVolume_.reset(select->GetNumberOfArrays(), 1);
     select->AddArray("internalMesh");
 
     if (debug)
     {
-        Info<< "<end> updateInfoInternalMesh" << endl;
+        Info<< "<end> " << FUNCTION_NAME << nl;
+    }
+}
+
+
+void Foam::vtkPVFoam::updateInfoAreaMesh
+(
+    vtkDataArraySelection* select
+)
+{
+    if (debug)
+    {
+        Info<< "<beg> " << FUNCTION_NAME << nl;
+    }
+
+    rangeArea_.reset(select->GetNumberOfArrays(), 0);
+
+    // Use the db directly since this might be called without a mesh,
+    // but the region must get added back in
+    fileName faMeshPrefix(faMesh::meshSubDir);
+    if (meshRegion_ != polyMesh::defaultRegion)
+    {
+        faMeshPrefix = meshRegion_/faMeshPrefix;
+    }
+
+    // Check for finiteArea mesh
+    if
+    (
+        isFile
+        (
+            dbPtr_->constantPath()/faMeshPrefix/"faceLabels"
+        )
+    )
+    {
+        rangeArea_ += 1;
+        select->AddArray("areaMesh");
+    }
+
+    if (debug)
+    {
+        Info<< "<end> " << FUNCTION_NAME << nl;
     }
 }
 
@@ -165,7 +211,7 @@ void Foam::vtkPVFoam::updateInfoLagrangian
     if (debug)
     {
         Info<< "<beg> " << FUNCTION_NAME << nl
-            << "    " << dbPtr_->timePath()/cloud::prefix << endl;
+            << "    " << dbPtr_->timePath()/cloud::prefix << nl;
     }
 
     // Use the db directly since this might be called without a mesh,
@@ -191,12 +237,12 @@ void Foam::vtkPVFoam::updateInfoLagrangian
         );
     }
 
-    rangeLagrangian_.reset(select->GetNumberOfArrays());
-    rangeLagrangian_ += addToArray(select, "lagrangian/", names.sortedToc());
+    rangeClouds_.reset(select->GetNumberOfArrays());
+    rangeClouds_ += addToArray(select, "lagrangian/", names.sortedToc());
 
     if (debug)
     {
-        Info<< "<end> " << FUNCTION_NAME << endl;
+        Info<< "<end> " << FUNCTION_NAME << nl;
     }
 }
 
@@ -210,14 +256,14 @@ void Foam::vtkPVFoam::updateInfoPatches
     if (debug)
     {
         Info<< "<beg> " << FUNCTION_NAME
-            << " [meshPtr=" << (meshPtr_ ? "set" : "null") << "]" << endl;
+            << " [volMeshPtr=" << (volMeshPtr_ ? "set" : "null") << "]" << nl;
     }
 
     rangePatches_.reset(select->GetNumberOfArrays());
 
-    if (meshPtr_)
+    if (volMeshPtr_)
     {
-        const polyBoundaryMesh& patches = meshPtr_->boundaryMesh();
+        const polyBoundaryMesh& patches = volMeshPtr_->boundaryMesh();
         const HashTable<labelList>& groups = patches.groupPatchIDs();
         DynamicList<string> displayNames(groups.size());
 
@@ -387,7 +433,7 @@ void Foam::vtkPVFoam::updateInfoPatches
 
     if (debug)
     {
-        Info<< "<end> " << FUNCTION_NAME << endl;
+        Info<< "<end> " << FUNCTION_NAME << nl;
     }
 }
 
@@ -404,17 +450,16 @@ void Foam::vtkPVFoam::updateInfoZones
 
     if (debug)
     {
-        Info<< "<beg> updateInfoZones"
-            << " [meshPtr=" << (meshPtr_ ? "set" : "null") << "]" << endl;
+        Info<< "<beg> " << FUNCTION_NAME
+            << " [volMeshPtr=" << (volMeshPtr_ ? "set" : "null") << "]" << nl;
     }
 
-
     // cellZones
     {
         const wordList names =
         (
-            meshPtr_
-          ? getZoneNames(meshPtr_->cellZones())
+            volMeshPtr_
+          ? getZoneNames(volMeshPtr_->cellZones())
           : getZoneNames("cellZones")
         );
 
@@ -426,8 +471,8 @@ void Foam::vtkPVFoam::updateInfoZones
     {
         const wordList names =
         (
-            meshPtr_
-          ? getZoneNames(meshPtr_->faceZones())
+            volMeshPtr_
+          ? getZoneNames(volMeshPtr_->faceZones())
           : getZoneNames("faceZones")
         );
 
@@ -439,8 +484,8 @@ void Foam::vtkPVFoam::updateInfoZones
     {
         const wordList names =
         (
-            meshPtr_
-          ? getZoneNames(meshPtr_->pointZones())
+            volMeshPtr_
+          ? getZoneNames(volMeshPtr_->pointZones())
           : getZoneNames("pointZones")
         );
 
@@ -450,7 +495,7 @@ void Foam::vtkPVFoam::updateInfoZones
 
     if (debug)
     {
-        Info<< "<end> updateInfoZones" << endl;
+        Info<< "<end> " << FUNCTION_NAME << nl;
     }
 }
 
@@ -467,7 +512,7 @@ void Foam::vtkPVFoam::updateInfoSets
 
     if (debug)
     {
-        Info<< "<beg> updateInfoSets" << endl;
+        Info<< "<beg> " << FUNCTION_NAME << nl;
     }
 
     // Add names of sets. Search for last time directory with a sets
@@ -493,7 +538,7 @@ void Foam::vtkPVFoam::updateInfoSets
     if (debug)
     {
         Info<< "     updateInfoSets read "
-            << objects.names() << " from " << setsInstance << endl;
+            << objects.names() << " from " << setsInstance << nl;
     }
 
 
@@ -523,7 +568,91 @@ void Foam::vtkPVFoam::updateInfoSets
 
     if (debug)
     {
-        Info<< "<end> updateInfoSets" << endl;
+        Info<< "<end> " << FUNCTION_NAME << nl;
+    }
+}
+
+
+void Foam::vtkPVFoam::updateInfoContinuumFields
+(
+    vtkDataArraySelection* select
+)
+{
+    if (debug)
+    {
+        Info<< "<beg> " << FUNCTION_NAME << nl;
+    }
+
+    // Preserve the enabled selections
+    HashSet<string> enabled;
+    if (!select->GetNumberOfArrays() && !volMeshPtr_)
+    {
+        // First call: automatically enable 'p' and 'U'
+        enabled = { "p", "U" };
+    }
+    else
+    {
+        enabled = getSelectedArraySet(select);
+    }
+
+    select->RemoveAllArrays();   // Clear existing list
+
+    // Use the db directly since this might be called without a mesh,
+    // but the region name is still required
+    word regionPrefix;
+    if (meshRegion_ != polyMesh::defaultRegion)
+    {
+        regionPrefix = meshRegion_;
+    }
+
+    // Objects for this time and mesh region
+    IOobjectList objects(dbPtr_(), dbPtr_().timeName(), regionPrefix);
+
+    updateInfoFields<fvPatchField, volMesh>(select, objects);
+    updateInfoFields<faPatchField, areaMesh>(select, objects);
+
+    setSelectedArrayEntries(select, enabled);  // Adjust/restore selected
+
+    if (debug)
+    {
+        Info<< "<end> " << FUNCTION_NAME << nl;
+    }
+}
+
+
+void Foam::vtkPVFoam::updateInfoPointFields
+(
+    vtkDataArraySelection* select
+)
+{
+    if (debug)
+    {
+        Info<< "<beg> " << FUNCTION_NAME << nl;
+    }
+
+    // Preserve the enabled selections
+    HashSet<string> enabled = getSelectedArraySet(select);
+
+    select->RemoveAllArrays();   // Clear existing list
+
+    // Use the db directly since this might be called without a mesh,
+    // but the region name is still required
+    word regionPrefix;
+    if (meshRegion_ != polyMesh::defaultRegion)
+    {
+        regionPrefix = meshRegion_;
+    }
+
+    // Objects for this time and mesh region
+    IOobjectList objects(dbPtr_(), dbPtr_().timeName(), regionPrefix);
+
+    updateInfoFields<pointPatchField, pointMesh>(select, objects);
+
+    setSelectedArrayEntries(select, enabled);  // Adjust/restore selected
+
+    if (debug)
+    {
+        Info<< "<end> " << FUNCTION_NAME << nl;
     }
 }
 
@@ -535,14 +664,14 @@ void Foam::vtkPVFoam::updateInfoLagrangianFields
 {
     if (debug)
     {
-        Info<< "<beg> updateInfoLagrangianFields" << endl;
+        Info<< "<beg> " << FUNCTION_NAME << nl;
     }
 
     // Preserve the enabled selections
     HashSet<string> enabled = getSelectedArraySet(select);
     select->RemoveAllArrays();
 
-    const arrayRange& range = rangeLagrangian_;
+    const arrayRange& range = rangeClouds_;
     if (range.empty())
     {
         return;
@@ -600,7 +729,7 @@ void Foam::vtkPVFoam::updateInfoLagrangianFields
 
     if (debug)
     {
-        Info<< "<end> " << FUNCTION_NAME << endl;
+        Info<< "<end> " << FUNCTION_NAME << nl;
     }
 }
 
diff --git a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamUpdateTemplates.C b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamUpdateTemplates.C
index 4f11fdb1d6fb5caa5d75554d7efdcc420626366a..796d39f3771f5108786a1ea90c7ff2da5620e0be 100644
--- a/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamUpdateTemplates.C
+++ b/applications/utilities/postProcessing/graphics/PVReaders/vtkPVFoam/vtkPVFoamUpdateTemplates.C
@@ -34,43 +34,19 @@ InClass
 template<template<class> class patchType, class meshType>
 void Foam::vtkPVFoam::updateInfoFields
 (
-    vtkDataArraySelection* select
+    vtkDataArraySelection* select,
+    const IOobjectList& objects
 )
 {
     if (debug)
     {
         Info<< "<beg> updateInfoFields <"
             << meshType::Mesh::typeName
-            << "> [meshPtr=" << (meshPtr_ ? "set" : "null") << "]"
-            << endl;
+            << "> [volMeshPtr=" << (volMeshPtr_ ? "set" : "null") << "]"
+            << nl;
     }
 
-    HashSet<string> enabled;
-    if (!select->GetNumberOfArrays() && !meshPtr_)
-    {
-        // Enable 'p' and 'U' only on the first call
-        enabled = { "p", "U" };
-    }
-    else
-    {
-        // Preserve the enabled selections
-        enabled = getSelectedArraySet(select);
-    }
-
-    select->RemoveAllArrays();
-
-    // use the db directly since this might be called without a mesh,
-    // but the region must get added back in
-    word regionPrefix;
-    if (meshRegion_ != polyMesh::defaultRegion)
-    {
-        regionPrefix = meshRegion_;
-    }
-
-    // Search for list of objects for this time and mesh region
-    IOobjectList objects(dbPtr_(), dbPtr_().timeName(), regionPrefix);
-
-    // Add volume fields to GUI
+    // Add geometric fields (volume/area) to GUI
     addToSelection<GeometricField<scalar, patchType, meshType>>
     (
         select,
@@ -91,14 +67,13 @@ void Foam::vtkPVFoam::updateInfoFields
         select,
         objects
     );
-
     addToSelection<GeometricField<tensor, patchType, meshType>>
     (
         select,
         objects
     );
 
-    // Add dimensioned fields to GUI
+    // Add dimensioned fields (volume/area) to GUI
     addToSelection<DimensionedField<scalar, meshType>>
     (
         select,
@@ -125,16 +100,13 @@ void Foam::vtkPVFoam::updateInfoFields
         objects
     );
 
-
-    // Restore the enabled selections
-    setSelectedArrayEntries(select, enabled);
-
     if (debug)
     {
-        Info<< "<end> updateInfoFields" << endl;
+        Info<< "<end> updateInfoFields" << nl;
     }
 }
 
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 #endif