cutPoly: New polyhedral cutting routines and isoSurface algorithm

A set of routines for cutting polyhedra have been added. These can cut
polyhedral cells based on the adjacent point values and an iso-value
which defines the surface. The method operates directly on the
polyhedral cells; it does not decompose them into tetrahedra at any
point. The routines can compute the cut topology as well as integrals of
properties above and below the cut surface.

An iso-surface algorithm has been added based on these polyhedral
cutting routines. It is significantly more robust than the previous
algorithm, and produces compact surfaces equivalent to the previous
algorithm's maximum filtering level. It is also approximately 3 times
faster than the previous algorithm, and 10 times faster when run
repeatedly on the same set of cells (this is because some addressing is
cached and reused).

This algorithm is used by the 'isoSurface', 'distanceSurface' and
'cutPlane' sampled surfaces.

The 'cutPlane' sampled surface is a renaming of 'cuttingPlane' to make
it consistent with the corresponding packaged function. The name
'cuttingPlane' has been retained for backwards compatibility and can
still be used to select a 'cutPlane' surface. The legacy 'plane' surface
has been removed.

The 'average' keyword has been removed from specification of these
sampled surfaces as cell-centred values are no longer used in the
generation of or interpolation to an iso-surface. The 'filtering'
keyword has also been removed as it relates to options within the
previous algorithm. Zone support has been reinstated into the
'isoSurface' sampled surface. Interpolation to all these sampled
surfaces has been corrected to exactly match the user-selected
interpolation scheme, and the interpolation procedure no longer
unnecessarily re-generates data that is already available.

applications/test/tetTetOverlap/Test-tetTetOverlap.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2012-2018 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -32,7 +32,7 @@ Description
 #include "tetPointRef.H"
 #include "OFstream.H"
 #include "meshTools.H"
-#include "cut.H"
+#include "cutTriTet.H"
 
 using namespace Foam;
 
@@ -87,9 +87,9 @@ int main(int argc, char *argv[])
     // Do intersection
     typedef DynamicList<FixedList<point, 4>> tetList;
     tetList tetsIn1, tetsIn2, tetsOut;
-    cut::appendOp<tetList> tetOpIn1(tetsIn1);
-    cut::appendOp<tetList> tetOpIn2(tetsIn2);
-    cut::appendOp<tetList> tetOpOut(tetsOut);
+    cutTriTet::appendOp<tetList> tetOpIn1(tetsIn1);
+    cutTriTet::appendOp<tetList> tetOpIn2(tetsIn2);
+    cutTriTet::appendOp<tetList> tetOpOut(tetsOut);
 
     const plane p0(tetB[1], tetB[3], tetB[2]);
     tetsIn1.clear();

applications/utilities/mesh/generation/foamyMesh/foamyHexMeshBackgroundMesh/foamyHexMeshBackgroundMesh.C

@@ -41,7 +41,7 @@ Description
 #include "cellShape.H"
 #include "cellModeller.H"
 #include "DynamicField.H"
-#include "isoSurface.H"
+#include "cutPolyIsoSurface.H"
 #include "vtkSurfaceWriter.H"
 #include "syncTools.H"
 
@@ -707,13 +707,7 @@ int main(int argc, char *argv[])
             }
         }
 
-        isoSurface iso
-        (
-            mesh,
-            cellDistance,
-            pointDistance,
-            0
-        );
+        cutPolyIsoSurface iso(mesh, pointDistance, 0);
 
         isoFaces.setSize(iso.size());
         forAll(isoFaces, i)

etc/caseDicts/postProcessing/surface/cutPlaneSurface.cfg

@@ -12,7 +12,7 @@ surfaces
 (
     cutPlane
     {
-        type        cuttingPlane;
+        type        cutPlane;
         interpolate $interpolate;
         planeType   pointAndNormal;
         point       $point;

src/OpenFOAM/containers/HashTables/HashList/HashList.C

@@ -0,0 +1,110 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "HashList.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    template<class Type, class Key, class Hash>
+    const Key HashList<Type, Key, Hash>::nullKey = Key::null;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class Type, class Key, class Hash>
+Foam::HashList<Type, Key, Hash>::HashList(const label size)
+:
+    List<Tuple2<Key, Type>>(size, Tuple2<Key, Type>(nullKey, Type()))
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+template<class Type, class Key, class Hash>
+bool Foam::HashList<Type, Key, Hash>::insert(const Key& k, const Type& t)
+{
+    List<Tuple2<Key, Type>>& map = *this;
+
+    const label n = map.size();
+
+    const unsigned h = Hash()(k);
+
+    for (label i = 0; i < n; i ++)
+    {
+        const label hi = (h + i) % n;
+
+        if (map[hi].first() == nullKey)
+        {
+            map[hi] = Tuple2<Key, Type>(k, t);
+            return true;
+        }
+
+        if (map[hi].first() == k)
+        {
+            return false;
+        }
+    }
+
+    FatalErrorInFunction
+        << "Hash list is full"
+        << exit(FatalError);
+
+    return false;
+}
+
+
+// * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * * //
+
+template<class Type, class Key, class Hash>
+const Type& Foam::HashList<Type, Key, Hash>::operator[](const Key& k) const
+{
+    const List<Tuple2<Key, Type>>& map = *this;
+
+    const label n = map.size();
+
+    const unsigned h = Hash()(k);
+
+    for (label i = 0; i < n; i ++)
+    {
+        const label hi = (h + i) % n;
+
+        if (map[hi].first() == k)
+        {
+            return map[hi].second();
+        }
+    }
+
+    FatalErrorInFunction
+        << "Hash list does not contain key \"" << k << "\""
+        << exit(FatalError);
+
+    return NullObjectRef<Type>();
+}
+
+
+// ************************************************************************* //

src/OpenFOAM/containers/HashTables/HashList/HashList.H

@@ -0,0 +1,104 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::HashList
+
+Description
+    HashList class. Like HashTable, but much less dynamic memory-y. Should be
+    faster for small sets of non-dynamic primitive types (labels, edges,
+    points, etc...). It is also much less functional at present. There is no
+    re-sizing, so you have to make sure it is constructed sufficiently large to
+    hold all the data that will ever be inserted into it.
+
+SourceFiles
+    HashList.C
+
+\*---------------------------------------------------------------------------*/
+
+#include "List.H"
+#include "Tuple2.H"
+#include "word.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifndef HashList_H
+#define HashList_H
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                          Class HashList Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class Type, class Key=word, class Hash=string::hash>
+class HashList
+:
+    private List<Tuple2<Key, Type>>
+{
+public:
+
+    // Public Static Member Data
+
+        //- Null key value for unset elements in the list
+        static const Key nullKey;
+
+
+    // Constructors
+
+        //- Construct given a size
+        HashList(const label size);
+
+
+    // Member Functions
+
+        //- Insert into the hash list. Return true if the value was newly
+        //  inserted, or false if it was already there.
+        bool insert(const Key& k, const Type& t);
+
+
+    // Member Operators
+
+        //- Retrieve from the hash list
+        const Type& operator[](const Key& k) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+    #include "HashList.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/OpenFOAM/meshes/primitiveShapes/plane/plane.C

@@ -309,13 +309,19 @@ Foam::point Foam::plane::aPoint() const
 
 Foam::point Foam::plane::nearestPoint(const point& p) const
 {
-    return p - normal_*((p - point_) & normal_);
+    return p - normal_*signedDistance(p);
 }
 
 
 Foam::scalar Foam::plane::distance(const point& p) const
 {
-    return mag((p - point_) & normal_);
+    return mag(signedDistance(p));
+}
+
+
+Foam::scalar Foam::plane::signedDistance(const point& p) const
+{
+    return (p - point_) & normal_;
 }
 
 

src/OpenFOAM/meshes/primitiveShapes/plane/plane.H

@@ -194,6 +194,9 @@ public:
         //- Return distance from the given point to the plane
         scalar distance(const point& p) const;
 
+        //- Return signed distance from the given point to the plane
+        scalar signedDistance(const point& p) const;
+
         //- Return cut coefficient for plane and line defined by
         //  origin and direction
         scalar normalIntersect(const point& pnt0, const vector& dir) const;

src/finiteVolume/Make/files

@@ -290,6 +290,7 @@ $(interpolation)/interpolation/interpolations.C
 
 $(interpolation)/interpolationCell/makeInterpolationCell.C
 $(interpolation)/interpolationCellPatchConstrained/makeInterpolationCellPatchConstrained.C
+$(interpolation)/interpolationVolPointInterpolation/makeInterpolationVolPointInterpolation.C
 $(interpolation)/interpolationCellPoint/cellPointWeight/cellPointWeight.C
 $(interpolation)/interpolationCellPoint/makeInterpolationCellPoint.C
 $(interpolation)/interpolationCellPointFace/makeInterpolationCellPointFace.C

src/finiteVolume/cfdTools/general/levelSet/levelSet.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2017-2019 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2017-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -24,7 +24,7 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "levelSet.H"
-#include "cut.H"
+#include "cutTriTet.H"
 #include "polyMeshTetDecomposition.H"
 #include "tetIndices.H"
 
@@ -38,6 +38,9 @@ Foam::tmp<Foam::scalarField> Foam::levelSetFraction
     const bool above
 )
 {
+    typedef cutTriTet::noOp noOp;
+    typedef cutTriTet::volumeOp volumeOp;
+
     tmp<scalarField> tResult(new scalarField(mesh.nCells(), Zero));
     scalarField& result = tResult.ref();
 
@@ -69,15 +72,15 @@ Foam::tmp<Foam::scalarField> Foam::levelSetFraction
                     levelP[triIs[2]]
                 };
 
-            v += cut::volumeOp()(tet);
+            v += volumeOp()(tet);
 
             if (above)
             {
-                r += tetCut(tet, level, cut::volumeOp(), cut::noOp());
+                r += tetCut(tet, level, volumeOp(), noOp());
             }
             else
             {
-                r += tetCut(tet, level, cut::noOp(), cut::volumeOp());
+                r += tetCut(tet, level, noOp(), volumeOp());
             }
         }
 
@@ -96,6 +99,9 @@ Foam::tmp<Foam::scalarField> Foam::levelSetFraction
     const bool above
 )
 {
+    typedef cutTriTet::noOp noOp;
+    typedef cutTriTet::areaMagOp areaMagOp;
+
     tmp<scalarField> tResult(new scalarField(patch.size(), 0));
     scalarField& result = tResult.ref();
 
@@ -124,15 +130,15 @@ Foam::tmp<Foam::scalarField> Foam::levelSetFraction
                     levelP[e[1]]
                 };
 
-            a += cut::areaMagOp()(tri);
+            a += areaMagOp()(tri);
 
             if (above)
             {
-                r += triCut(tri, level, cut::areaMagOp(), cut::noOp());
+                r += triCut(tri, level, areaMagOp(), noOp());
             }
             else
             {
-                r += triCut(tri, level, cut::noOp(), cut::areaMagOp());
+                r += triCut(tri, level, noOp(), areaMagOp());
             }
         }
 

src/finiteVolume/cfdTools/general/levelSet/levelSetTemplates.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2017-2019 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2017-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -24,7 +24,7 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "levelSet.H"
-#include "cut.H"
+#include "cutTriTet.H"
 #include "polyMeshTetDecomposition.H"
 #include "tetIndices.H"
 
@@ -73,7 +73,7 @@ Foam::tmp<Foam::Field<Type>> Foam::levelSetAverage
                     levelP[triIs[1]],
                     levelP[triIs[2]]
                 };
-            const cut::volumeIntegrateOp<Type>
+            const cutTriTet::volumeIntegrateOp<Type>
                 positive = FixedList<Type, 4>
                 ({
                     positiveC[cI],
@@ -81,7 +81,7 @@ Foam::tmp<Foam::Field<Type>> Foam::levelSetAverage
                     positiveP[triIs[1]],
                     positiveP[triIs[2]]
                 });
-            const cut::volumeIntegrateOp<Type>
+            const cutTriTet::volumeIntegrateOp<Type>
                 negative = FixedList<Type, 4>
                 ({
                     negativeC[cI],
@@ -90,7 +90,7 @@ Foam::tmp<Foam::Field<Type>> Foam::levelSetAverage
                     negativeP[triIs[2]]
                 });
 
-            v += cut::volumeOp()(tet);
+            v += cutTriTet::volumeOp()(tet);
 
             r += tetCut(tet, level, positive, negative);
         }
@@ -142,14 +142,14 @@ Foam::tmp<Foam::Field<Type>> Foam::levelSetAverage
                     levelP[e[0]],
                     levelP[e[1]]
                 };
-            const cut::areaMagIntegrateOp<Type>
+            const cutTriTet::areaMagIntegrateOp<Type>
                 positive = FixedList<Type, 3>
                 ({
                     positiveF[fI],
                     positiveP[e[0]],
                     positiveP[e[1]]
                 });
-            const cut::areaMagIntegrateOp<Type>
+            const cutTriTet::areaMagIntegrateOp<Type>
                 negative = FixedList<Type, 3>
                 ({
                     negativeF[fI],
@@ -157,7 +157,7 @@ Foam::tmp<Foam::Field<Type>> Foam::levelSetAverage
                     negativeP[e[1]]
                 });
 
-            a += cut::areaMagOp()(tri);
+            a += cutTriTet::areaMagOp()(tri);
 
             r += triCut(tri, level, positive, negative);
         }

src/finiteVolume/interpolation/interpolation/interpolation/interpolation.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2020 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -103,8 +103,8 @@ Foam::tmp<Foam::Field<Type>>
 Foam::fieldInterpolation<Type, InterpolationType>::interpolate
 (
     const vectorField& position,
-    const labelField& celli,
-    const labelField& facei
+    const labelList& celli,
+    const labelList& facei
 ) const
 {
     tmp<Field<Type>> tField(new Field<Type>(position.size()));
@@ -131,10 +131,10 @@ Foam::tmp<Foam::Field<Type>>
 Foam::fieldInterpolation<Type, InterpolationType>::interpolate
 (
     const Field<barycentric>& coordinates,
-    const labelField& celli,
-    const labelField& tetFacei,
-    const labelField& tetPti,
-    const labelField& facei
+    const labelList& celli,
+    const labelList& tetFacei,
+    const labelList& tetPti,
+    const labelList& facei
 ) const
 {
     tmp<Field<Type>> tField(new Field<Type>(coordinates.size()));

src/finiteVolume/interpolation/interpolation/interpolation/interpolation.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2020 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -136,8 +136,8 @@ public:
         virtual tmp<Field<Type>> interpolate
         (
             const vectorField& position,
-            const labelField& celli,
-            const labelField& facei = NullObjectRef<labelField>()
+            const labelList& celli,
+            const labelList& facei = NullObjectRef<labelList>()
         ) const = 0;
 
         //- Interpolate field to the given coordinates in the tetrahedron
@@ -155,10 +155,10 @@ public:
         virtual tmp<Field<Type>> interpolate
         (
             const Field<barycentric>& coordinates,
-            const labelField& celli,
-            const labelField& tetFacei,
-            const labelField& tetPti,
-            const labelField& facei = NullObjectRef<labelField>()
+            const labelList& celli,
+            const labelList& tetFacei,
+            const labelList& tetPti,
+            const labelList& facei = NullObjectRef<labelList>()
         ) const = 0;
 };
 
@@ -186,18 +186,18 @@ public:
         virtual tmp<Field<Type>> interpolate
         (
             const vectorField& position,
-            const labelField& celli,
-            const labelField& facei = NullObjectRef<labelField>()
+            const labelList& celli,
+            const labelList& facei = NullObjectRef<labelList>()
         ) const;
 
         //- Interpolate field to the given coordinates in the given tetrahedra
         virtual tmp<Field<Type>> interpolate
         (
             const Field<barycentric>& coordinates,
-            const labelField& celli,
-            const labelField& tetFacei,
-            const labelField& tetPti,
-            const labelField& facei = NullObjectRef<labelField>()
+            const labelList& celli,
+            const labelList& tetFacei,
+            const labelList& tetPti,
+            const labelList& facei = NullObjectRef<labelList>()
         ) const;
 };
 
@@ -208,22 +208,13 @@ public:
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-#define makeInterpolationType(SS, Type)                                        \
-                                                                               \
-defineNamedTemplateTypeNameAndDebug(SS<Type>, 0);                              \
-                                                                               \
-interpolation<Type>::adddictionaryConstructorToTable<SS<Type>>                 \
-    add##SS##Type##ConstructorToTable_;
-
-
-#define makeInterpolation(SS)                                                  \
-                                                                               \
-makeInterpolationType(SS, scalar)                                              \
-makeInterpolationType(SS, vector)                                              \
-makeInterpolationType(SS, sphericalTensor)                                     \
-makeInterpolationType(SS, symmTensor)                                          \
-makeInterpolationType(SS, tensor)
+#define defineInterpolation(Type, InterpolationType)                           \
+    defineNamedTemplateTypeNameAndDebug(InterpolationType<Type>, 0);
 
+#define makeInterpolation(Type, InterpolationType)                             \
+    interpolation<Type>::                                                      \
+        adddictionaryConstructorToTable<InterpolationType<Type>>               \
+        add##InterpolationType##Type##ConstructorToTable_;
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 

src/finiteVolume/interpolation/interpolation/interpolationCell/makeInterpolationCell.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -29,7 +29,8 @@ License
 
 namespace Foam
 {
-    makeInterpolation(interpolationCell);
+    FOR_ALL_FIELD_TYPES(defineInterpolation, interpolationCell);
+    FOR_ALL_FIELD_TYPES(makeInterpolation, interpolationCell);
 }
 
 // ************************************************************************* //

src/finiteVolume/interpolation/interpolation/interpolationCellPatchConstrained/makeInterpolationCellPatchConstrained.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -29,7 +29,16 @@ License
 
 namespace Foam
 {
-    makeInterpolation(interpolationCellPatchConstrained);
+    FOR_ALL_FIELD_TYPES
+    (
+        defineInterpolation,
+        interpolationCellPatchConstrained
+    );
+    FOR_ALL_FIELD_TYPES
+    (
+        makeInterpolation,
+        interpolationCellPatchConstrained
+    );
 }
 
 // ************************************************************************* //

src/finiteVolume/interpolation/interpolation/interpolationCellPoint/interpolationCellPoint.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2020 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -35,15 +35,7 @@ Foam::interpolationCellPoint<Type>::interpolationCellPoint
 )
 :
     fieldInterpolation<Type, interpolationCellPoint<Type>>(psi),
-    psip_
-    (
-        volPointInterpolation::New(psi.mesh()).interpolate
-        (
-            psi,
-            "volPointInterpolate(" + psi.name() + ')',
-            true        // use cache
-        )
-    )
+    interpolationVolPointInterpolation<Type>(psi)
 {
     // Uses cellPointWeight to do interpolation which needs tet decomposition
     (void)psi.mesh().tetBasePtIs();
@@ -58,7 +50,7 @@ Foam::interpolationCellPoint<Type>::interpolationCellPoint
 )
 :
     fieldInterpolation<Type, interpolationCellPoint<Type>>(psi),
-    psip_(psip)
+    interpolationVolPointInterpolation<Type>(psi, psip)
 {
     // Uses cellPointWeight to do interpolation which needs tet decomposition
     (void)psi.mesh().tetBasePtIs();

src/finiteVolume/interpolation/interpolation/interpolationCellPoint/interpolationCellPoint.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2020 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -33,7 +33,7 @@ Description
 #ifndef interpolationCellPoint_H
 #define interpolationCellPoint_H
 
-#include "interpolation.H"
+#include "interpolationVolPointInterpolation.H"
 #include "cellPointWeight.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -48,16 +48,9 @@ namespace Foam
 template<class Type>
 class interpolationCellPoint
 :
-    public fieldInterpolation<Type, interpolationCellPoint<Type>>
+    public fieldInterpolation<Type, interpolationCellPoint<Type>>,
+    public interpolationVolPointInterpolation<Type>
 {
-protected:
-
-    // Protected data
-
-        //- Interpolated volfield
-        const GeometricField<Type, pointPatchField, pointMesh> psip_;
-
-
 public:
 
     //- Runtime type information

src/finiteVolume/interpolation/interpolation/interpolationCellPoint/interpolationCellPointI.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2020 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -35,9 +35,9 @@ inline Type Foam::interpolationCellPoint<Type>::interpolate
     const triFace& faceVertices = cpw.faceVertices();
 
     Type t = this->psi_[cpw.cell()]*weights[0];
-    t += psip_[faceVertices[0]]*weights[1];
-    t += psip_[faceVertices[1]]*weights[2];
-    t += psip_[faceVertices[2]]*weights[3];
+    t += this->psip_[faceVertices[0]]*weights[1];
+    t += this->psip_[faceVertices[1]]*weights[2];
+    t += this->psip_[faceVertices[2]]*weights[3];
 
     return t;
 }
@@ -83,9 +83,9 @@ inline Type Foam::interpolationCellPoint<Type>::interpolate
 
     return
         this->psi_[tetIs.cell()]*coordinates[0]
-      + psip_[triIs[0]]*coordinates[1]
-      + psip_[triIs[1]]*coordinates[2]
-      + psip_[triIs[2]]*coordinates[3];
+      + this->psip_[triIs[0]]*coordinates[1]
+      + this->psip_[triIs[1]]*coordinates[2]
+      + this->psip_[triIs[2]]*coordinates[3];
 }
 
 

src/finiteVolume/interpolation/interpolation/interpolationCellPoint/makeInterpolationCellPoint.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -29,7 +29,8 @@ License
 
 namespace Foam
 {
-    makeInterpolation(interpolationCellPoint);
+    FOR_ALL_FIELD_TYPES(defineInterpolation, interpolationCellPoint);
+    FOR_ALL_FIELD_TYPES(makeInterpolation, interpolationCellPoint);
 }
 
 // ************************************************************************* //

src/finiteVolume/interpolation/interpolation/interpolationCellPointFace/interpolationCellPointFace.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2020 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -40,15 +40,7 @@ Foam::interpolationCellPointFace<Type>::interpolationCellPointFace
 )
 :
     fieldInterpolation<Type, interpolationCellPointFace<Type>>(psi),
-    psip_
-    (
-        volPointInterpolation::New(psi.mesh()).interpolate
-        (
-            psi,
-            "volPointInterpolate(" + psi.name() + ')',
-            true        // use cache
-        )
-    ),
+    interpolationVolPointInterpolation<Type>(psi),
     psis_(linearInterpolate(psi))
 {}
 
@@ -211,7 +203,7 @@ Type Foam::interpolationCellPointFace<Type>::interpolate
         {
             for (label i=0; i<2; i++)
             {
-                ts[i] = psip_[tetPointLabels[i]];
+                ts[i] = this->psip_[tetPointLabels[i]];
             }
 
             if (closestFace < psis_.size())
@@ -301,7 +293,7 @@ Type Foam::interpolationCellPointFace<Type>::interpolate
             // add up the point values ...
             for (label i=0; i<2; i++)
             {
-                Type vel = psip_[tetPointLabels[i]];
+                Type vel = this->psip_[tetPointLabels[i]];
                 t += phi[i]*vel;
             }
 

src/finiteVolume/interpolation/interpolation/interpolationCellPointFace/interpolationCellPointFace.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2020 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -32,7 +32,7 @@ Description
 #ifndef interpolationCellPointFace_H
 #define interpolationCellPointFace_H
 
-#include "interpolation.H"
+#include "interpolationVolPointInterpolation.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -46,13 +46,11 @@ namespace Foam
 template<class Type>
 class interpolationCellPointFace
 :
-    public fieldInterpolation<Type, interpolationCellPointFace<Type>>
+    public fieldInterpolation<Type, interpolationCellPointFace<Type>>,
+    public interpolationVolPointInterpolation<Type>
 {
     // Private Data
 
-        //- Interpolated volfield
-        const GeometricField<Type, pointPatchField, pointMesh> psip_;
-
         //- Linearly interpolated volfield
         const GeometricField<Type, fvsPatchField, surfaceMesh> psis_;
 

src/finiteVolume/interpolation/interpolation/interpolationCellPointFace/makeInterpolationCellPointFace.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -29,7 +29,8 @@ License
 
 namespace Foam
 {
-    makeInterpolation(interpolationCellPointFace);
+    FOR_ALL_FIELD_TYPES(defineInterpolation, interpolationCellPointFace);
+    FOR_ALL_FIELD_TYPES(makeInterpolation, interpolationCellPointFace);
 }
 
 // ************************************************************************* //

src/finiteVolume/interpolation/interpolation/interpolationCellPointWallModified/makeInterpolationCellPointWallModified.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -29,7 +29,16 @@ License
 
 namespace Foam
 {
-    makeInterpolation(interpolationCellPointWallModified);
+    FOR_ALL_FIELD_TYPES
+    (
+        defineInterpolation,
+        interpolationCellPointWallModified
+    );
+    FOR_ALL_FIELD_TYPES
+    (
+        makeInterpolation,
+        interpolationCellPointWallModified
+    );
 }
 
 // ************************************************************************* //

src/finiteVolume/interpolation/interpolation/interpolationPointMVC/interpolationPointMVC.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2020 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -35,15 +35,7 @@ Foam::interpolationPointMVC<Type>::interpolationPointMVC
 )
 :
     fieldInterpolation<Type, interpolationPointMVC<Type>>(psi),
-    psip_
-    (
-        volPointInterpolation::New(psi.mesh()).interpolate
-        (
-            psi,
-            "volPointInterpolate(" + psi.name() + ')',
-            true            // use cache
-        )
-    )
+    interpolationVolPointInterpolation<Type>(psi)
 {}
 
 

src/finiteVolume/interpolation/interpolation/interpolationPointMVC/interpolationPointMVC.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2020 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -33,7 +33,7 @@ Description
 #ifndef interpolationPointMVC_H
 #define interpolationPointMVC_H
 
-#include "interpolation.H"
+#include "interpolationVolPointInterpolation.H"
 #include "pointMVCWeight.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -48,16 +48,9 @@ namespace Foam
 template<class Type>
 class interpolationPointMVC
 :
-    public fieldInterpolation<Type, interpolationPointMVC<Type>>
+    public fieldInterpolation<Type, interpolationPointMVC<Type>>,
+    public interpolationVolPointInterpolation<Type>
 {
-protected:
-
-    // Protected data
-
-        //- Interpolated volfield
-        const GeometricField<Type, pointPatchField, pointMesh> psip_;
-
-
 public:
 
     //- Runtime type information

src/finiteVolume/interpolation/interpolation/interpolationPointMVC/interpolationPointMVCI.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2020 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -31,7 +31,7 @@ inline Type Foam::interpolationPointMVC<Type>::interpolate
     const pointMVCWeight& cpw
 ) const
 {
-    return cpw.interpolate(psip_);
+    return cpw.interpolate(this->psip_);
 }
 
 

src/finiteVolume/interpolation/interpolation/interpolationPointMVC/makeInterpolationPointMVC.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -29,7 +29,8 @@ License
 
 namespace Foam
 {
-    makeInterpolation(interpolationPointMVC);
+    FOR_ALL_FIELD_TYPES(defineInterpolation, interpolationPointMVC);
+    FOR_ALL_FIELD_TYPES(makeInterpolation, interpolationPointMVC);
 }
 
 // ************************************************************************* //

src/finiteVolume/interpolation/interpolation/interpolationVolPointInterpolation/interpolationVolPointInterpolation.C

@@ -0,0 +1,62 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "interpolationVolPointInterpolation.H"
+#include "volPointInterpolation.H"
+
+// * * * * * * * * * * * * * * * * Constructor * * * * * * * * * * * * * * * //
+
+template<class Type>
+Foam::interpolationVolPointInterpolation<Type>::
+interpolationVolPointInterpolation
+(
+    const GeometricField<Type, fvPatchField, volMesh>& psi
+)
+:
+    psip_
+    (
+        volPointInterpolation::New(psi.mesh()).interpolate
+        (
+            psi,
+            "volPointInterpolate(" + psi.name() + ')',
+            true        // use cache
+        )
+    )
+{}
+
+
+template<class Type>
+Foam::interpolationVolPointInterpolation<Type>::
+interpolationVolPointInterpolation
+(
+    const GeometricField<Type, fvPatchField, volMesh>& psi,
+    tmp<GeometricField<Type, pointPatchField, pointMesh>> psip
+)
+:
+    psip_(psip)
+{}
+
+
+// ************************************************************************* //

src/finiteVolume/interpolation/interpolation/interpolationVolPointInterpolation/interpolationVolPointInterpolation.H

@@ -0,0 +1,111 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::interpolationVolPointInterpolation
+
+Description
+    Base class for interpolations that require a vol-point interpolated field
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef interpolationVolPointInterpolation_H
+#define interpolationVolPointInterpolation_H
+
+#include "interpolation.H"
+#include "cellPointWeight.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                   Class interpolationVolPointInterpolation Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class Type>
+class interpolationVolPointInterpolation
+{
+protected:
+
+    // Protected data
+
+        //- Interpolated volfield
+        const GeometricField<Type, pointPatchField, pointMesh> psip_;
+
+
+public:
+
+    //- Runtime type information
+    TypeName("interpolationVolPointInterpolation");
+
+
+
+    // Constructors
+
+        //- Construct from components
+        interpolationVolPointInterpolation
+        (
+            const GeometricField<Type, fvPatchField, volMesh>& psi
+        );
+
+        //- Construct from components
+        interpolationVolPointInterpolation
+        (
+            const GeometricField<Type, fvPatchField, volMesh>& psi,
+            tmp<GeometricField<Type, pointPatchField, pointMesh>> psip
+        );
+
+
+    // Destructor
+    virtual ~interpolationVolPointInterpolation()
+    {}
+
+
+    // Member Functions
+
+        //- Access the point field
+        inline const GeometricField<Type, pointPatchField, pointMesh>&
+            psip() const
+        {
+            return psip_;
+        }
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+    #include "interpolationVolPointInterpolation.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/finiteVolume/interpolation/interpolation/interpolationVolPointInterpolation/makeInterpolationVolPointInterpolation.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -21,36 +21,19 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
-Description
-    Cutting plane sampling functionality
-
 \*---------------------------------------------------------------------------*/
 
-#include "cuttingPlane.H"
+#include "interpolationVolPointInterpolation.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-template<class Type>
-Foam::tmp<Foam::Field<Type>> Foam::cuttingPlane::sample
-(
-    const Field<Type>& fld
-) const
+namespace Foam
 {
-    return tmp<Field<Type>>(new Field<Type>(fld, cutCells()));
-}
-
-
-template<class Type>
-Foam::tmp<Foam::Field<Type>> Foam::cuttingPlane::sample
-(
-    const tmp<Field<Type>>& tfld
-) const
-{
-    tmp<Field<Type>> tsf = sample(tfld());
-    tfld.clear();
-    return tsf;
+    FOR_ALL_FIELD_TYPES
+    (
+        defineInterpolation,
+        interpolationVolPointInterpolation
+    );
 }
 
 // ************************************************************************* //

src/meshTools/Make/files

@@ -255,6 +255,10 @@ tetOverlapVolume/tetOverlapVolume.C
 triIntersect/triIntersect.C
 triIntersect/triIntersectLocationIO.C
 
+cutPoly/cellEdgeAddressing.C
+cutPoly/cutPoly.C
+cutPoly/cutPolyIsoSurface.C
+
 nonConformal/polyPatches/nonConformal/nonConformalPolyPatch.C
 nonConformal/polyPatches/nonConformalCoupled/nonConformalCoupledPolyPatch.C
 nonConformal/polyPatches/nonConformalCyclic/nonConformalCyclicPolyPatch.C

src/meshTools/cutPoly/cellEdgeAddressing.C

@@ -0,0 +1,223 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "cellEdgeAddressing.H"
+#include "polyDistributionMap.H"
+#include "polyTopoChangeMap.H"
+#include "polyMeshMap.H"
+#include "HashList.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// An edge hash functor that is much faster than Hash<edge>, but is also more
+// prone to collisions
+struct QuickHashEdge
+{
+    unsigned operator()(const edge& e)
+    {
+        return e[0]*e[1];
+    }
+};
+
+// An invalid null edge for unset entries in the edge map
+template<>
+const edge HashList<label, edge, QuickHashEdge>::nullKey(-labelMax, -labelMax);
+
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::cellEdgeAddressing::cellEdgeAddressing
+(
+    const cell& c,
+    const faceList& fs,
+    const bool cOwnsFirst
+)
+{
+    // Compute the number of cell edges
+    label nCellEdges = 0;
+    forAll(c, cfi)
+    {
+        const face& f = fs[c[cfi]];
+
+        nCellEdges += f.size();
+    }
+    nCellEdges /= 2;
+
+    // Construct a map to enumerate the cell edges
+    HashList<label, edge, QuickHashEdge> edgeToCei(nCellEdges*6);
+    {
+        label cellEdgei = 0;
+
+        forAll(c, cfi)
+        {
+            const face& f = fs[c[cfi]];
+
+            forAll(f, fei)
+            {
+                cellEdgei += edgeToCei.insert(f.faceEdge(fei), cellEdgei);
+            }
+        }
+    }
+
+    // Allocate and initialise the addressing
+    cfiAndFeiToCei_ = labelListList(c.size());
+    ceiToCfiAndFei_ =
+        List<Pair<labelPair>>
+        (
+            nCellEdges,
+            Pair<labelPair>(labelPair(-1, -1), labelPair(-1, -1))
+        );
+
+    // Copy data out of the map into the addressing
+    forAll(c, cfi)
+    {
+        const face& f = fs[c[cfi]];
+
+        cfiAndFeiToCei_[cfi] = labelList(f.size(), -1);
+
+        forAll(f, fei)
+        {
+            const label cei = edgeToCei[f.faceEdge(fei)];
+
+            cfiAndFeiToCei_[cfi][fei] = cei;
+            ceiToCfiAndFei_[cei]
+            [
+                ceiToCfiAndFei_[cei][0] != labelPair(-1, -1)
+            ] = {cfi, fei};
+        }
+    }
+
+    // Allocate and initialise the face signs
+    cOwns_ = boolList(c.size(), false);
+    cOwns_[0] = cOwnsFirst;
+    boolList cOwnsIsSet(c.size(), false);
+    cOwnsIsSet[0] = true;
+
+    // Compare cell-face-edges to determine face signs
+    forAll(c, cfi)
+    {
+        const face& f = fs[c[cfi]];
+
+        if (!cOwnsIsSet[cfi]) continue;
+
+        forAll(f, fei)
+        {
+            const label cei = cfiAndFeiToCei_[cfi][fei];
+
+            const labelPair& other =
+                ceiToCfiAndFei_[cei]
+                [
+                    ceiToCfiAndFei_[cei][0] == labelPair(cfi, fei)
+                ];
+            const label cfj = other[0], fej = other[1];
+
+            if (cOwnsIsSet[cfj]) continue;
+
+            const label sign =
+                edge::compare
+                (
+                    f.faceEdge(fei),
+                    fs[c[cfj]].faceEdge(fej)
+                );
+
+            cOwns_[cfj] = sign < 0 ? cOwns_[cfi] : !cOwns_[cfi];
+            cOwnsIsSet[cfj] = true;
+        }
+    }
+}
+
+
+Foam::cellEdgeAddressingList::cellEdgeAddressingList(const polyMesh& mesh)
+:
+    MeshObject<polyMesh, UpdateableMeshObject, cellEdgeAddressingList>(mesh),
+    list_(mesh.nCells())
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+bool Foam::cellEdgeAddressingList::movePoints()
+{
+    return true;
+}
+
+
+void Foam::cellEdgeAddressingList::distribute(const polyDistributionMap& map)
+{
+    // We could be more selective here and try to keep addressing for cells
+    // that haven't changed
+
+    list_.clear();
+    list_.resize(map.mesh().nCells());
+}
+
+
+void Foam::cellEdgeAddressingList::topoChange(const polyTopoChangeMap& map)
+{
+    // We could be more selective here and try to keep addressing for cells
+    // that haven't changed
+
+    list_.clear();
+    list_.resize(map.mesh().nCells());
+}
+
+
+void Foam::cellEdgeAddressingList::mapMesh(const polyMeshMap& map)
+{
+    list_.clear();
+    list_.resize(map.mesh().nCells());
+}
+
+
+// * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * * //
+
+const Foam::cellEdgeAddressing& Foam::cellEdgeAddressingList::operator[]
+(
+    const label celli
+) const
+{
+    if (!list_.set(celli))
+    {
+        const cell& c = mesh().cells()[celli];
+        const faceList& fs = mesh().faces();
+        const labelList& fOwners = mesh().faceOwner();
+
+        list_.set
+        (
+            celli,
+            new cellEdgeAddressing(c, fs, fOwners[c[0]] == celli)
+        );
+    }
+
+    return list_[celli];
+}
+
+
+// ************************************************************************* //

src/meshTools/cutPoly/cellEdgeAddressing.H

@@ -0,0 +1,163 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::cellEdgeAddressing
+    Foam::cellEdgeAddressingList
+
+Description
+    Engine for providing cell-local cell-edge to face-edge addressing
+
+SourceFiles
+    cellEdgeAddressing.C
+    cellEdgeAddressingI.H
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef cellEdgeAddressing_H
+#define cellEdgeAddressing_H
+
+#include "MeshObject.H"
+#include "polyMesh.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+class polyDistributionMap;
+class polyTopoChangeMap;
+class polyMeshMap;
+
+/*---------------------------------------------------------------------------*\
+                      Class cellEdgeAddressing Declaration
+\*---------------------------------------------------------------------------*/
+
+class cellEdgeAddressing
+{
+private:
+
+    // Private Data
+
+        //- Map from cell-face-index and face-edge-index to cell-edge-index
+        labelListList cfiAndFeiToCei_;
+
+        //- Map from cell-edge-index to cell-face-index and face-edge-index.
+        //  Note that there are two sets of indices for each cell-edge,
+        //  corresponding to the two faces that are connected to that edge.
+        List<Pair<labelPair>> ceiToCfiAndFei_;
+
+        //- For each cell-face, whether or not the cell owns it. The same as
+        //  the concept of owner in the polyMesh, except this is calculated
+        //  using the direction of the face-edges. I.e., if the first face is
+        //  know to be "owned" by the face then any face that is edge connected
+        //  to the first face is also owned by the cell if the connecting
+        //  face-edge is numbered in an opposing direction.
+        boolList cOwns_;
+
+
+public:
+
+    // Constructors
+
+        //- Construct from a cell and its faces
+        cellEdgeAddressing
+        (
+            const cell& c,
+            const faceList& fs,
+            const bool cOwnsFirst
+        );
+
+
+    // Member Functions
+
+        //- Map from cell-face-index and face-edge-index to cell-edge-index
+        inline const labelListList& cfiAndFeiToCei() const;
+
+        //- Map from cell-edge-index to cell-face-index and face-edge-index
+        inline const List<Pair<labelPair>>& ceiToCfiAndFei() const;
+
+        //- For each cell-face, whether or not the cell owns it
+        inline const boolList& cOwns() const;
+};
+
+
+/*---------------------------------------------------------------------------*\
+                    Class cellEdgeAddressingList Declaration
+\*---------------------------------------------------------------------------*/
+
+class cellEdgeAddressingList
+:
+    public MeshObject<polyMesh, UpdateableMeshObject, cellEdgeAddressingList>
+{
+private:
+
+    // Private Data
+
+        //- Edge addressing for each cell
+        mutable PtrList<cellEdgeAddressing> list_;
+
+
+public:
+
+    // Constructors
+
+        //- Construct for a mesh
+        cellEdgeAddressingList(const polyMesh& mesh);
+
+
+    // Member Functions
+
+        //- Update following mesh motion
+        virtual bool movePoints();
+
+        //- Update following mesh distribution
+        virtual void distribute(const polyDistributionMap& map);
+
+        //- Update following topology change
+        virtual void topoChange(const polyTopoChangeMap& map);
+
+        //- Update following mapping
+        virtual void mapMesh(const polyMeshMap& map);
+
+
+    // Member Operators
+
+        //- Get the addressing for a given cell
+        const cellEdgeAddressing& operator[](const label celli) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#include "cellEdgeAddressingI.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/meshTools/cutPoly/cellEdgeAddressingI.H

@@ -0,0 +1,50 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "cellEdgeAddressing.H"
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+inline const Foam::labelListList&
+Foam::cellEdgeAddressing::cfiAndFeiToCei() const
+{
+    return cfiAndFeiToCei_;
+}
+
+
+inline const Foam::List<Foam::Pair<Foam::labelPair>>&
+Foam::cellEdgeAddressing::ceiToCfiAndFei() const
+{
+    return ceiToCfiAndFei_;
+}
+
+
+inline const Foam::boolList& Foam::cellEdgeAddressing::cOwns() const
+{
+    return cOwns_;
+}
+
+
+// ************************************************************************* //

src/meshTools/cutPoly/cutPoly.C

@@ -0,0 +1,314 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "cutPolyValue.H"
+#include "OBJstream.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+template<class Type>
+Type min(const UIndirectList<Type>& l)
+{
+    Type result = pTraits<Type>::max;
+    forAll(l, i)
+    {
+        result = min(result, l[i]);
+    }
+    return result;
+}
+
+template<class Type>
+Type max(const UIndirectList<Type>& l)
+{
+    Type result = pTraits<Type>::min;
+    forAll(l, i)
+    {
+        result = max(result, l[i]);
+    }
+    return result;
+}
+
+}
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+Foam::List<Foam::labelPair> Foam::cutPoly::faceCuts
+(
+    const face& f,
+    const scalarField& pAlphas,
+    const scalar isoAlpha
+)
+{
+    UIndirectList<scalar> fAlphas(pAlphas, f);
+
+    // Quick reject if all alpha values are above or below the iso value
+    if (min(fAlphas) > isoAlpha || isoAlpha > max(fAlphas))
+    {
+        return List<labelPair>();
+    }
+
+    // Find the starting point
+    label fpi0 = 0;
+    while ((fAlphas[fpi0] < isoAlpha) == separatedBelow)
+    {
+        ++ fpi0;
+    }
+
+    // Create cuts on every edge for which there is a sign change
+    DynamicList<labelPair> cuts(1);
+    label cuti = 0;
+    forAll(f, i)
+    {
+        const label fpi1 = f.fcIndex(fpi0);
+
+        if ((fAlphas[fpi0] < isoAlpha) != (fAlphas[fpi1] < isoAlpha))
+        {
+            if (cuti == 0)
+            {
+                cuts.append({fpi0, -1});
+            }
+            else
+            {
+                cuts.last()[1] = fpi0;
+            }
+
+            cuti = cuts.last().fcIndex(cuti);
+        }
+
+        fpi0 = fpi1;
+    }
+
+    // There should have been an even number of cuts
+    if (cuti != 0)
+    {
+        FatalErrorInFunction
+            << "Cutting values " << fAlphas << " with iso-value " << isoAlpha
+            << " resulted in an odd number of cuts"
+            << exit(FatalError);
+    }
+
+    cuts.shrink();
+    return cuts;
+}
+
+
+Foam::labelListList Foam::cutPoly::cellCuts
+(
+    const cell& c,
+    const cellEdgeAddressing& cAddr,
+    const faceUList& fs,
+    const List<List<labelPair>>& fCuts,
+    const scalarField& pAlphas,
+    const scalar isoAlpha
+)
+{
+    // Quick return if not cut
+    label nCellFaceCuts = 0;
+    forAll(c, cfi)
+    {
+        nCellFaceCuts += fCuts[c[cfi]].size();
+    }
+    if (nCellFaceCuts == 0)
+    {
+        return labelListList();
+    }
+
+    // Get local cell-face-edge addressing
+    const labelListList& cfiAndFeiToCei = cAddr.cfiAndFeiToCei();
+    const List<Pair<labelPair>>& ceiToCfiAndFei = cAddr.ceiToCfiAndFei();
+    const boolList& cOwns = cAddr.cOwns();
+
+    // For each cell-face and face-edge, get the face-edge at the other end of
+    // the edge's cut, or -1 if the edge is not cut. This allows us to walk
+    // along face cuts. It also doubles as the "visited" array during the walk.
+    // As cuts are visited, the relevant label in this array is set to -1 to
+    // ensure that the cut is not considered twice.
+    labelListList cellFaceAndFaceEdgeCutToFaceEdge(c.size());
+    forAll(c, cfi)
+    {
+        cellFaceAndFaceEdgeCutToFaceEdge[cfi] =
+            labelList(fs[c[cfi]].size(), -1);
+
+        forAll(fCuts[c[cfi]], fci)
+        {
+            const label fei0 = fCuts[c[cfi]][fci].first();
+            const label fei1 = fCuts[c[cfi]][fci].second();
+
+            if (cOwns[cfi])
+            {
+                cellFaceAndFaceEdgeCutToFaceEdge[cfi][fei0] = fei1;
+            }
+            else
+            {
+                cellFaceAndFaceEdgeCutToFaceEdge[cfi][fei1] = fei0;
+            }
+        }
+    }
+
+    // Walk around the face cuts to generate the cell cuts
+    DynamicList<labelList> cuts;
+    {
+        label nCellFaceCutsVisited = 0;
+        label cfi0 = 0, fei0 = 0;
+        while (nCellFaceCutsVisited < nCellFaceCuts)
+        {
+            // Find the next unvisited connection
+            bool found = false;
+            for (label cfj0 = cfi0; cfj0 < c.size(); ++ cfj0)
+            {
+                for (label fej0 = fei0; fej0 < fs[c[cfj0]].size(); ++ fej0)
+                {
+                    const label fej1 =
+                        cellFaceAndFaceEdgeCutToFaceEdge[cfj0][fej0];
+
+                    if (fej1 != -1)
+                    {
+                        cfi0 = cfj0;
+                        fei0 = fej0;
+                        found = true;
+                    }
+
+                    if (found) break;
+                }
+
+                if (found) break;
+
+                fei0 = 0;
+            }
+            if (!found)
+            {
+                FatalErrorInFunction
+                    << "Could not find next unvisited connection for cell cut"
+                    << exit(FatalError);
+            }
+
+            // Walk around the cell from face to face to form a cut
+            label cfi = cfi0, fei = fei0;
+            DynamicList<label> cut(8);
+            while (cellFaceAndFaceEdgeCutToFaceEdge[cfi][fei] != -1)
+            {
+                ++ nCellFaceCutsVisited;
+
+                const label otherFei =
+                    cellFaceAndFaceEdgeCutToFaceEdge[cfi][fei];
+
+                cellFaceAndFaceEdgeCutToFaceEdge[cfi][fei] = -1;
+
+                const label cei = cfiAndFeiToCei[cfi][otherFei];
+
+                cut.append(cei);
+
+                const labelPair next =
+                    ceiToCfiAndFei[cei].other({cfi, otherFei});
+
+                cfi = next.first();
+                fei = next.second();
+            }
+
+            // Copy the data over into the result list
+            cuts.append(labelList());
+            cuts.last().transfer(cut);
+        }
+    }
+
+    cuts.shrink();
+    return cuts;
+}
+
+
+void Foam::cutPoly::writeFaceCuts
+(
+    const face& f,
+    const List<labelPair>& fCuts,
+    const pointField& ps,
+    const scalarField& pAlphas,
+    const scalar isoAlpha,
+    OBJstream& obj
+)
+{
+    forAll(fCuts, cuti)
+    {
+        pointField cutPs(2);
+
+        forAll(fCuts[cuti], i)
+        {
+            const edge e = f.faceEdge(fCuts[cuti][i]);
+
+            cutPs[i] = edgeCutValue(e, pAlphas, isoAlpha, ps);
+        }
+
+        obj.write(edge(0, 1), cutPs);
+    }
+}
+
+
+void Foam::cutPoly::writeCellCuts
+(
+    const cell& c,
+    const cellEdgeAddressing& cAddr,
+    const List<List<label>>& cCuts,
+    const faceList& fs,
+    const pointField& ps,
+    const scalarField& pAlphas,
+    const scalar isoAlpha,
+    OBJstream& obj
+)
+{
+    // Quick return if not cut
+    if (cCuts.size() == 0)
+    {
+        return;
+    }
+
+    // Get local cell-face-edge addressing
+    const List<Pair<labelPair>>& ceiToCfiAndFei = cAddr.ceiToCfiAndFei();
+
+    // Write out each cut as a single face
+    forAll(cCuts, cuti)
+    {
+        if (cCuts[cuti].size() < 3) continue;
+
+        pointField cutPs(cCuts[cuti].size());
+
+        forAll(cCuts[cuti], i)
+        {
+            const label cei = cCuts[cuti][i];
+            const label cfi = ceiToCfiAndFei[cei][0][0];
+            const label fei = ceiToCfiAndFei[cei][0][1];
+
+            const edge e = fs[c[cfi]].faceEdge(fei);
+
+            cutPs[i] = edgeCutValue(e, pAlphas, isoAlpha, ps);
+        }
+
+        obj.write(face(identity(cCuts[cuti].size())), cutPs, false);
+    }
+}
+
+
+// ************************************************************************* //

src/meshTools/cutPoly/cutPoly.H

@@ -0,0 +1,122 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Namespace
+    Foam::cutPoly
+
+Description
+    Low level functions for cutting poly faces and cells
+
+    Cutting a face using an iso-surface defined at the points is done by
+    walking around the face. A start point is chosen and then the walk proceeds
+    until a "cut" edge is found with points on either side of the iso-surface.
+    This edge is noted and the walk continues until another cut edge is found.
+    These two edge cuts, along with all points in-between, form a sub-face on
+    one side of the iso-surface. This algorithm continues around the entire
+    face, potentially cutting off multiple sub-faces. The remainder is a single
+    contiguous face on the other side of the iso-surface.
+
+    Cutting a cell is then done by walking from cut-edge to cut-edge around the
+    cell's faces.
+
+SourceFiles
+    cutPoly.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef cutPoly_H
+#define cutPoly_H
+
+#include "cellEdgeAddressing.H"
+#include "PrimitivePatch.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+class OBJstream;
+
+namespace cutPoly
+{
+
+//- This flag determines which side of the iso-surface is separated into
+//  multiple sub-faces, and which side is kept as the central contiguous face
+static const bool separatedBelow = true;
+
+//- Return the cuts for a given face. This returns a list of pairs of
+//  face-edge indices. Both edges in the pair are cut, and connecting these
+//  cuts creates an edge of the iso-surface.
+List<labelPair> faceCuts
+(
+    const face& f,
+    const scalarField& pAlphas,
+    const scalar isoAlpha
+);
+
+//- Return the cuts for a given cell. This returns a list of lists of cell-edge
+//  indices. Each sub-list represents a single closed loop of cell-edges. The
+//  edges listed are all cut, and connecting these cuts up creates a face of
+//  the iso-surface.
+labelListList cellCuts
+(
+    const cell& c,
+    const cellEdgeAddressing& cAddr,
+    const faceUList& fs,
+    const List<List<labelPair>>& fCuts,
+    const scalarField& pAlphas,
+    const scalar isoAlpha
+);
+
+//- Write the cuts for a given face to an OBJ file
+void writeFaceCuts
+(
+    const face& f,
+    const List<labelPair>& fCuts,
+    const pointField& ps,
+    const scalarField& pAlphas,
+    const scalar isoAlpha,
+    OBJstream& obj
+);
+
+//- Write the cuts for a given cell to an OBJ file
+void writeCellCuts
+(
+    const cell& c,
+    const cellEdgeAddressing& cAddr,
+    const List<List<label>>& cCuts,
+    const faceList& fs,
+    const pointField& ps,
+    const scalarField& pAlphas,
+    const scalar isoAlpha,
+    OBJstream& obj
+);
+
+} // End namespace cutPoly
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/meshTools/cutPoly/cutPolyIntegral.H

@@ -0,0 +1,115 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::cutPolyIntegral
+
+Description
+    Functions for computing integrals over cut faces and cells
+
+SourceFiles
+    cutPolyIntegralTemplates.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef cutPolyIntegral_H
+#define cutPolyIntegral_H
+
+#include "cutPolyValue.H"
+
+#include <tuple>
+#include <type_traits>
+#include <utility>
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace cutPoly
+{
+
+//- Compute the face-area and face-area-integrals of the given properties over
+//  a face
+template<template<class> class FaceValues, class ... Types>
+std::tuple<vector, typename outerProduct<vector, Types>::type ...>
+faceAreaIntegral
+(
+    const FaceValues<point>& fPs,
+    const point& fPAvg,
+    const std::tuple<FaceValues<Types> ...>& fPsis,
+    const std::tuple<Types ...>& fPsiAvg
+);
+
+//- Compute the face-area and face-area-integral of the given property over
+//  a cut-face. ToDo: Make variadic as faceAreaIntegral.
+template<class Type>
+Tuple2<vector, typename outerProduct<vector, Type>::type> faceCutAreaIntegral
+(
+    const face& f,
+    const vector& fArea,
+    const Type& fPsi,
+    const List<labelPair>& fCuts,
+    const pointField& ps,
+    const Field<Type>& pPsis,
+    const scalarField& pAlphas,
+    const scalar isoAlpha,
+    const bool below
+);
+
+//- Compute the cell-volume and cell-volume-integral of the given property over
+//  a cut-cell. ToDo: Make variadic as faceAreaIntegral.
+template<class Type>
+Tuple2<scalar, Type> cellCutVolumeIntegral
+(
+    const cell& c,
+    const cellEdgeAddressing& cAddr,
+    const scalar cVolume,
+    const Type& cPsi,
+    const labelListList& cCuts,
+    const faceUList& fs,
+    const vectorField& fAreas,
+    const vectorField& fCentres,
+    const vectorField& fPsis,
+    const vectorField& fCutAreas,
+    const vectorField& fCutPsis,
+    const pointField& ps,
+    const Field<Type>& pPsis,
+    const scalarField& pAlphas,
+    const scalar isoAlpha,
+    const bool below
+);
+
+} // End namespace cutPoly
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+    #include "cutPolyIntegralTemplates.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/meshTools/cutPoly/cutPolyIntegralTemplates.C

@@ -0,0 +1,498 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2018-2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "cutPolyIntegral.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace cutPoly
+{
+
+template<class Type>
+Type average(const UIndirectList<Type>& l)
+{
+    Type result = pTraits<Type>::zero;
+    forAll(l, i)
+    {
+        result += l[i];
+    }
+    return result/l.size();
+}
+
+template<class Type>
+Type average(const List<Type>& xs, const labelUList& is)
+{
+    return average(UIndirectList<Type>(xs, is));
+}
+
+template<class Container>
+typename Container::value_type iterableAverage(const Container& xs)
+{
+    label n = 0;
+    typename Container::value_type nResult =
+        pTraits<typename Container::value_type>::zero;
+
+    forAllConstIter(typename Container, xs, iter)
+    {
+        ++ n;
+        nResult += *iter;
+    }
+
+    return nResult/n;
+}
+
+} // End namespace cutPoly
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace cutPoly
+{
+
+template<class Op, class Tuple, class Int, Int ... Is>
+auto tupleOp
+(
+    const Tuple& tuple,
+    const std::integer_sequence<Int, Is ...>&,
+    const Op& op
+)
+{
+    return std::make_tuple(op(std::get<Is>(tuple)) ...);
+}
+
+
+template<class Op, class ... Types>
+auto tupleOp(const std::tuple<Types...>& tuple, const Op& op)
+{
+    return tupleOp(tuple, std::make_index_sequence<sizeof ... (Types)>(), op);
+}
+
+
+struct OpBegin
+{
+    template<class Type>
+    auto operator()(const Type& x) const
+    {
+        return x.begin();
+    }
+};
+
+
+struct OpDereference
+{
+    template<class Type>
+    auto operator()(const Type& x) const
+    {
+        return *x;
+    }
+};
+
+
+struct OpNext
+{
+    template<class Type>
+    auto operator()(const Type& x) const
+    {
+        return x.next();
+    }
+};
+
+
+template<class ScaleType>
+struct OpScaled
+{
+    const ScaleType s_;
+
+    OpScaled(const ScaleType& s)
+    :
+        s_(s)
+    {}
+
+    template<class Type>
+    auto operator()(const Type& x) const
+    {
+        return s_*x;
+    }
+};
+
+
+template<class Op, class Tuple, class Int, Int ... Is>
+void tupleInPlaceOp
+(
+    Tuple& tuple,
+    const std::integer_sequence<Int, Is ...>&,
+    const Op& op
+)
+{
+    (void)std::initializer_list<nil>
+    {(
+        op(std::get<Is>(tuple)),
+        nil()
+    ) ... };
+}
+
+
+template<class Op, class ... Types>
+void tupleInPlaceOp(std::tuple<Types...>& tuple, const Op& op)
+{
+    tupleInPlaceOp(tuple, std::make_index_sequence<sizeof ... (Types)>(), op);
+}
+
+
+struct InPlaceOpAdvance
+{
+    template<class Type>
+    void operator()(Type& x) const
+    {
+        ++ x;
+    }
+};
+
+
+template<class BinaryOp, class Tuple, class Int, Int ... Is>
+auto tupleBinaryOp
+(
+    const Tuple& tupleA,
+    const Tuple& tupleB,
+    const std::integer_sequence<Int, Is ...>&,
+    const BinaryOp& bop
+)
+{
+    return std::make_tuple(bop(std::get<Is>(tupleA), std::get<Is>(tupleB)) ...);
+}
+
+
+template<class BinaryOp, class ... TypesA, class ... TypesB>
+auto tupleBinaryOp
+(
+    const std::tuple<TypesA ...>& tupleA,
+    const std::tuple<TypesB ...>& tupleB,
+    const BinaryOp& bop
+)
+{
+    return tupleBinaryOp
+    (
+        tupleA,
+        tupleB,
+        std::make_index_sequence<sizeof ... (TypesA)>(),
+        bop
+    );
+}
+
+
+struct BinaryOpAdd
+{
+    template<class Type>
+    auto operator()(const Type& a, const Type& b) const
+    {
+        return a + b;
+    }
+};
+
+} // End namespace cutPoly
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+template<template<class> class FaceValues, class ... Types>
+std::tuple
+<
+    Foam::vector,
+    typename Foam::outerProduct<Foam::vector, Types>::type ...
+>
+Foam::cutPoly::faceAreaIntegral
+(
+    const FaceValues<point>& fPs,
+    const point& fPAvg,
+    const std::tuple<FaceValues<Types> ...>& fPsis,
+    const std::tuple<Types ...>& fPsiAvg
+)
+{
+    vector fCutsArea = Zero;
+    auto fCutsAreaPsi =
+        std::make_tuple(typename outerProduct<vector, Types>::type(Zero) ...);
+
+    typename FaceValues<point>::const_iterator fPIter(fPs.begin());
+    auto fPsiIter = tupleOp(fPsis, OpBegin());
+
+    for(; fPIter != fPs.end(); ++ fPIter)
+    {
+        const point p0 = *fPIter;
+        const point p1 = fPIter.next();
+        auto psi0 = tupleOp(fPsiIter, OpDereference());
+        auto psi1 = tupleOp(fPsiIter, OpNext());
+
+        const vector a = ((p1 - p0)^(fPAvg - p0))/2;
+        auto v =
+            tupleBinaryOp
+            (
+                psi0,
+                tupleBinaryOp
+                (
+                    psi1,
+                    fPsiAvg,
+                    BinaryOpAdd()
+                ),
+                BinaryOpAdd()
+            );
+
+        auto av = tupleOp(v, OpScaled<vector>(a/3));
+
+        fCutsArea += a;
+        fCutsAreaPsi = tupleBinaryOp(fCutsAreaPsi, av, BinaryOpAdd());
+
+        tupleInPlaceOp(fPsiIter, InPlaceOpAdvance());
+    }
+
+    return std::tuple_cat(std::tuple<vector>(fCutsArea), fCutsAreaPsi);
+}
+
+
+template<class Type>
+Foam::Tuple2
+<
+    Foam::vector,
+    typename Foam::outerProduct<Foam::vector, Type>::type
+>
+Foam::cutPoly::faceCutAreaIntegral
+(
+    const face& f,
+    const vector& fArea,
+    const Type& fPsi,
+    const List<labelPair>& fCuts,
+    const pointField& ps,
+    const Field<Type>& pPsis,
+    const scalarField& pAlphas,
+    const scalar isoAlpha,
+    const bool below
+)
+{
+    typedef typename outerProduct<vector, Type>::type IntegralType;
+
+    // If there are no cuts return either the entire face or zero, depending on
+    // which side of the iso-surface the face is
+    if (fCuts.size() == 0)
+    {
+        if ((pAlphas[f[0]] < isoAlpha) == below)
+        {
+            return Tuple2<vector, IntegralType>(fArea, fArea*fPsi);
+        }
+        else
+        {
+            return Tuple2<vector, IntegralType>(Zero, Zero);
+        }
+    }
+
+    auto result =
+        faceAreaIntegral
+        (
+            cutPoly::FaceCutValues<Type>
+            (
+                f,
+                fCuts,
+                ps,
+                pAlphas,
+                isoAlpha,
+                below
+            ),
+            average(ps, f),
+            std::make_tuple
+            (
+                cutPoly::FaceCutValues<Type>
+                (
+                    f,
+                    fCuts,
+                    pPsis,
+                    pAlphas,
+                    isoAlpha,
+                    below
+                )
+            ),
+            std::make_tuple(average(pPsis, f))
+        );
+
+    return
+        Tuple2<vector, IntegralType>
+        (
+            std::get<0>(result),
+            std::get<1>(result)
+        );
+}
+
+
+template<class Type>
+Foam::Tuple2<Foam::scalar, Type> Foam::cutPoly::cellCutVolumeIntegral
+(
+    const cell& c,
+    const cellEdgeAddressing& cAddr,
+    const scalar cVolume,
+    const Type& cPsi,
+    const labelListList& cCuts,
+    const faceUList& fs,
+    const vectorField& fAreas,
+    const vectorField& fCentres,
+    const vectorField& fPsis,
+    const vectorField& fCutAreas,
+    const vectorField& fCutPsis,
+    const pointField& ps,
+    const Field<Type>& pPsis,
+    const scalarField& pAlphas,
+    const scalar isoAlpha,
+    const bool below
+)
+{
+    // If there are no cuts return either the entire cell or zero, depending on
+    // which side of the iso-surface the cell is
+    if (cCuts.size() == 0)
+    {
+        if ((pAlphas[fs[c[0]][0]] < isoAlpha) == below)
+        {
+            return Tuple2<scalar, Type>(cVolume, cVolume*cPsi);
+        }
+        else
+        {
+            return Tuple2<scalar, Type>(0, Zero);
+        }
+    }
+
+    // Averages
+    const point cPAvg = average(fCentres, c);
+    const Type cPsiAvg = average(fPsis, c);
+
+    // Initialise totals
+    scalar cCutsVolume = 0;
+    Type cCutsVolumePsi = Zero;
+
+    // Face contributions
+    forAll(c, cfi)
+    {
+        // We use the un-cut face's centroid as the base of the pyramid formed
+        // by the cut face. This is potentially less exact, but it is
+        // consistent with the un-cut cell volume calculation. This means if
+        // you run this function with both values of "below" then the result
+        // will exactly sum to the volume of the overall cell. If we used the
+        // cut-face's centroid this would not be the case.
+        const vector& fBaseP = fCentres[c[cfi]];
+
+        const scalar pyrVolume =
+            (cAddr.cOwns()[cfi] ? +1 : -1)
+           *(fCutAreas[c[cfi]] & (fBaseP - cPAvg))/3;
+
+        cCutsVolume += pyrVolume;
+        cCutsVolumePsi += pyrVolume*(3*fCutPsis[c[cfi]] + cPsiAvg)/4;
+    }
+
+    // Cut contributions
+    {
+        const cutPoly::CellCutValues<point> cPValues
+        (
+            c,
+            cAddr,
+            cCuts,
+            fs,
+            ps,
+            pAlphas,
+            isoAlpha
+        );
+
+        const cutPoly::CellCutValues<Type> cPsiValues
+        (
+            c,
+            cAddr,
+            cCuts,
+            fs,
+            pPsis,
+            pAlphas,
+            isoAlpha
+        );
+
+        // This method is more exact, as it uses the true centroid of the cell
+        // cut. However, to obtain that centroid we have to divide by the area
+        // magnitude, so this can't be generalised to types that do not support
+        // division.
+        /*
+        auto fSumPPsis =
+            faceAreaIntegral
+            (
+                cPValues,
+                iterableAverage(cPValues),
+                std::make_tuple
+                (
+                    cPValues,
+                    cPsiValues
+                ),
+                std::make_tuple
+                (
+                    iterableAverage(cPValues),
+                    iterableAverage(cPsiValues)
+                )
+            );
+
+        const vector& fCutArea = std::get<0>(fSumPPsis);
+        const scalar fMagSqrCutArea = magSqr(fCutArea);
+
+        const point fCutCentre =
+            fMagSqrCutArea > vSmall
+          ? (fCutArea & std::get<1>(fSumPPsis))/fMagSqrCutArea
+          : cPAvg;
+        const Type fCutPsi =
+            fMagSqrCutArea > vSmall
+          ? (fCutArea & std::get<2>(fSumPPsis))/fMagSqrCutArea
+          : cPsiAvg;
+        */
+
+        // This method is more approximate, as it uses point averages rather
+        // than centroids. This does not involve division, though, so this
+        // could be used with types like polynomials that only support addition
+        // and multiplication.
+        auto fSumPPsis =
+            faceAreaIntegral
+            (
+                cPValues,
+                iterableAverage(cPValues),
+                std::make_tuple(),
+                std::make_tuple()
+            );
+
+        const vector& fCutArea = std::get<0>(fSumPPsis);
+        const point fCutCentre = iterableAverage(cPValues);
+        const Type fCutPsi = iterableAverage(cPsiValues);
+
+        const scalar pyrVolume =
+            (below ? +1 : -1)*(fCutArea & (fCutCentre - cPAvg))/3;
+
+        cCutsVolume += pyrVolume;
+        cCutsVolumePsi += pyrVolume*(3*fCutPsi + cPsiAvg)/4;
+    }
+
+    return Tuple2<scalar, Type>(cCutsVolume, cCutsVolumePsi);
+}
+
+
+// ************************************************************************* //

src/meshTools/cutPoly/cutPolyIsoSurface.C

@@ -0,0 +1,238 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "cutPolyIsoSurface.H"
+#include "cellEdgeAddressing.H"
+#include "cutPolyValue.H"
+#include "EdgeMap.H"
+#include "cpuTime.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(cutPolyIsoSurface, 0);
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::cutPolyIsoSurface::cutPolyIsoSurface
+(
+    const polyMesh& mesh,
+    const scalarField& pAlphas,
+    const scalar isoAlpha,
+    const labelList& zoneIDs
+)
+:
+    points_(),
+    pointEdges_(),
+    pointEdgeLambdas_(),
+    faces_(),
+    faceCells_()
+{
+    cpuTime cpuTime;
+
+    // Request the cell-edge addressing engine
+    const cellEdgeAddressingList& cAddrs = cellEdgeAddressingList::New(mesh);
+
+    // Cut the faces
+    List<List<labelPair>> faceCuts(mesh.faces().size());
+    forAll(mesh.faces(), facei)
+    {
+        faceCuts[facei] =
+            cutPoly::faceCuts(mesh.faces()[facei], pAlphas, isoAlpha);
+    }
+
+    // Cut the cells
+    List<labelListList> cellCuts(mesh.cells().size());
+    label nCutCells = 0;
+    auto cutCell = [&](const label celli)
+    {
+        cellCuts[celli] =
+            cutPoly::cellCuts
+            (
+                mesh.cells()[celli],
+                cAddrs[celli],
+                mesh.faces(),
+                faceCuts,
+                pAlphas,
+                isoAlpha
+            );
+        nCutCells += !cellCuts[celli].empty();
+    };
+    if (!isNull<labelList>(zoneIDs))
+    {
+        forAll(zoneIDs, i)
+        {
+            forAll(mesh.cellZones()[zoneIDs[i]], zoneCelli)
+            {
+                cutCell(mesh.cellZones()[zoneIDs[i]][zoneCelli]);
+            }
+        }
+    }
+    else
+    {
+        forAll(mesh.cells(), celli)
+        {
+            cutCell(celli);
+        }
+    }
+
+    // Generate the cell cut polygons
+    const label nAllocate = nCutCells + nCutCells/10;
+    EdgeMap<label> meshEdgePoint(nAllocate*4);
+    DynamicList<point> pointsDyn(nAllocate);
+    DynamicList<edge> pointEdgesDyn(nAllocate);
+    DynamicList<scalar> pointEdgeLambdasDyn(nAllocate);
+    DynamicList<face> facesDyn(nAllocate);
+    DynamicList<label> faceCellsDyn(nAllocate);
+    forAll(mesh.cells(), celli)
+    {
+        const cell& c = mesh.cells()[celli];
+
+        const labelListList& cCuts = cellCuts[celli];
+
+        const List<Pair<labelPair>>& ceiToCfiAndFei =
+            cAddrs[celli].ceiToCfiAndFei();
+
+        forAll(cCuts, cuti)
+        {
+            if (cCuts[cuti].size() < 3) continue;
+
+            facesDyn.append(face(cCuts[cuti].size()));
+
+            forAll(cCuts[cuti], i)
+            {
+                const label cei = cCuts[cuti][i];
+                const label cfi = ceiToCfiAndFei[cei][0][0];
+                const label fei = ceiToCfiAndFei[cei][0][1];
+
+                const edge e = mesh.faces()[c[cfi]].faceEdge(fei);
+
+                EdgeMap<label>::iterator iter = meshEdgePoint.find(e);
+
+                const label surfacePointi =
+                    iter != meshEdgePoint.end() ? *iter : pointsDyn.size();
+
+                if (surfacePointi == pointsDyn.size())
+                {
+                    meshEdgePoint.insert(e, surfacePointi);
+
+                    const scalar lambda =
+                        cutPoly::edgeCutLambda(e, pAlphas, isoAlpha);
+
+                    pointsDyn.append
+                    (
+                        cutPoly::edgeCutValue(e, lambda, mesh.points())
+                    );
+                    pointEdgesDyn.append(e);
+                    pointEdgeLambdasDyn.append(lambda);
+                }
+
+                facesDyn.last()[i] = surfacePointi;
+            }
+
+            faceCellsDyn.append(celli);
+        }
+    }
+
+    // Transfer to non-dynamic storage
+    points_.transfer(pointsDyn);
+    pointEdges_.transfer(pointEdgesDyn);
+    pointEdgeLambdas_.transfer(pointEdgeLambdasDyn);
+    faces_.transfer(facesDyn);
+    faceCells_.transfer(faceCellsDyn);
+
+    if (debug)
+    {
+        const label nLabels =
+            sum(ListListOps::subSizes(faces(), accessOp<face>()));
+
+        Pout<< typeName << " : constructed surface of size "
+            << points().size()*sizeof(point) + nLabels*sizeof(label)
+            << " in " << cpuTime.cpuTimeIncrement() << "s" << nl << endl;
+    }
+}
+
+
+Foam::cutPolyIsoSurface::cutPolyIsoSurface
+(
+    const PtrList<cutPolyIsoSurface>& isos
+)
+:
+    points_(),
+    pointEdges_(),
+    pointEdgeLambdas_(),
+    faces_(),
+    faceCells_()
+{
+    label nPoints = 0, nFaces = 0;
+    forAll(isos, i)
+    {
+        nPoints += isos[i].points_.size();
+        nFaces += isos[i].faces_.size();
+    }
+
+    points_.resize(nPoints);
+    pointEdges_.resize(nPoints);
+    pointEdgeLambdas_.resize(nPoints);
+    faces_.resize(nFaces);
+    faceCells_.resize(nFaces);
+
+    label pointi0 = 0, facei0 = 0;
+    forAll(isos, i)
+    {
+        const label nPs = isos[i].points_.size();
+        const label nFs = isos[i].faces_.size();
+
+        SubList<point>(points_, nPs, pointi0) = isos[i].points_;
+        SubList<edge>(pointEdges_, nPs, pointi0) = isos[i].pointEdges_;
+        SubList<scalar>(pointEdgeLambdas_, nPs, pointi0) =
+            isos[i].pointEdgeLambdas_;
+        SubList<face>(faces_, nFs, facei0) = isos[i].faces_;
+        SubList<label>(faceCells_, nFs, facei0) = isos[i].faceCells_;
+
+        forAll(isos[i].faces_, fi)
+        {
+            forAll(faces_[facei0 + fi], fpi)
+            {
+                faces_[facei0 + fi][fpi] += pointi0;
+            }
+        }
+
+        pointi0 += nPs;
+        facei0 += nFs;
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::cutPolyIsoSurface::~cutPolyIsoSurface()
+{}
+
+
+// ************************************************************************* //

src/meshTools/cutPoly/cutPolyIsoSurface.H

@@ -0,0 +1,147 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::cutPolyIsoSurface
+
+Description
+
+SourceFiles
+    cutPolyIsoSurfaceI.H
+    cutPolyIsoSurface.C
+    cutPolyIsoSurfaceIO.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef cutPolyIsoSurface_H
+#define cutPolyIsoSurface_H
+
+#include "polyMesh.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                          Class cutPolyIsoSurface Declaration
+\*---------------------------------------------------------------------------*/
+
+class cutPolyIsoSurface
+{
+    // Private Data
+
+        //- Points of the iso-surface
+        pointField points_;
+
+        //- The mesh edge associated with each point
+        edgeList pointEdges_;
+
+        //- The local coordinate of each point within its associated edge
+        scalarList pointEdgeLambdas_;
+
+        //- Faces of the iso-surface
+        faceList faces_;
+
+        //- The mesh cell index associated with each face
+        labelList faceCells_;
+
+
+public:
+
+    //- Runtime type information
+    TypeName("cutPolyIsoSurface");
+
+
+    // Constructors
+
+        //- Construct from a mesh, point values and an iso-value
+        cutPolyIsoSurface
+        (
+            const polyMesh& mesh,
+            const scalarField& pAlphas,
+            const scalar isoAlpha,
+            const labelList& zoneIDs = NullObjectRef<labelList>()
+        );
+
+        //- Construct by combining a list of iso-surfaces
+        cutPolyIsoSurface(const PtrList<cutPolyIsoSurface>& isos);
+
+        //- Disallow default bitwise copy construction
+        cutPolyIsoSurface(const cutPolyIsoSurface&) = delete;
+
+
+    //- Destructor
+    virtual ~cutPolyIsoSurface();
+
+
+    // Member Functions
+
+        // Access
+
+            //- Points of the iso-surface
+            inline const pointField& points() const
+            {
+                return points_;
+            }
+
+            //- Faces of the iso-surface
+            inline const faceList& faces() const
+            {
+                return faces_;
+            }
+
+            //- The mesh cell index associated with each face
+            inline const labelList& faceCells() const
+            {
+                return faceCells_;
+            }
+
+
+        //- Return a field of values for the surface faces, sampled from their
+        //  associated cells
+        template<class Type>
+        tmp<Field<Type>> sample(const Field<Type> cPsis) const;
+
+        //- Return a field of values for the surface points, interpolated
+        //  within their associated edges
+        template<class Type>
+        tmp<Field<Type>> interpolate(const Field<Type>& pPsis) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+    #include "cutPolyIsoSurfaceTemplates.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/meshTools/cutPoly/cutPolyIsoSurfaceTemplates.C

@@ -0,0 +1,65 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "cutPolyIsoSurface.H"
+#include "cutPolyValue.H"
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+template<class Type>
+Foam::tmp<Foam::Field<Type>> Foam::cutPolyIsoSurface::sample
+(
+    const Field<Type> cPsis
+) const
+{
+    return tmp<Field<Type>>(new Field<Type>(cPsis, faceCells_));
+}
+
+
+template<class Type>
+Foam::tmp<Foam::Field<Type>> Foam::cutPolyIsoSurface::interpolate
+(
+    const Field<Type>& pPsis
+) const
+{
+    tmp<Field<Type>> tValues(new Field<Type>(points_.size()));
+    Field<Type>& values = tValues.ref();
+
+    forAll(points_, pointi)
+    {
+        values[pointi] =
+            cutPoly::edgeCutValue
+            (
+                pointEdges_[pointi],
+                pointEdgeLambdas_[pointi],
+                pPsis
+            );
+    }
+
+    return tValues;
+}
+
+
+// ************************************************************************* //

src/meshTools/cutPoly/cutPolyValue.H

@@ -0,0 +1,348 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Namespace
+    Foam::cutPoly
+
+Description
+    Functions and classes for extracting values from cut edges, faces and cells
+
+SourceFiles
+    cutPolyValueTemplates.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef cutPolyValue_H
+#define cutPolyValue_H
+
+#include "cutPoly.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace cutPoly
+{
+
+//- Get the local coordinate within an edge, given end point values and an
+//  iso-value
+inline scalar edgeCutLambda
+(
+    const edge& e,
+    const scalarField& pAlphas,
+    const scalar isoAlpha
+);
+
+//- Linearly interpolate a value from the end points to the cut point of an
+//  edge, given a local coordinate within the edge
+template<class Type>
+Type edgeCutValue
+(
+    const edge& e,
+    const scalar lambda,
+    const Field<Type>& pPsis
+);
+
+//- Linearly interpolate a value from the end points to the cut point of an
+//  edge, given end point values and an iso-value
+template<class Type>
+Type edgeCutValue
+(
+    const edge& e,
+    const scalarField& pAlphas,
+    const scalar isoAlpha,
+    const Field<Type>& pPsis
+);
+
+
+/*---------------------------------------------------------------------------*\
+                       Class FaceCutValues Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class Type>
+class FaceCutValues
+{
+    // Private Data
+
+        //- The face
+        const face& f_;
+
+        //- The face's cuts
+        const List<labelPair>& fCuts_;
+
+        //- The point values that we want to iterate over
+        const Field<Type>& pPsis_;
+
+        //- The point values that define the iso-surface
+        const scalarField& pAlphas_;
+
+        //- The value that defines the iso-surface
+        const scalar isoAlpha_;
+
+        //- Do we want to iterate over the cut below or above the iso surface?
+        const bool below_;
+
+
+public:
+
+    // Public Typedefs
+
+        //- The value type
+        typedef Type value_type;
+
+
+    // Public Classes
+
+        //- Forward iterator
+        class const_iterator
+        {
+            // Private Data
+
+                //- Reference to the face-values
+                const FaceCutValues<Type>& fValues_;
+
+                //- Index of the sub-face and sub-face-point
+                label i_, j_;
+
+                //- Cache of values on cut edges for this sub-face segment
+                mutable Tuple2<label, Pair<Type>> iAndCutPsis_;
+
+
+            // Private Member Functions
+
+                //- Get the values on the cut edges for this sub-face segment
+                const Pair<Type>& cutPsis(const label i) const;
+
+                //- Return the size of the given sub-face
+                label size(const label i) const;
+
+                //- Return the value for the given sub-face and sub-face-point
+                const Type psi(const label i, const label j) const;
+
+
+        public:
+
+            // Constructors
+
+                //- Construct from components
+                const_iterator
+                (
+                    const FaceCutValues<Type>& fValues,
+                    const label i,
+                    const label j
+                );
+
+
+            // Member Functions
+
+                //- Get the next value around the sub-face
+                Type next() const;
+
+
+            // Member Operators
+
+                //- Equality comparison
+                bool operator==(const const_iterator& it) const;
+
+                //- Inequality comparison
+                bool operator!=(const const_iterator& it) const;
+
+                //- Dereference
+                Type operator*() const;
+
+                //- Increment
+                inline const_iterator& operator++();
+
+                //- Increment
+                inline const_iterator operator++(int);
+        };
+
+
+    // Constructors
+
+        //- Construct from components
+        FaceCutValues
+        (
+            const face& f,
+            const List<labelPair>& fCuts,
+            const Field<Type>& pPsis,
+            const scalarField& pAlphas,
+            const scalar isoAlpha,
+            const bool below
+        );
+
+
+    // Member Functions
+
+        //- Get the beginning of the iteration over the values
+        const_iterator begin() const;
+
+        //- Get the end of the iteration over the values
+        const_iterator end() const;
+};
+
+
+/*---------------------------------------------------------------------------*\
+                       Class CellCutValues Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class Type>
+class CellCutValues
+{
+    // Private Data
+
+        //- The cell
+        const cell& c_;
+
+        //- The cell's edge addressing
+        const cellEdgeAddressing& cAddr_;
+
+        //- The cell's cuts
+        const labelListList& cCuts_;
+
+        //- The faces
+        const faceList& fs_;
+
+        //- The point values that we want to iterate over
+        const Field<Type>& pPsis_;
+
+        //- The point values that define the iso-surface
+        const scalarField& pAlphas_;
+
+        //- The value that defines the iso-surface
+        const scalar isoAlpha_;
+
+
+public:
+
+    // Public Typedefs
+
+        //- The value type
+        typedef Type value_type;
+
+
+    // Public Classes
+
+        //- Forward iterator
+        class const_iterator
+        {
+            // Private Data
+
+                //- Reference to the cell-values
+                const CellCutValues<Type>& cValues_;
+
+                //- Index of the loop and loop-point
+                label i_, j_;
+
+
+            // Private Member Functions
+
+                //- The values at this loop-point and the next
+                Pair<Type> psis_;
+
+                //- Return the value for the given loop and loop-point
+                Type psi(const label i, const label j);
+
+
+        public:
+
+            // Constructors
+
+                //- Construct from components
+                const_iterator
+                (
+                    const CellCutValues<Type>& cValues,
+                    const label i,
+                    const label j
+                );
+
+
+            // Member Functions
+
+                //- Get the next value around the loop
+                const Type& next() const;
+
+
+            // Member Operators
+
+                //- Equality comparison
+                bool operator==(const const_iterator& it) const;
+
+                //- Inequality comparison
+                bool operator!=(const const_iterator& it) const;
+
+                //- Dereference
+                const Type& operator*() const;
+
+                //- Increment
+                inline const_iterator& operator++();
+
+                //- Increment
+                inline const_iterator operator++(int);
+        };
+
+
+    // Constructors
+
+        //- Construct from components
+        CellCutValues
+        (
+            const cell& c,
+            const cellEdgeAddressing& cAddr,
+            const labelListList& cCuts,
+            const faceList& fs,
+            const Field<Type>& pPsis,
+            const scalarField& pAlphas,
+            const scalar isoAlpha
+        );
+
+
+    // Member Functions
+
+        //- Get the beginning of the iteration over the values
+        const_iterator begin() const;
+
+        //- Get the end of the iteration over the values
+        const_iterator end() const;
+};
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace cutPoly
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#include "cutPolyValueI.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+    #include "cutPolyValueTemplates.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/meshTools/cutPoly/cutPolyValueI.H

@@ -0,0 +1,41 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "cutPolyValue.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+inline Foam::scalar Foam::cutPoly::edgeCutLambda
+(
+    const edge& e,
+    const scalarField& pAlphas,
+    const scalar isoAlpha
+)
+{
+    return (isoAlpha - pAlphas[e[0]])/(pAlphas[e[1]] - pAlphas[e[0]]);
+}
+
+
+// ************************************************************************* //

src/meshTools/cutPoly/cutPolyValueTemplates.C

@@ -0,0 +1,384 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "cutPolyValue.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+template<class Type>
+Type Foam::cutPoly::edgeCutValue
+(
+    const edge& e,
+    const scalar lambda,
+    const Field<Type>& pPsis
+)
+{
+    return (1 - lambda)*pPsis[e[0]] + lambda*pPsis[e[1]];
+}
+
+
+template<class Type>
+Type Foam::cutPoly::edgeCutValue
+(
+    const edge& e,
+    const scalarField& pAlphas,
+    const scalar isoAlpha,
+    const Field<Type>& pPsis
+)
+{
+    return edgeCutValue(e, edgeCutLambda(e, pAlphas, isoAlpha), pPsis);
+}
+
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+template<class Type>
+const Foam::Pair<Type>&
+Foam::cutPoly::FaceCutValues<Type>::const_iterator::cutPsis(const label i) const
+{
+    if (i != iAndCutPsis_.first())
+    {
+        const face& f = fValues_.f_;
+        const labelPair& fCut = fValues_.fCuts_[i];
+
+        iAndCutPsis_.first() = i;
+
+        forAll(iAndCutPsis_.second(), j)
+        {
+            iAndCutPsis_.second()[j] =
+                edgeCutValue
+                (
+                    f.faceEdge(fCut[j]),
+                    fValues_.pAlphas_,
+                    fValues_.isoAlpha_,
+                    fValues_.pPsis_
+                );
+        }
+    }
+
+    return iAndCutPsis_.second();
+}
+
+
+template<class Type>
+Foam::label
+Foam::cutPoly::FaceCutValues<Type>::const_iterator::size(const label i) const
+{
+    const face& f = fValues_.f_;
+    const List<labelPair>& fCuts = fValues_.fCuts_;
+    const labelPair& fCut0 = fCuts[i];
+    const labelPair& fCut1 = fCuts[(i + 1) % fCuts.size()];
+
+    const label dfCut =
+        fValues_.below_ == separatedBelow
+      ? fCut0[1] - fCut0[0]
+      : fCut1[0] - fCut0[1];
+
+    return 2 + ((dfCut + f.size()) % f.size());
+}
+
+
+template<class Type>
+const Type Foam::cutPoly::FaceCutValues<Type>::const_iterator::psi
+(
+    const label i,
+    const label j
+) const
+{
+    const face& f = fValues_.f_;
+    const labelPair& fCut = fValues_.fCuts_[i];
+
+    if (j < 2)
+    {
+        const label fCutj = fValues_.below_ == separatedBelow ? 1 - j : j;
+        return cutPsis(i)[fCutj];
+    }
+    else
+    {
+        const label fCutj = fValues_.below_ == separatedBelow ? 0 : 1;
+        return fValues_.pPsis_[f[(j - 2 + fCut[fCutj] + 1) % f.size()]];
+    }
+}
+
+
+template<class Type>
+Type Foam::cutPoly::CellCutValues<Type>::const_iterator::psi
+(
+    const label i,
+    const label j
+)
+{
+    const label iStar = i % cValues_.cCuts_.size();
+    const label jStar = j % cValues_.cCuts_[iStar].size();
+
+    const label cei = cValues_.cCuts_[iStar][jStar];
+    const label cfi = cValues_.cAddr_.ceiToCfiAndFei()[cei][0][0];
+    const label fei = cValues_.cAddr_.ceiToCfiAndFei()[cei][0][1];
+
+    return
+        edgeCutValue
+        (
+            cValues_.fs_[cValues_.c_[cfi]].faceEdge(fei),
+            cValues_.pAlphas_,
+            cValues_.isoAlpha_,
+            cValues_.pPsis_
+        );
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class Type>
+Foam::cutPoly::FaceCutValues<Type>::const_iterator::const_iterator
+(
+    const FaceCutValues<Type>& fValues,
+    const label i,
+    const label j
+)
+:
+    fValues_(fValues),
+    i_(i),
+    j_(j),
+    iAndCutPsis_(-1, Pair<Type>(Zero, Zero))
+{}
+
+
+template<class Type>
+Foam::cutPoly::FaceCutValues<Type>::FaceCutValues
+(
+    const face& f,
+    const List<labelPair>& fCuts,
+    const Field<Type>& pPsis,
+    const scalarField& pAlphas,
+    const scalar isoAlpha,
+    const bool below
+)
+:
+    f_(f),
+    fCuts_(fCuts),
+    pPsis_(pPsis),
+    pAlphas_(pAlphas),
+    isoAlpha_(isoAlpha),
+    below_(below)
+{}
+
+
+template<class Type>
+Foam::cutPoly::CellCutValues<Type>::const_iterator::const_iterator
+(
+    const CellCutValues<Type>& cValues,
+    const label i,
+    const label j
+)
+:
+    cValues_(cValues),
+    i_(i),
+    j_(j),
+    psis_(psi(i, j), psi(i, j + 1))
+{}
+
+
+template<class Type>
+Foam::cutPoly::CellCutValues<Type>::CellCutValues
+(
+    const cell& c,
+    const cellEdgeAddressing& cAddr,
+    const labelListList& cCuts,
+    const faceList& fs,
+    const Field<Type>& pPsis,
+    const scalarField& pAlphas,
+    const scalar isoAlpha
+)
+:
+    c_(c),
+    cAddr_(cAddr),
+    cCuts_(cCuts),
+    fs_(fs),
+    pPsis_(pPsis),
+    pAlphas_(pAlphas),
+    isoAlpha_(isoAlpha)
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+template<class Type>
+Type Foam::cutPoly::FaceCutValues<Type>::const_iterator::next() const
+{
+    const label j = (j_ + 1) % size(i_);
+    const label i = i_ + (fValues_.below_ != separatedBelow && j == 0);
+    return i == fValues_.fCuts_.size() ? psi(0, 0) : psi(i, j);
+}
+
+
+template<class Type>
+typename Foam::cutPoly::FaceCutValues<Type>::const_iterator
+Foam::cutPoly::FaceCutValues<Type>::begin() const
+{
+    return const_iterator(*this, 0, 0);
+}
+
+
+template<class Type>
+typename Foam::cutPoly::FaceCutValues<Type>::const_iterator
+Foam::cutPoly::FaceCutValues<Type>::end() const
+{
+    return const_iterator(*this, this->fCuts_.size(), 0);
+}
+
+
+template<class Type>
+const Type& Foam::cutPoly::CellCutValues<Type>::const_iterator::next() const
+{
+    return psis_.second();
+}
+
+
+template<class Type>
+typename Foam::cutPoly::CellCutValues<Type>::const_iterator
+Foam::cutPoly::CellCutValues<Type>::begin() const
+{
+    return const_iterator(*this, 0, 0);
+}
+
+
+template<class Type>
+typename Foam::cutPoly::CellCutValues<Type>::const_iterator
+Foam::cutPoly::CellCutValues<Type>::end() const
+{
+    return const_iterator(*this, this->cCuts_.size(), 0);
+}
+
+
+// * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * * //
+
+template<class Type>
+bool Foam::cutPoly::FaceCutValues<Type>::const_iterator::operator==
+(
+    const const_iterator& it
+) const
+{
+    return it.i_ == i_ && it.j_ == j_;
+}
+
+
+template<class Type>
+bool Foam::cutPoly::FaceCutValues<Type>::const_iterator::operator!=
+(
+    const const_iterator& it
+) const
+{
+    return !(it == *this);
+}
+
+
+template<class Type>
+Type Foam::cutPoly::FaceCutValues<Type>::const_iterator::operator*() const
+{
+    return psi(i_, j_);
+}
+
+
+template<class Type>
+inline typename Foam::cutPoly::FaceCutValues<Type>::const_iterator&
+Foam::cutPoly::FaceCutValues<Type>::const_iterator::operator++()
+{
+    j_ = (j_ + 1) % size(i_);
+    i_ += j_ == 0;
+    return *this;
+}
+
+
+template<class Type>
+inline typename Foam::cutPoly::FaceCutValues<Type>::const_iterator
+Foam::cutPoly::FaceCutValues<Type>::const_iterator::operator++(int)
+{
+    const const_iterator it(*this);
+    ++ *this;
+    return it;
+}
+
+
+template<class Type>
+bool Foam::cutPoly::CellCutValues<Type>::const_iterator::operator==
+(
+    const const_iterator& it
+) const
+{
+    return it.i_ == i_ && it.j_ == j_;
+}
+
+template<class Type>
+bool Foam::cutPoly::CellCutValues<Type>::const_iterator::operator!=
+(
+    const const_iterator& it
+) const
+{
+    return !(it == *this);
+}
+
+
+template<class Type>
+const Type&
+Foam::cutPoly::CellCutValues<Type>::const_iterator::operator*() const
+{
+    return psis_.first();
+}
+
+
+template<class Type>
+inline typename Foam::cutPoly::CellCutValues<Type>::const_iterator&
+Foam::cutPoly::CellCutValues<Type>::const_iterator::operator++()
+{
+    ++ j_;
+
+    if (j_ == cValues_.cCuts_[i_].size())
+    {
+        j_ = 0;
+        ++ i_;
+        psis_.first() = psi(i_, j_);
+    }
+    else
+    {
+        psis_.first() = psis_.second();
+    }
+
+    psis_.second() = psi(i_, j_ + 1);
+
+    return *this;
+}
+
+
+template<class Type>
+inline typename Foam::cutPoly::CellCutValues<Type>::const_iterator
+Foam::cutPoly::CellCutValues<Type>::const_iterator::operator++(int)
+{
+    const const_iterator it(*this);
+    ++ *this;
+    return it;
+}
+
+
+// ************************************************************************* //

src/meshTools/cutTriTet/cutTriTet.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2017-2020 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2017-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -26,13 +26,13 @@ Description:
     applied to each half of a cut.
 
 SourceFiles:
-    cutI.H
-    cutTemplates.C
+    cutTriTetI.H
+    cutTriTetTemplates.C
 
 \*---------------------------------------------------------------------------*/
 
-#ifndef cut_H
-#define cut_H
+#ifndef cutTriTet_H
+#define cutTriTet_H
 
 #include "FixedList.H"
 #include "nil.H"
@@ -45,7 +45,7 @@ SourceFiles:
 
 namespace Foam
 {
-namespace cut
+namespace cutTriTet
 {
 
 /*---------------------------------------------------------------------------*\
@@ -560,7 +560,7 @@ public:
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-} // End namespace cut
+} // End namespace cutTriTet
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -568,7 +568,7 @@ public:
 //  on aboveOp and belowOp; the operations applied to the regions on either side
 //  of the cut.
 template<class Point, class AboveOp, class BelowOp>
-typename cut::opAddResult<AboveOp, BelowOp>::type triCut
+typename cutTriTet::opAddResult<AboveOp, BelowOp>::type triCut
 (
     const FixedList<Point, 3>& tri,
     const FixedList<scalar, 3>& level,
@@ -578,7 +578,7 @@ typename cut::opAddResult<AboveOp, BelowOp>::type triCut
 
 //- As above, but with a plane specifying the location of the cut
 template<class AboveOp, class BelowOp>
-typename cut::opAddResult<AboveOp, BelowOp>::type triCut
+typename cutTriTet::opAddResult<AboveOp, BelowOp>::type triCut
 (
     const FixedList<point, 3>& tri,
     const plane& s,
@@ -588,7 +588,7 @@ typename cut::opAddResult<AboveOp, BelowOp>::type triCut
 
 //- As triCut, but for a tetrahedron.
 template<class Point, class AboveOp, class BelowOp>
-typename cut::opAddResult<AboveOp, BelowOp>::type tetCut
+typename cutTriTet::opAddResult<AboveOp, BelowOp>::type tetCut
 (
     const FixedList<Point, 4>& tet,
     const FixedList<scalar, 4>& level,
@@ -598,7 +598,7 @@ typename cut::opAddResult<AboveOp, BelowOp>::type tetCut
 
 //- As above, but with a plane specifying the location of the cut
 template<class AboveOp, class BelowOp>
-typename cut::opAddResult<AboveOp, BelowOp>::type tetCut
+typename cutTriTet::opAddResult<AboveOp, BelowOp>::type tetCut
 (
     const FixedList<point, 4>& tet,
     const plane& s,
@@ -612,10 +612,10 @@ typename cut::opAddResult<AboveOp, BelowOp>::type tetCut
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-#include "cutI.H"
+#include "cutTriTetI.H"
 
 #ifdef NoRepository
-    #include "cutTemplates.C"
+    #include "cutTriTetTemplates.C"
 #endif
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

src/meshTools/cutTriTet/cutTriTetI.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2017-2019 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2017-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -23,7 +23,7 @@ License
 
 \*---------------------------------------------------------------------------*/
 
-#include "cut.H"
+#include "cutTriTet.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -34,9 +34,9 @@ namespace Foam
 
 //- Modify a uniform operation for reordering a tri (does nothing)
 template<class Type>
-inline const cut::uniformOp<Type>& triReorder
+inline const cutTriTet::uniformOp<Type>& triReorder
 (
-    const cut::uniformOp<Type>& x,
+    const cutTriTet::uniformOp<Type>& x,
     const FixedList<label, 3>&
 )
 {
@@ -46,9 +46,9 @@ inline const cut::uniformOp<Type>& triReorder
 
 //- Modify a uniform operation for cutting a tri from a tri (does nothing)
 template<class Type>
-inline const cut::uniformOp<Type>& triCutTri
+inline const cutTriTet::uniformOp<Type>& triCutTri
 (
-    const cut::uniformOp<Type>& x,
+    const cutTriTet::uniformOp<Type>& x,
     const Pair<scalar>&
 )
 {
@@ -58,9 +58,9 @@ inline const cut::uniformOp<Type>& triCutTri
 
 //- Modify a uniform operation for cutting a quad from a tri (does nothing)
 template<class Type>
-inline const cut::uniformOp<Type>& triCutQuad
+inline const cutTriTet::uniformOp<Type>& triCutQuad
 (
-    const cut::uniformOp<Type>& x,
+    const cutTriTet::uniformOp<Type>& x,
     const Pair<scalar>&
 )
 {
@@ -70,9 +70,9 @@ inline const cut::uniformOp<Type>& triCutQuad
 
 //- Modify a uniform operation for reordering a tet (does nothing)
 template<class Type>
-inline const cut::uniformOp<Type>& tetReorder
+inline const cutTriTet::uniformOp<Type>& tetReorder
 (
-    const cut::uniformOp<Type>& x,
+    const cutTriTet::uniformOp<Type>& x,
     const FixedList<label, 4>&
 )
 {
@@ -82,9 +82,9 @@ inline const cut::uniformOp<Type>& tetReorder
 
 //- Modify a uniform operation for cutting a tet from a tet (does nothing)
 template<class Type>
-inline const cut::uniformOp<Type>& tetCutTet
+inline const cutTriTet::uniformOp<Type>& tetCutTet
 (
-    const cut::uniformOp<Type>& x,
+    const cutTriTet::uniformOp<Type>& x,
     const FixedList<scalar, 3>&
 )
 {
@@ -94,9 +94,9 @@ inline const cut::uniformOp<Type>& tetCutTet
 
 //- Modify a uniform operation for cutting prism0 from a tet (does nothing)
 template<class Type>
-inline const cut::uniformOp<Type>& tetCutPrism0
+inline const cutTriTet::uniformOp<Type>& tetCutPrism0
 (
-    const cut::uniformOp<Type>& x,
+    const cutTriTet::uniformOp<Type>& x,
     const FixedList<scalar, 3>&
 )
 {
@@ -106,9 +106,9 @@ inline const cut::uniformOp<Type>& tetCutPrism0
 
 //- Modify a uniform operation for cutting prism01 from a tet (does nothing)
 template<class Type>
-inline const cut::uniformOp<Type>& tetCutPrism01
+inline const cutTriTet::uniformOp<Type>& tetCutPrism01
 (
-    const cut::uniformOp<Type>& x,
+    const cutTriTet::uniformOp<Type>& x,
     const FixedList<scalar, 4>&
 )
 {
@@ -118,9 +118,9 @@ inline const cut::uniformOp<Type>& tetCutPrism01
 
 //- Modify a uniform operation for cutting prism23 from a tet (does nothing)
 template<class Type>
-inline const cut::uniformOp<Type>& tetCutPrism23
+inline const cutTriTet::uniformOp<Type>& tetCutPrism23
 (
-    const cut::uniformOp<Type>& x,
+    const cutTriTet::uniformOp<Type>& x,
     const FixedList<scalar, 4>&
 )
 {

src/meshTools/cutTriTet/cutTriTetTemplates.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2017-2019 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2017-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -23,12 +23,12 @@ License
 
 \*---------------------------------------------------------------------------*/
 
-#include "cut.H"
+#include "cutTriTet.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 template<class Point, class AboveOp, class BelowOp>
-typename Foam::cut::opAddResult<AboveOp, BelowOp>::type Foam::triCut
+typename Foam::cutTriTet::opAddResult<AboveOp, BelowOp>::type Foam::triCut
 (
     const FixedList<Point, 3>& tri,
     const FixedList<scalar, 3>& level,
@@ -99,7 +99,7 @@ typename Foam::cut::opAddResult<AboveOp, BelowOp>::type Foam::triCut
 
 
 template<class AboveOp, class BelowOp>
-typename Foam::cut::opAddResult<AboveOp, BelowOp>::type Foam::triCut
+typename Foam::cutTriTet::opAddResult<AboveOp, BelowOp>::type Foam::triCut
 (
     const FixedList<point, 3>& tri,
     const plane& p,
@@ -120,7 +120,7 @@ typename Foam::cut::opAddResult<AboveOp, BelowOp>::type Foam::triCut
 
 
 template<class Point, class AboveOp, class BelowOp>
-typename Foam::cut::opAddResult<AboveOp, BelowOp>::type Foam::tetCut
+typename Foam::cutTriTet::opAddResult<AboveOp, BelowOp>::type Foam::tetCut
 (
     const FixedList<Point, 4>& tet,
     const FixedList<scalar, 4>& level,
@@ -248,7 +248,7 @@ typename Foam::cut::opAddResult<AboveOp, BelowOp>::type Foam::tetCut
 
 
 template<class AboveOp, class BelowOp>
-typename Foam::cut::opAddResult<AboveOp, BelowOp>::type Foam::tetCut
+typename Foam::cutTriTet::opAddResult<AboveOp, BelowOp>::type Foam::tetCut
 (
     const FixedList<point, 4>& tet,
     const plane& p,

src/meshTools/tetOverlapVolume/tetOverlapVolume.C

@@ -29,7 +29,7 @@ License
 #include "treeBoundBox.H"
 #include "indexedOctree.H"
 #include "treeDataCell.H"
-#include "cut.H"
+#include "cutTriTet.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -53,6 +53,8 @@ Foam::scalar Foam::tetOverlapVolume::tetTetOverlapVol
     const tetPointRef& tetB
 ) const
 {
+    typedef cutTriTet::noOp noOp;
+
     // A maximum of three cuts are made (the tets that result from the final cut
     // are not stored), and each cut can create at most three tets. The
     // temporary storage must therefore extend to 3^3 = 27 tets.
@@ -67,7 +69,7 @@ Foam::scalar Foam::tetOverlapVolume::tetTetOverlapVol
     }
     const FixedList<point, 4> t({tetA.a(), tetA.b(), tetA.c(), tetA.d()});
     cutTetList1.clear();
-    tetCut(t, pl0, cut::appendOp<tetListType>(cutTetList1), cut::noOp());
+    tetCut(t, pl0, cutTriTet::appendOp<tetListType>(cutTetList1), noOp());
     if (cutTetList1.size() == 0)
     {
         return 0;
@@ -83,7 +85,7 @@ Foam::scalar Foam::tetOverlapVolume::tetTetOverlapVol
     for (label i = 0; i < cutTetList1.size(); i++)
     {
         const FixedList<point, 4>& t = cutTetList1[i];
-        tetCut(t, pl1, cut::appendOp<tetListType>(cutTetList2), cut::noOp());
+        tetCut(t, pl1, cutTriTet::appendOp<tetListType>(cutTetList2), noOp());
     }
     if (cutTetList2.size() == 0)
     {
@@ -100,7 +102,7 @@ Foam::scalar Foam::tetOverlapVolume::tetTetOverlapVol
     for (label i = 0; i < cutTetList2.size(); i++)
     {
         const FixedList<point, 4>& t = cutTetList2[i];
-        tetCut(t, pl2, cut::appendOp<tetListType>(cutTetList1), cut::noOp());
+        tetCut(t, pl2, cutTriTet::appendOp<tetListType>(cutTetList1), noOp());
     }
     if (cutTetList1.size() == 0)
     {
@@ -117,7 +119,7 @@ Foam::scalar Foam::tetOverlapVolume::tetTetOverlapVol
     for (label i = 0; i < cutTetList1.size(); i++)
     {
         const FixedList<point, 4>& t = cutTetList1[i];
-        v += tetCut(t, pl3, cut::volumeOp(), cut::noOp());
+        v += tetCut(t, pl3, cutTriTet::volumeOp(), noOp());
     }
 
     return v;

src/sampling/Make/files

@@ -34,20 +34,17 @@ $(setWriters)/csv/csvSetWriter.C
 $(setWriters)/gnuplot/gnuplotSetWriter.C
 $(setWriters)/none/noSetWriter.C
 
-cuttingPlane/cuttingPlane.C
-
 sampledSurface/sampledPatch/sampledPatch.C
 sampledSurface/sampledPatchInternalField/sampledPatchInternalField.C
-sampledSurface/sampledPlane/sampledPlane.C
-sampledSurface/isoSurface/isoSurface.C
-sampledSurface/isoSurface/sampledIsoSurface.C
-sampledSurface/distanceSurface/distanceSurface.C
-sampledSurface/sampledCuttingPlane/sampledCuttingPlane.C
+sampledSurface/sampledIsoSurfaceSurface/sampledIsoSurfaceSurface.C
+sampledSurface/sampledIsoSurface/sampledIsoSurface.C
+sampledSurface/sampledDistanceSurface/sampledDistanceSurface.C
+sampledSurface/sampledCutPlane/sampledCutPlane.C
 sampledSurface/sampledSurface/sampledSurface.C
 sampledSurface/sampledSurfaces/sampledSurfaces.C
 sampledSurface/sampledTriSurfaceMesh/sampledTriSurfaceMesh.C
-sampledSurface/thresholdCellFaces/thresholdCellFaces.C
-sampledSurface/thresholdCellFaces/sampledThresholdCellFaces.C
+sampledSurface/sampledThresholdCellFaces/thresholdCellFaces.C
+sampledSurface/sampledThresholdCellFaces/sampledThresholdCellFaces.C
 
 surfWriters = sampledSurface/writers
 

src/sampling/cuttingPlane/cuttingPlane.C

@@ -1,430 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2021 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "cuttingPlane.H"
-#include "primitiveMesh.H"
-#include "linePointRef.H"
-#include "meshTools.H"
-#include "polygonTriangulate.H"
-
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-// Set values for what is close to zero and what is considered to
-// be positive (and not just rounding noise)
-//! \cond localScope
-const Foam::scalar zeroish  = Foam::small;
-const Foam::scalar positive = Foam::small * 1E3;
-//! \endcond
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-void Foam::cuttingPlane::calcCutCells
-(
-    const primitiveMesh& mesh,
-    const scalarField& dotProducts,
-    const labelUList& cellIdLabels
-)
-{
-    const labelListList& cellEdges = mesh.cellEdges();
-    const edgeList& edges = mesh.edges();
-
-    label listSize = cellEdges.size();
-    if (notNull(cellIdLabels))
-    {
-        listSize = cellIdLabels.size();
-    }
-
-    cutCells_.setSize(listSize);
-    label cutcelli(0);
-
-    // fp: check if the list (cellZone) is not empty.
-    const bool isZoneEmpty
-    (
-        (returnReduce(cellIdLabels.size(), sumOp<label>()) > 0) ? false : true
-    );
-
-    // Find the cut cells by detecting any cell that uses points with
-    // opposing dotProducts.
-    for (label listI = 0; listI < listSize; ++listI)
-    {
-        label celli = listI;
-
-        if (notNull(cellIdLabels))
-        {
-            celli = cellIdLabels[listI];
-        }
-        else
-        {
-            // fp: in parallel computation, if the cellZone exists globally
-            // but not locally, the postprocessing must be still be limited to
-            // the crossing plane.
-            if (!isZoneEmpty)
-            {
-                cutCells_.setSize(0);
-                return;
-            }
-        }
-        const labelList& cEdges = cellEdges[celli];
-
-        label nCutEdges = 0;
-
-        forAll(cEdges, i)
-        {
-            const edge& e = edges[cEdges[i]];
-
-            if
-            (
-                (dotProducts[e[0]] < zeroish && dotProducts[e[1]] > positive)
-             || (dotProducts[e[1]] < zeroish && dotProducts[e[0]] > positive)
-            )
-            {
-                nCutEdges++;
-
-                if (nCutEdges > 2)
-                {
-                    cutCells_[cutcelli++] = celli;
-
-                    break;
-                }
-            }
-        }
-    }
-
-    // Set correct list size
-    cutCells_.setSize(cutcelli);
-}
-
-
-void Foam::cuttingPlane::intersectEdges
-(
-    const primitiveMesh& mesh,
-    const scalarField& dotProducts,
-    List<label>& edgePoint
-)
-{
-    // Use the dotProducts to find out the cut edges.
-    const edgeList& edges = mesh.edges();
-    const pointField& points = mesh.points();
-
-    // Per edge -1 or the label of the intersection point
-    edgePoint.setSize(edges.size());
-
-    DynamicList<point> dynCuttingPoints(4*cutCells_.size());
-
-    forAll(edges, edgeI)
-    {
-        const edge& e = edges[edgeI];
-
-        if
-        (
-            (dotProducts[e[0]] < zeroish && dotProducts[e[1]] > positive)
-         || (dotProducts[e[1]] < zeroish && dotProducts[e[0]] > positive)
-        )
-        {
-            // Edge is cut
-            edgePoint[edgeI] = dynCuttingPoints.size();
-
-            const point& p0 = points[e[0]];
-            const point& p1 = points[e[1]];
-
-            scalar alpha = lineIntersect(linePointRef(p0, p1));
-
-            if (alpha < zeroish)
-            {
-                dynCuttingPoints.append(p0);
-            }
-            else if (alpha >= 1.0)
-            {
-                dynCuttingPoints.append(p1);
-            }
-            else
-            {
-                dynCuttingPoints.append((1-alpha)*p0 + alpha*p1);
-            }
-        }
-        else
-        {
-            edgePoint[edgeI] = -1;
-        }
-    }
-
-    this->storedPoints().transfer(dynCuttingPoints);
-}
-
-
-bool Foam::cuttingPlane::walkCell
-(
-    const primitiveMesh& mesh,
-    const labelUList& edgePoint,
-    const label celli,
-    const label startEdgeI,
-    DynamicList<label>& faceVerts
-)
-{
-    label facei = -1;
-    label edgeI = startEdgeI;
-
-    label nIter = 0;
-
-    faceVerts.clear();
-    do
-    {
-        faceVerts.append(edgePoint[edgeI]);
-
-        // Cross edge to other face
-        facei = meshTools::otherFace(mesh, celli, facei, edgeI);
-
-        // Find next cut edge on face.
-        const labelList& fEdges = mesh.faceEdges()[facei];
-
-        label nextEdgeI = -1;
-
-        // Note: here is where we should check for whether there are more
-        // than 2 intersections with the face (warped/non-convex face).
-        // If so should e.g. decompose the cells on both faces and redo
-        // the calculation.
-
-        forAll(fEdges, i)
-        {
-            label edge2I = fEdges[i];
-
-            if (edge2I != edgeI && edgePoint[edge2I] != -1)
-            {
-                nextEdgeI = edge2I;
-                break;
-            }
-        }
-
-        if (nextEdgeI == -1)
-        {
-            // Did not find another cut edge on facei. Do what?
-            WarningInFunction
-                << "Did not find closed walk along surface of cell " << celli
-                << " starting from edge " << startEdgeI
-                << " in " << nIter << " iterations." << nl
-                << "Collected cutPoints so far:" << faceVerts
-                << endl;
-
-            return false;
-        }
-
-        edgeI = nextEdgeI;
-
-        nIter++;
-
-        if (nIter > 1000)
-        {
-            WarningInFunction
-                << "Did not find closed walk along surface of cell " << celli
-                << " starting from edge " << startEdgeI
-                << " in " << nIter << " iterations." << nl
-                << "Collected cutPoints so far:" << faceVerts
-                << endl;
-            return false;
-        }
-
-    } while (edgeI != startEdgeI);
-
-
-    if (faceVerts.size() >= 3)
-    {
-        return true;
-    }
-    else
-    {
-        WarningInFunction
-            << "Did not find closed walk along surface of cell " << celli
-            << " starting from edge " << startEdgeI << nl
-            << "Collected cutPoints so far:" << faceVerts
-            << endl;
-
-        return false;
-    }
-}
-
-
-void Foam::cuttingPlane::walkCellCuts
-(
-    const primitiveMesh& mesh,
-    const bool triangulate,
-    const labelUList& edgePoint
-)
-{
-    const pointField& cutPoints = this->points();
-
-    // use dynamic lists to handle triangulation and/or missed cuts
-    DynamicList<face>  dynCutFaces(cutCells_.size());
-    DynamicList<label> dynCutCells(cutCells_.size());
-
-    // scratch space for calculating the face vertices
-    DynamicList<label> faceVerts(10);
-
-    // Create a triangulation engine
-    polygonTriangulate triEngine;
-
-    forAll(cutCells_, i)
-    {
-        label celli = cutCells_[i];
-
-        // Find the starting edge to walk from.
-        const labelList& cEdges = mesh.cellEdges()[celli];
-
-        label startEdgeI = -1;
-
-        forAll(cEdges, cEdgeI)
-        {
-            label edgeI = cEdges[cEdgeI];
-
-            if (edgePoint[edgeI] != -1)
-            {
-                startEdgeI = edgeI;
-                break;
-            }
-        }
-
-        // Check for the unexpected ...
-        if (startEdgeI == -1)
-        {
-            FatalErrorInFunction
-                << "Cannot find cut edge for cut cell " << celli
-                << abort(FatalError);
-        }
-
-        // Walk from starting edge around the circumference of the cell.
-        bool okCut = walkCell
-        (
-            mesh,
-            edgePoint,
-            celli,
-            startEdgeI,
-            faceVerts
-        );
-
-        if (okCut)
-        {
-            face f(faceVerts);
-
-            // Orient face to point in the same direction as the plane normal
-            if ((f.area(cutPoints) & normal()) < 0)
-            {
-                f.flip();
-            }
-
-            // the cut faces are usually quite ugly, so optionally triangulate
-            if (triangulate)
-            {
-                triEngine.triangulate
-                (
-                    UIndirectList<point>(cutPoints, f)
-                );
-                forAll(triEngine.triPoints(), i)
-                {
-                    dynCutFaces.append(triEngine.triPoints(i, f));
-                    dynCutCells.append(celli);
-                }
-            }
-            else
-            {
-                dynCutFaces.append(f);
-                dynCutCells.append(celli);
-            }
-        }
-    }
-
-    this->storedFaces().transfer(dynCutFaces);
-    cutCells_.transfer(dynCutCells);
-}
-
-
-// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
-
-void Foam::cuttingPlane::reCut
-(
-    const primitiveMesh& mesh,
-    const bool triangulate,
-    const labelUList& cellIdLabels
-)
-{
-    MeshStorage::clear();
-    cutCells_.clear();
-
-    const scalarField dotProducts((mesh.points() - refPoint()) & normal());
-
-    // Determine cells that are (probably) cut.
-    calcCutCells(mesh, dotProducts, cellIdLabels);
-
-    // Determine cutPoints and return list of edge cuts.
-    // per edge -1 or the label of the intersection point
-    labelList edgePoint;
-    intersectEdges(mesh, dotProducts, edgePoint);
-
-    // Do topological walk around cell to find closed loop.
-    walkCellCuts(mesh, triangulate, edgePoint);
-}
-
-
-void Foam::cuttingPlane::remapFaces
-(
-    const labelUList& faceMap
-)
-{
-    // recalculate the cells cut
-    if (notNull(faceMap) && faceMap.size())
-    {
-        MeshStorage::remapFaces(faceMap);
-
-        List<label> newCutCells(faceMap.size());
-        forAll(faceMap, facei)
-        {
-            newCutCells[facei] = cutCells_[faceMap[facei]];
-        }
-        cutCells_.transfer(newCutCells);
-    }
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::cuttingPlane::cuttingPlane(const plane& pln)
-:
-    plane(pln)
-{}
-
-
-Foam::cuttingPlane::cuttingPlane
-(
-    const plane& pln,
-    const primitiveMesh& mesh,
-    const bool triangulate,
-    const labelUList& cellIdLabels
-)
-:
-    plane(pln)
-{
-    reCut(mesh, triangulate, cellIdLabels);
-}
-
-
-// ************************************************************************* //

src/sampling/cuttingPlane/cuttingPlane.H

@@ -1,195 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-Class
-    Foam::cuttingPlane
-
-Description
-    Constructs plane through mesh.
-
-    No attempt at resolving degenerate cases. Since the cut faces are
-    usually quite ugly, they will always be triangulated.
-
-    Note:
-        When the cutting plane coincides with a mesh face, the cell edge on the
-        positive side of the plane is taken.
-
-SourceFiles
-    cuttingPlane.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef cuttingPlane_H
-#define cuttingPlane_H
-
-#include "plane.H"
-#include "pointField.H"
-#include "faceList.H"
-#include "MeshedSurface.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-class primitiveMesh;
-
-/*---------------------------------------------------------------------------*\
-                        Class cuttingPlane Declaration
-\*---------------------------------------------------------------------------*/
-
-class cuttingPlane
-:
-    public plane,
-    public MeshedSurface<face>
-{
-    // Private Typedef
-
-        typedef MeshedSurface<face> MeshStorage;
-
-
-    // Private Data
-
-        //- List of cells cut by the plane
-        labelList cutCells_;
-
-
-    // Private Member Functions
-
-        //- Determine cut cells, possibly restricted to a list of cells
-        void calcCutCells
-        (
-            const primitiveMesh&,
-            const scalarField& dotProducts,
-            const labelUList& cellIdLabels = labelUList::null()
-        );
-
-        //- Determine intersection points (cutPoints).
-        void intersectEdges
-        (
-            const primitiveMesh&,
-            const scalarField& dotProducts,
-            List<label>& edgePoint
-        );
-
-        //- Walk circumference of cell, starting from startEdgeI crossing
-        //  only cut edges. Record cutPoint labels in faceVerts.
-        static bool walkCell
-        (
-            const primitiveMesh&,
-            const labelUList& edgePoint,
-            const label celli,
-            const label startEdgeI,
-            DynamicList<label>& faceVerts
-        );
-
-        //- Determine cuts for all cut cells.
-        void walkCellCuts
-        (
-            const primitiveMesh& mesh,
-            const bool triangulate,
-            const labelUList& edgePoint
-        );
-
-
-protected:
-
-    // Constructors
-
-        //- Construct plane description without cutting
-        cuttingPlane(const plane&);
-
-
-    // Protected Member Functions
-
-        //- Recut mesh with existing planeDesc, restricted to a list of cells
-        void reCut
-        (
-            const primitiveMesh&,
-            const bool triangulate,
-            const labelUList& cellIdLabels = labelUList::null()
-        );
-
-        //- Remap action on triangulation or cleanup
-        virtual void remapFaces(const labelUList& faceMap);
-
-
-public:
-
-    // Constructors
-
-        //- Construct from plane and mesh reference,
-        //  possibly restricted to a list of cells
-        cuttingPlane
-        (
-            const plane&,
-            const primitiveMesh&,
-            const bool triangulate,
-            const labelUList& cellIdLabels = labelUList::null()
-        );
-
-
-    // Member Functions
-
-        //- Return plane used
-        const plane& planeDesc() const
-        {
-            return static_cast<const plane&>(*this);
-        }
-
-        //- Return List of cells cut by the plane
-        const labelList& cutCells() const
-        {
-            return cutCells_;
-        }
-
-        //- Return true or false to question: have any cells been cut?
-        bool cut() const
-        {
-            return cutCells_.size();
-        }
-
-        //- Sample the cell field
-        template<class Type>
-        tmp<Field<Type>> sample(const Field<Type>&) const;
-
-        template<class Type>
-        tmp<Field<Type>> sample(const tmp<Field<Type>>&) const;
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#ifdef NoRepository
-    #include "cuttingPlaneTemplates.C"
-#endif
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //

src/sampling/sampledSurface/distanceSurface/distanceSurface.C

@@ -1,551 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2021 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "distanceSurface.H"
-#include "addToRunTimeSelectionTable.H"
-
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-namespace Foam
-{
-namespace sampledSurfaces
-{
-    defineTypeNameAndDebug(distanceSurface, 0);
-    addToRunTimeSelectionTable(sampledSurface, distanceSurface, word);
-}
-}
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-void Foam::sampledSurfaces::distanceSurface::createGeometry()
-{
-    if (debug)
-    {
-        Pout<< "distanceSurface::createGeometry :updating geometry." << endl;
-    }
-
-    // Clear any stored topologies
-    isoSurfPtr_.clear();
-
-    // Clear derived data
-    clearGeom();
-
-    // Get any subMesh
-    if (zoneKey_.size() && !subMeshPtr_.valid())
-    {
-        const polyBoundaryMesh& patches = mesh().boundaryMesh();
-
-        // Patch to put exposed internal faces into
-        const label exposedPatchi = patches.findPatchID(exposedPatchName_);
-
-        subMeshPtr_.reset
-        (
-            new fvMeshSubset(static_cast<const fvMesh&>(mesh()))
-        );
-        subMeshPtr_().setLargeCellSubset
-        (
-            labelHashSet(mesh().cellZones().findMatching(zoneKey_).used()),
-            exposedPatchi
-        );
-
-        DebugInfo
-            << "Allocating subset of size " << subMeshPtr_().subMesh().nCells()
-            << " with exposed faces into patch "
-            << patches[exposedPatchi].name() << endl;
-    }
-
-    // Select either the submesh or the underlying mesh
-    const fvMesh& mesh =
-    (
-        subMeshPtr_.valid()
-      ? subMeshPtr_().subMesh()
-      : static_cast<const fvMesh&>(this->mesh())
-    );
-
-
-    // Distance to cell centres
-    // ~~~~~~~~~~~~~~~~~~~~~~~~
-
-    cellDistancePtr_.reset
-    (
-        new volScalarField
-        (
-            IOobject
-            (
-                "cellDistance",
-                mesh.time().timeName(),
-                mesh.time(),
-                IOobject::NO_READ,
-                IOobject::NO_WRITE,
-                false
-            ),
-            mesh,
-            dimensionedScalar(dimLength, 0)
-        )
-    );
-    volScalarField& cellDistance = cellDistancePtr_();
-
-    // Internal field
-    {
-        const pointField& cc = mesh.C();
-        scalarField& fld = cellDistance.primitiveFieldRef();
-
-        List<pointIndexHit> nearest;
-        surfPtr_().findNearest
-        (
-            cc,
-            scalarField(cc.size(), great),
-            nearest
-        );
-
-        if (signed_)
-        {
-            List<volumeType> volType;
-
-            surfPtr_().getVolumeType(cc, volType);
-
-            forAll(volType, i)
-            {
-                volumeType vT = volType[i];
-
-                if (vT == volumeType::outside)
-                {
-                    fld[i] = Foam::mag(cc[i] - nearest[i].hitPoint());
-                }
-                else if (vT == volumeType::inside)
-                {
-                    fld[i] = -Foam::mag(cc[i] - nearest[i].hitPoint());
-                }
-                else
-                {
-                    FatalErrorInFunction
-                        << "getVolumeType failure, neither INSIDE or OUTSIDE"
-                        << exit(FatalError);
-                }
-            }
-        }
-        else
-        {
-            forAll(nearest, i)
-            {
-                fld[i] = Foam::mag(cc[i] - nearest[i].hitPoint());
-            }
-        }
-    }
-
-    volScalarField::Boundary& cellDistanceBf =
-        cellDistance.boundaryFieldRef();
-
-    // Patch fields
-    {
-        forAll(mesh.C().boundaryField(), patchi)
-        {
-            const pointField& cc = mesh.C().boundaryField()[patchi];
-            fvPatchScalarField& fld = cellDistanceBf[patchi];
-
-            List<pointIndexHit> nearest;
-            surfPtr_().findNearest
-            (
-                cc,
-                scalarField(cc.size(), great),
-                nearest
-            );
-
-            if (signed_)
-            {
-                List<volumeType> volType;
-
-                surfPtr_().getVolumeType(cc, volType);
-
-                forAll(volType, i)
-                {
-                    volumeType vT = volType[i];
-
-                    if (vT == volumeType::outside)
-                    {
-                        fld[i] = Foam::mag(cc[i] - nearest[i].hitPoint());
-                    }
-                    else if (vT == volumeType::inside)
-                    {
-                        fld[i] = -Foam::mag(cc[i] - nearest[i].hitPoint());
-                    }
-                    else
-                    {
-                        FatalErrorInFunction
-                            << "getVolumeType failure, "
-                            << "neither INSIDE or OUTSIDE"
-                            << exit(FatalError);
-                    }
-                }
-            }
-            else
-            {
-                forAll(nearest, i)
-                {
-                    fld[i] = Foam::mag(cc[i] - nearest[i].hitPoint());
-                }
-            }
-        }
-    }
-
-
-    // On processor patches the mesh.C() will already be the cell centre
-    // on the opposite side so no need to swap cellDistance.
-
-
-    // Distance to points
-    pointDistance_.setSize(mesh.nPoints());
-    {
-        const pointField& pts = mesh.points();
-
-        List<pointIndexHit> nearest;
-        surfPtr_().findNearest
-        (
-            pts,
-            scalarField(pts.size(), great),
-            nearest
-        );
-
-        if (signed_)
-        {
-            List<volumeType> volType;
-
-            surfPtr_().getVolumeType(pts, volType);
-
-            forAll(volType, i)
-            {
-                volumeType vT = volType[i];
-
-                if (vT == volumeType::outside)
-                {
-                    pointDistance_[i] =
-                        Foam::mag(pts[i] - nearest[i].hitPoint());
-                }
-                else if (vT == volumeType::inside)
-                {
-                    pointDistance_[i] =
-                        -Foam::mag(pts[i] - nearest[i].hitPoint());
-                }
-                else
-                {
-                    FatalErrorInFunction
-                        << "getVolumeType failure, neither INSIDE or OUTSIDE"
-                        << exit(FatalError);
-                }
-            }
-        }
-        else
-        {
-            forAll(nearest, i)
-            {
-                pointDistance_[i] = Foam::mag(pts[i]-nearest[i].hitPoint());
-            }
-        }
-    }
-
-
-    if (debug)
-    {
-        Pout<< "Writing cell distance:" << cellDistance.objectPath() << endl;
-        cellDistance.write();
-        pointScalarField pDist
-        (
-            IOobject
-            (
-                "pointDistance",
-                mesh.time().timeName(),
-                mesh.time(),
-                IOobject::NO_READ,
-                IOobject::NO_WRITE,
-                false
-            ),
-            pointMesh::New(mesh),
-            dimensionedScalar(dimLength, 0)
-        );
-        pDist.primitiveFieldRef() = pointDistance_;
-
-        Pout<< "Writing point distance:" << pDist.objectPath() << endl;
-        pDist.write();
-    }
-
-
-    //- Direct from cell field and point field.
-    isoSurfPtr_.reset
-    (
-        new isoSurface
-        (
-            mesh,
-            cellDistance,
-            pointDistance_,
-            distance_,
-            filter_
-        )
-    );
-
-    if (debug)
-    {
-        print(Pout);
-        Pout<< endl;
-    }
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::sampledSurfaces::distanceSurface::distanceSurface
-(
-    const word& name,
-    const polyMesh& mesh,
-    const dictionary& dict
-)
-:
-    sampledSurface(name, mesh, dict),
-    surfPtr_
-    (
-        searchableSurface::New
-        (
-            dict.lookup("surfaceType"),
-            IOobject
-            (
-                dict.lookupOrDefault("surfaceName", name),
-                mesh.time().constant(),
-                searchableSurface::geometryDir(mesh.time()),
-                mesh.time(),
-                IOobject::MUST_READ,
-                IOobject::NO_WRITE
-            ),
-            dict
-        )
-    ),
-    distance_(dict.lookup<scalar>("distance")),
-    signed_(readBool(dict.lookup("signed"))),
-    filter_
-    (
-        dict.found("filtering")
-      ? isoSurface::filterTypeNames_.read(dict.lookup("filtering"))
-      : isoSurface::filterType::full
-    ),
-    average_(dict.lookupOrDefault("average", false)),
-    zoneKey_(dict.lookupOrDefault("zone", wordRe::null)),
-    needsUpdate_(true),
-    subMeshPtr_(nullptr),
-    isoSurfPtr_(nullptr)
-{
-    if (zoneKey_.size())
-    {
-        dict.lookup("exposedPatchName") >> exposedPatchName_;
-
-        if (mesh.boundaryMesh().findPatchID(exposedPatchName_) == -1)
-        {
-            FatalErrorInFunction
-                << "Cannot find patch " << exposedPatchName_
-                << " in which to put exposed faces." << endl
-                << "Valid patches are " << mesh.boundaryMesh().names()
-                << exit(FatalError);
-        }
-
-        if (debug)
-        {
-            Info<< "Restricting to cellZone " << zoneKey_
-                << " with exposed internal faces into patch "
-                << exposedPatchName_ << endl;
-        }
-
-        if (mesh.cellZones().findIndex(zoneKey_) < 0)
-        {
-            WarningInFunction
-                << "cellZone " << zoneKey_
-                << " not found - using entire mesh" << endl;
-        }
-    }
-}
-
-
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-Foam::sampledSurfaces::distanceSurface::~distanceSurface()
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-bool Foam::sampledSurfaces::distanceSurface::needsUpdate() const
-{
-    return needsUpdate_;
-}
-
-
-bool Foam::sampledSurfaces::distanceSurface::expire()
-{
-    if (debug)
-    {
-        Pout<< "distanceSurface::expire :"
-            << " needsUpdate_:" << needsUpdate_ << endl;
-    }
-
-    // Clear derived data
-    clearGeom();
-
-    // already marked as expired
-    if (needsUpdate_)
-    {
-        return false;
-    }
-
-    needsUpdate_ = true;
-    return true;
-}
-
-
-bool Foam::sampledSurfaces::distanceSurface::update()
-{
-    if (debug)
-    {
-        Pout<< "distanceSurface::update :"
-            << " needsUpdate_:" << needsUpdate_ << endl;
-    }
-
-    if (!needsUpdate_)
-    {
-        return false;
-    }
-
-    createGeometry();
-
-    needsUpdate_ = false;
-    return true;
-}
-
-
-Foam::tmp<Foam::scalarField>
-Foam::sampledSurfaces::distanceSurface::sample
-(
-    const volScalarField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::vectorField>
-Foam::sampledSurfaces::distanceSurface::sample
-(
-    const volVectorField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::sphericalTensorField>
-Foam::sampledSurfaces::distanceSurface::sample
-(
-    const volSphericalTensorField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::symmTensorField>
-Foam::sampledSurfaces::distanceSurface::sample
-(
-    const volSymmTensorField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::tensorField>
-Foam::sampledSurfaces::distanceSurface::sample
-(
-    const volTensorField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::scalarField>
-Foam::sampledSurfaces::distanceSurface::interpolate
-(
-    const interpolation<scalar>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-
-Foam::tmp<Foam::vectorField>
-Foam::sampledSurfaces::distanceSurface::interpolate
-(
-    const interpolation<vector>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-Foam::tmp<Foam::sphericalTensorField>
-Foam::sampledSurfaces::distanceSurface::interpolate
-(
-    const interpolation<sphericalTensor>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-
-Foam::tmp<Foam::symmTensorField>
-Foam::sampledSurfaces::distanceSurface::interpolate
-(
-    const interpolation<symmTensor>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-
-Foam::tmp<Foam::tensorField>
-Foam::sampledSurfaces::distanceSurface::interpolate
-(
-    const interpolation<tensor>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-
-void Foam::sampledSurfaces::distanceSurface::print(Ostream& os) const
-{
-    os  << "distanceSurface: " << name() << " :"
-        << "  surface:" << surfPtr_().name()
-        << "  distance:" << distance_
-        << "  faces:" << faces().size()
-        << "  points:" << points().size();
-}
-
-
-// ************************************************************************* //

src/sampling/sampledSurface/distanceSurface/distanceSurface.H

@@ -1,277 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2021 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-Class
-     Foam::sampledSurfaces::distanceSurface
-
-Description
-    A sampledSurface defined by a distance to a surface.
-
-    Example:
-    \verbatim
-    {
-        type        distanceSurface;
-        surfaceType searchableBox;
-        min         (-1 -1 -1);
-        max         (1 1 1);
-        distance    0.1;
-        signed      yes;
-        filtering   full;
-        interpolate yes;
-    }
-    \endverbatim
-
-Usage
-    \table
-        Property     | Description             | Required    | Default value
-        surfaceType  | the type of surface to sample from | yes |
-        distance     | the distance from which to sample the surface | yes |
-        signed       | sample only on one side of the surface, specified by \\
-                       the sign of the distance | yes |
-        filtering    | the level of filtering to perform on the iso-surface \\
-                                               | no          | full
-        interpolate  | interpolate values to the surface points | no | no
-    \endtable
-
-See also
-    Foam::isoSurface
-
-SourceFiles
-    distanceSurface.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef sampledSurfaces_distanceSurface_H
-#define sampledSurfaces_distanceSurface_H
-
-#include "sampledSurface.H"
-#include "searchableSurface.H"
-#include "isoSurface.H"
-#include "fvMeshSubset.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-namespace sampledSurfaces
-{
-
-/*---------------------------------------------------------------------------*\
-                       Class distanceSurface Declaration
-\*---------------------------------------------------------------------------*/
-
-class distanceSurface
-:
-    public sampledSurface
-{
-    // Private Data
-
-        //- Surface
-        const autoPtr<searchableSurface> surfPtr_;
-
-        //- Distance value
-        const scalar distance_;
-
-        //- Signed distance
-        const bool signed_;
-
-        //- Whether to coarsen
-        const isoSurface::filterType filter_;
-
-        //- Whether to recalculate cell values as average of point values
-        const Switch average_;
-
-        //- If restricted to zones, name of this zone or a regular expression
-        wordRe zoneKey_;
-
-        //- For zones: patch to put exposed faces into
-        mutable word exposedPatchName_;
-
-        //- Track if the surface needs an update
-        mutable bool needsUpdate_;
-
-        //- Optional subsetted mesh
-        autoPtr<fvMeshSubset> subMeshPtr_;
-
-        //- Distance to cell centres
-        autoPtr<volScalarField> cellDistancePtr_;
-
-        //- Distance to points
-        scalarField pointDistance_;
-
-        //- Constructed iso surface
-        autoPtr<isoSurface> isoSurfPtr_;
-
-
-    // Private Member Functions
-
-        //- Create iso surface
-        void createGeometry();
-
-        //- Sample field on faces
-        template<class Type>
-        tmp<Field<Type>> sampleField
-        (
-            const GeometricField<Type, fvPatchField, volMesh>& vField
-        ) const;
-
-
-        template<class Type>
-        tmp<Field<Type>>
-        interpolateField(const interpolation<Type>&) const;
-
-
-public:
-
-    //- Runtime type information
-    TypeName("distanceSurface");
-
-
-    // Constructors
-
-        //- Construct from dictionary
-        distanceSurface
-        (
-            const word& name,
-            const polyMesh& mesh,
-            const dictionary& dict
-        );
-
-
-    //- Destructor
-    virtual ~distanceSurface();
-
-
-    // Member Functions
-
-        //- Does the surface need an update?
-        virtual bool needsUpdate() const;
-
-        //- Mark the surface as needing an update.
-        //  May also free up unneeded data.
-        //  Return false if surface was already marked as expired.
-        virtual bool expire();
-
-        //- Update the surface as required.
-        //  Do nothing (and return false) if no update was needed
-        virtual bool update();
-
-        //- Points of surface
-        virtual const pointField& points() const
-        {
-            return surface().points();
-        }
-
-        //- Faces of surface
-        virtual const faceList& faces() const
-        {
-            return surface().faces();
-        }
-
-        const isoSurface& surface() const
-        {
-            return isoSurfPtr_();
-        }
-
-        //- Sample field on surface
-        virtual tmp<scalarField> sample
-        (
-            const volScalarField&
-        ) const;
-
-        //- Sample field on surface
-        virtual tmp<vectorField> sample
-        (
-            const volVectorField&
-        ) const;
-
-        //- Sample field on surface
-        virtual tmp<sphericalTensorField> sample
-        (
-            const volSphericalTensorField&
-        ) const;
-
-        //- Sample field on surface
-        virtual tmp<symmTensorField> sample
-        (
-            const volSymmTensorField&
-        ) const;
-
-        //- Sample field on surface
-        virtual tmp<tensorField> sample
-        (
-            const volTensorField&
-        ) const;
-
-
-        //- Interpolate field on surface
-        virtual tmp<scalarField> interpolate
-        (
-            const interpolation<scalar>&
-        ) const;
-
-        //- Interpolate field on surface
-        virtual tmp<vectorField> interpolate
-        (
-            const interpolation<vector>&
-        ) const;
-
-        //- Interpolate field on surface
-        virtual tmp<sphericalTensorField> interpolate
-        (
-            const interpolation<sphericalTensor>&
-        ) const;
-
-        //- Interpolate field on surface
-        virtual tmp<symmTensorField> interpolate
-        (
-            const interpolation<symmTensor>&
-        ) const;
-
-        //- Interpolate field on surface
-        virtual tmp<tensorField> interpolate
-        (
-            const interpolation<tensor>&
-        ) const;
-
-        //- Write
-        virtual void print(Ostream&) const;
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace sampledSurfaces
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#ifdef NoRepository
-    #include "distanceSurfaceTemplates.C"
-#endif
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //

src/sampling/sampledSurface/distanceSurface/distanceSurfaceTemplates.C

@@ -1,97 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2020 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "distanceSurface.H"
-#include "volPointInterpolation.H"
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-template<class Type>
-Foam::tmp<Foam::Field<Type>>
-Foam::sampledSurfaces::distanceSurface::sampleField
-(
-    const GeometricField<Type, fvPatchField, volMesh>& vField
-) const
-{
-    return tmp<Field<Type>>
-    (
-        new Field<Type>(vField, isoSurfPtr_().meshCells())
-    );
-}
-
-
-template<class Type>
-Foam::tmp<Foam::Field<Type>>
-Foam::sampledSurfaces::distanceSurface::interpolateField
-(
-    const interpolation<Type>& interpolator
-) const
-{
-    // Get fields to sample. Assume volPointInterpolation!
-    const GeometricField<Type, fvPatchField, volMesh>& volFld =
-        interpolator.psi();
-
-    if (subMeshPtr_.valid())
-    {
-        tmp<GeometricField<Type, fvPatchField, volMesh>> tvolSubFld =
-            subMeshPtr_().interpolate(volFld);
-
-        const GeometricField<Type, fvPatchField, volMesh>& volSubFld =
-            tvolSubFld();
-
-        tmp<GeometricField<Type, pointPatchField, pointMesh>> tpointSubFld =
-            volPointInterpolation::New(volSubFld.mesh()).interpolate(volSubFld);
-
-        return isoSurfPtr_().interpolate
-        (
-            (
-                average_
-              ? pointAverage(tpointSubFld())()
-              : volSubFld
-            ),
-            tpointSubFld()
-        );
-    }
-    else
-    {
-        tmp<GeometricField<Type, pointPatchField, pointMesh>> pointFld
-        (
-            volPointInterpolation::New(volFld.mesh()).interpolate(volFld)
-        );
-
-        return isoSurfPtr_().interpolate
-        (
-            (
-                average_
-              ? pointAverage(pointFld())()
-              : volFld
-            ),
-            pointFld()
-        );
-    }
-}
-
-
-// ************************************************************************* //

src/sampling/sampledSurface/isoSurface/isoSurface.H

@@ -1,297 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2021 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-Class
-    Foam::isoSurface
-
-Description
-    Marching tet iso surface algorithm with filtering to remove unnecessary
-    topology.
-
-    Three levels of filtering are possible:
-
-    - none:
-        Don't filter; the full intersection between the iso surface and the
-        tetrahedral decomposition of the mesh is generated
-
-    - partial:
-        Remove points from vertex to cell-centre edges and merge triangles that
-        form a contiguous cut through a single cell
-
-    - full:
-        As partial, but also remove points from face-diagonals and merge edges
-        that originate from the same face
-
-    - clean:
-        As full, but also remove faces that were made non-manifold by the point
-        filtering process. Note that this process is not entirely robust. The
-        removal process can result in entire contiguous sections of the surface
-        being removed. If small imperfections in the surface can be tolerated,
-        then it is advised to use "full" filtering instead.
-
-SourceFiles
-    isoSurface.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef isoSurface_H
-#define isoSurface_H
-
-#include "labelPair.H"
-#include "pointIndexHit.H"
-#include "PackedBoolList.H"
-#include "MeshedSurface.H"
-#include "edgeList.H"
-#include "NamedEnum.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-class polyMesh;
-class tetMatcher;
-
-/*---------------------------------------------------------------------------*\
-                         Class isoSurface Declaration
-\*---------------------------------------------------------------------------*/
-
-class isoSurface
-:
-    public MeshedSurface<face>
-{
-public:
-
-        enum class filterType
-        {
-            none,
-            partial,
-            full,
-            clean
-        };
-
-        static const NamedEnum<filterType, 4> filterTypeNames_;
-
-
-private:
-
-    // Private Data
-
-        enum class cellCutType
-        {
-            notCut,     // Not cut
-            sphere,     // All edges to cell centre cut
-            cut         // Normal cut
-        };
-
-
-        //- Reference to mesh
-        const polyMesh& mesh_;
-
-        const scalarField& cVals_;
-
-        const scalarField& pVals_;
-
-        //- Iso value
-        const scalar iso_;
-
-        //- Corrected version of tetBasePtIs
-        labelList tetBasePtIs_;
-
-        //- Per point: originating mesh vertex/cc. See encoding above
-        edgeList pointToVerts_;
-
-        //- For every face the original cell in mesh
-        labelList meshCells_;
-
-        //- For every point the originating face in mesh
-        labelList pointToFace_;
-
-
-    // Private Member Functions
-
-        scalar minTetQ
-        (
-            const label facei,
-            const label faceBasePtI
-        ) const;
-
-        void fixTetBasePtIs();
-
-        //- Does any edge of triangle cross iso value?
-        bool isTriCut
-        (
-            const triFace& tri,
-            const scalarField& pointValues
-        ) const;
-
-        //- Determine whether cell is cut
-        cellCutType calcCutType
-        (
-            const bool isTet,
-            const label
-        ) const;
-
-        //- Determine for all mesh whether cell is cut
-        label calcCutTypes
-        (
-            tetMatcher& tet,
-            List<cellCutType>& cellCutTypes
-        );
-
-        //- Generate single point on edge
-        label generatePoint
-        (
-            const label facei,
-            const bool edgeIsDiag,
-            const edge& vertices,
-
-            DynamicList<edge>& pointToVerts,
-            DynamicList<label>& pointToFace,
-            DynamicList<bool>& pointFromDiag,
-            EdgeMap<label>& vertsToPoint
-        ) const;
-
-        //- Generate triangles from tet
-        void generateTriPoints
-        (
-            const label facei,
-            const FixedList<scalar, 4>& s,
-            const FixedList<point, 4>& p,
-            const FixedList<label, 4>& pIndex,
-            const FixedList<bool, 6>& edgeIsDiag,
-
-            DynamicList<edge>& pointToVerts,
-            DynamicList<label>& pointToFace,
-            DynamicList<bool>& pointFromDiag,
-
-            EdgeMap<label>& vertsToPoint,
-            DynamicList<label>& verts
-        ) const;
-
-        //- Generate triangles from cell
-        void generateTriPoints
-        (
-            const label celli,
-            const bool isTet,
-
-            DynamicList<edge>& pointToVerts,
-            DynamicList<label>& pointToFace,
-            DynamicList<bool>& pointFromDiag,
-
-            EdgeMap<label>& vertsToPoint,
-            DynamicList<label>& verts,
-            DynamicList<label>& faceLabels
-        ) const;
-
-
-        // Simplification
-
-            void triangulateOutside
-            (
-                const bool filterDiag,
-                const PrimitivePatch<SubList<face>, const pointField&>& pp,
-                const boolList& pointFromDiag,
-                const labelList& pointToFace,
-                const label cellID,
-
-                DynamicList<face>& compactFaces,
-                DynamicList<label>& compactCellIDs
-            ) const;
-
-            MeshedSurface<face> removeInsidePoints
-            (
-                const bool filterDiag,
-                const MeshedSurface<face>& s,
-                const boolList& pointFromDiag,
-                const labelList& pointToFace,
-                const labelList& start,              // Per cell:starting tri
-                DynamicList<label>& pointCompactMap, // Per point the original
-                DynamicList<label>& compactCellIDs   // Per face the cellID
-            ) const;
-
-
-public:
-
-    //- Runtime type information
-    TypeName("isoSurface");
-
-
-    // Constructors
-
-        //- Construct from dictionary
-        isoSurface
-        (
-            const polyMesh& mesh,
-            const scalarField& cellValues,
-            const scalarField& pointValues,
-            const scalar iso,
-            const filterType filter
-        );
-
-
-    // Member Functions
-
-        //- For every face original cell in mesh
-        const labelList& meshCells() const
-        {
-            return meshCells_;
-        }
-
-        //- For every point originating face (pyramid) in mesh
-        const labelList& pointToFace() const
-        {
-            return pointToFace_;
-        }
-
-        //- Per point: originating mesh vertex/cc. See encoding above
-        const edgeList& pointToVerts() const
-        {
-            return pointToVerts_;
-        }
-
-        //- Interpolates cCoords,pCoords.
-        template<class Type>
-        tmp<Field<Type>> interpolate
-        (
-            const Field<Type>& cCoords,
-            const Field<Type>& pCoords
-        ) const;
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#ifdef NoRepository
-    #include "isoSurfaceTemplates.C"
-#endif
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //

src/sampling/sampledSurface/isoSurface/isoSurfaceTemplates.C

@@ -1,91 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-template<class Type>
-Foam::tmp<Foam::Field<Type>>
-Foam::isoSurface::interpolate
-(
-    const Field<Type>& cellCoords,
-    const Field<Type>& pointCoords
-) const
-{
-    tmp<Field<Type>> tfld(new Field<Type>(pointToVerts_.size()));
-    Field<Type>& fld = tfld.ref();
-
-    forAll(pointToVerts_, i)
-    {
-        scalar s0;
-        Type p0;
-        {
-            label v0 = pointToVerts_[i][0];
-            if (v0 < mesh_.nPoints())
-            {
-                s0 = pVals_[v0];
-                p0 = pointCoords[v0];
-            }
-            else
-            {
-                label celli = v0-mesh_.nPoints();
-                s0 = cVals_[celli];
-                p0 = cellCoords[celli];
-            }
-        }
-
-        scalar s1;
-        Type p1;
-        {
-            label v1 = pointToVerts_[i][1];
-            if (v1 < mesh_.nPoints())
-            {
-                s1 = pVals_[v1];
-                p1 = pointCoords[v1];
-            }
-            else
-            {
-                label celli = v1-mesh_.nPoints();
-                s1 = cVals_[celli];
-                p1 = cellCoords[celli];
-            }
-        }
-
-        scalar d = s1-s0;
-        if (mag(d) > VSMALL)
-        {
-            scalar s = (iso_-s0)/d;
-            fld[i] = s*p1+(1.0-s)*p0;
-        }
-        else
-        {
-            fld[i] = 0.5*(p0+p1);
-        }
-    }
-
-    return tfld;
-}
-
-
-// ************************************************************************* //

src/sampling/sampledSurface/isoSurface/sampledIsoSurface.C

@@ -1,425 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2021 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "sampledIsoSurface.H"
-#include "addToRunTimeSelectionTable.H"
-
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-namespace Foam
-{
-namespace sampledSurfaces
-{
-    defineTypeNameAndDebug(isoSurface, 0);
-    addToRunTimeSelectionTable(sampledSurface, isoSurface, word);
-}
-}
-
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-bool Foam::sampledSurfaces::isoSurface::updateGeometry() const
-{
-    const fvMesh& fvm = static_cast<const fvMesh&>(mesh());
-
-    // no update needed
-    if (fvm.time().timeIndex() == prevTimeIndex_)
-    {
-        return false;
-    }
-
-    prevTimeIndex_ = fvm.time().timeIndex();
-
-    // Clear derived data
-    sampledSurface::clearGeom();
-
-    // Optionally read volScalarField
-    autoPtr<volScalarField> readFieldPtr_;
-
-    // 1. see if field in database
-    // 2. see if field can be read
-    const volScalarField* cellFldPtr = nullptr;
-    if (fvm.foundObject<volScalarField>(isoField_))
-    {
-        if (debug)
-        {
-            InfoInFunction << "Lookup " << isoField_ << endl;
-        }
-
-        cellFldPtr = &fvm.lookupObject<volScalarField>(isoField_);
-    }
-    else
-    {
-        // Bit of a hack. Read field and store.
-
-        if (debug)
-        {
-            InfoInFunction
-                << "Reading " << isoField_
-                << " from time " <<fvm.time().timeName()
-                << endl;
-        }
-
-        readFieldPtr_.reset
-        (
-            new volScalarField
-            (
-                IOobject
-                (
-                    isoField_,
-                    fvm.time().timeName(),
-                    fvm,
-                    IOobject::MUST_READ,
-                    IOobject::NO_WRITE,
-                    false
-                ),
-                fvm
-            )
-        );
-
-        cellFldPtr = readFieldPtr_.operator->();
-    }
-    const volScalarField& cellFld = *cellFldPtr;
-
-    tmp<pointScalarField> pointFld
-    (
-        volPointInterpolation::New(fvm).interpolate(cellFld)
-    );
-
-    PtrList<Foam::isoSurface> isos(isoVals_.size());
-    forAll(isos, isoi)
-    {
-        isos.set
-        (
-            isoi,
-            new Foam::isoSurface
-            (
-                fvm,
-                cellFld.primitiveField(),
-                pointFld().primitiveField(),
-                isoVals_[isoi],
-                filter_
-            )
-        );
-    }
-
-    if (isos.size() == 1)
-    {
-        // Straight transfer
-        const_cast<isoSurface&>
-        (
-            *this
-        ).MeshedSurface<face>::transfer(isos[0]);
-        meshCells_ = isos[0].meshCells();
-    }
-    else
-    {
-        label nFaces = 0;
-        label nPoints = 0;
-        forAll(isos, isoi)
-        {
-            nFaces += isos[isoi].size();
-            nPoints += isos[isoi].points().size();
-        }
-
-        faceList allFaces(nFaces);
-        labelList allMeshCells(nFaces);
-        pointField allPoints(nPoints);
-
-        nFaces = 0;
-        nPoints = 0;
-        forAll(isos, isoi)
-        {
-            Foam::isoSurface& iso = isos[isoi];
-
-            SubList<face> subAll(allFaces, iso.size(), nFaces);
-            subAll = iso;
-            // Offset point indices
-            if (nPoints > 0)
-            {
-                forAll(subAll, i)
-                {
-                    face& f = subAll[i];
-                    forAll(f, fp)
-                    {
-                        f[fp] += nPoints;
-                    }
-                }
-            }
-            SubList<label>(allMeshCells, iso.size(), nFaces) = iso.meshCells();
-            nFaces += iso.size();
-
-            const pointField& pts = iso.points();
-            SubList<point>(allPoints, pts.size(), nPoints) = pts;
-            nPoints += pts.size();
-
-            // Clear out asap
-            isos.set(isoi, nullptr);
-        }
-
-        if (nFaces != allFaces.size() || nPoints != allPoints.size())
-        {
-            FatalErrorInFunction << "nFaces:" << nFaces
-                << " nPoints:" << nPoints << exit(FatalError);
-        }
-
-
-        surfZoneList allZones(1);
-        allZones[0] = surfZone
-        (
-            "allFaces",
-            allFaces.size(),    // size
-            0,                  // start
-            0                   // index
-        );
-
-        // Transfer
-        const_cast<isoSurface&>
-        (
-            *this
-        ).MeshedSurface<face>::reset
-        (
-            move(allPoints),
-            move(allFaces),
-            move(allZones)
-        );
-        meshCells_.transfer(allMeshCells);
-    }
-    if (debug)
-    {
-        Pout<< "sampledSurfaces::isoSurface::updateGeometry() : "
-               "constructed iso:"
-            << nl
-            << "    filtering      : "
-            << Foam::isoSurface::filterTypeNames_[filter_] << nl
-            << "    isoField       : " << isoField_ << nl;
-        if (isoVals_.size() == 1)
-        {
-            Pout<< "    isoValue       : " << isoVals_[0] << nl;
-        }
-        else
-        {
-            Pout<< "    isoValues      : " << isoVals_ << nl;
-        }
-        Pout<< "    points         : " << points().size() << nl
-            << "    faces          : " << faces().size() << nl
-            << "    cut cells      : " << meshCells_.size() << endl;
-    }
-
-    return true;
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::sampledSurfaces::isoSurface::isoSurface
-(
-    const word& name,
-    const polyMesh& mesh,
-    const dictionary& dict
-)
-:
-    sampledSurface(name, mesh, dict),
-    isoField_(dict.lookup("isoField")),
-    isoVals_
-    (
-        dict.found("isoValues")
-      ? scalarField(dict.lookup("isoValues"))
-      : scalarField(1, dict.lookup<scalar>("isoValue"))
-    ),
-    filter_
-    (
-        dict.found("filtering")
-      ? Foam::isoSurface::filterTypeNames_.read(dict.lookup("filtering"))
-      : Foam::isoSurface::filterType::full
-    ),
-    prevTimeIndex_(-1),
-    meshCells_(0)
-{}
-
-
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-Foam::sampledSurfaces::isoSurface::~isoSurface()
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-Foam::wordList Foam::sampledSurfaces::isoSurface::fields() const
-{
-    return wordList(isoField_);
-}
-
-
-bool Foam::sampledSurfaces::isoSurface::needsUpdate() const
-{
-    const fvMesh& fvm = static_cast<const fvMesh&>(mesh());
-
-    return fvm.time().timeIndex() != prevTimeIndex_;
-}
-
-
-bool Foam::sampledSurfaces::isoSurface::expire()
-{
-    // Clear derived data
-    sampledSurface::clearGeom();
-    MeshedSurface<face>::clearGeom();
-
-    // already marked as expired
-    if (prevTimeIndex_ == -1)
-    {
-        return false;
-    }
-
-    // force update
-    prevTimeIndex_ = -1;
-    return true;
-}
-
-
-bool Foam::sampledSurfaces::isoSurface::update()
-{
-    return updateGeometry();
-}
-
-
-Foam::tmp<Foam::scalarField>
-Foam::sampledSurfaces::isoSurface::sample
-(
-    const volScalarField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::vectorField>
-Foam::sampledSurfaces::isoSurface::sample
-(
-    const volVectorField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::sphericalTensorField>
-Foam::sampledSurfaces::isoSurface::sample
-(
-    const volSphericalTensorField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::symmTensorField>
-Foam::sampledSurfaces::isoSurface::sample
-(
-    const volSymmTensorField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::tensorField>
-Foam::sampledSurfaces::isoSurface::sample
-(
-    const volTensorField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::scalarField>
-Foam::sampledSurfaces::isoSurface::interpolate
-(
-    const interpolation<scalar>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-
-Foam::tmp<Foam::vectorField>
-Foam::sampledSurfaces::isoSurface::interpolate
-(
-    const interpolation<vector>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-Foam::tmp<Foam::sphericalTensorField>
-Foam::sampledSurfaces::isoSurface::interpolate
-(
-    const interpolation<sphericalTensor>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-
-Foam::tmp<Foam::symmTensorField>
-Foam::sampledSurfaces::isoSurface::interpolate
-(
-    const interpolation<symmTensor>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-
-Foam::tmp<Foam::tensorField>
-Foam::sampledSurfaces::isoSurface::interpolate
-(
-    const interpolation<tensor>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-
-void Foam::sampledSurfaces::isoSurface::print(Ostream& os) const
-{
-    os  << "isoSurface: " << name() << " :"
-        << "  field:" << isoField_;
-    if (isoVals_.size() == 1)
-    {
-        os << "  value:" << isoVals_[0];
-    }
-    else
-    {
-        os << "  values:" << isoVals_;
-    }
-    os  << "  faces:" << faces().size()
-        << "  points:" << points().size();
-}
-
-
-// ************************************************************************* //

src/sampling/sampledSurface/isoSurface/sampledIsoSurface.H

@@ -1,252 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2021 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-Class
-    Foam::sampledSurfaces::isoSurface
-
-Description
-    A sampledSurface defined by a surface of iso value.
-
-    Example:
-    \verbatim
-    {
-        type        isoSurface;
-        isoField    p;
-        isoValue    1e5;
-        filtering   full;
-        interpolate yes;
-    }
-    \endverbatim
-
-Usage
-    \table
-        Property     | Description             | Required    | Default value
-        isoField     | the field to get an iso-surface of | yes |
-        isoValue     | the iso-surface value   | yes         |
-        filtering    | the level of filtering to perform on the iso-surface \\
-                                               | no          | full
-        interpolate  | interpolate values to the surface points | no | no
-    \endtable
-
-See also
-    Foam::isoSurface
-
-SourceFiles
-    sampledIsoSurface.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef sampledIsoSurface_H
-#define sampledIsoSurface_H
-
-#include "sampledSurface.H"
-#include "isoSurface.H"
-#include "MeshedSurface.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-namespace sampledSurfaces
-{
-
-/*---------------------------------------------------------------------------*\
-                         Class isoSurface Declaration
-\*---------------------------------------------------------------------------*/
-
-class isoSurface
-:
-    public sampledSurface,
-    public MeshedSurface<face>
-{
-    // Private Data
-
-        //- Field to get isoSurface of
-        const word isoField_;
-
-        //- Iso value
-        const scalarField isoVals_;
-
-        //- Whether to coarsen
-        const Foam::isoSurface::filterType filter_;
-
-
-        // Recreated for every isoSurface
-
-            //- Time at last call, also track it surface needs an update
-            mutable label prevTimeIndex_;
-
-            //- For every triangle/face the original cell in mesh
-            mutable labelList meshCells_;
-
-
-    // Private Member Functions
-
-        //- Create iso surface (if time has changed)
-        //  Do nothing (and return false) if no update was needed
-        bool updateGeometry() const;
-
-        //- Sample field on faces
-        template<class Type>
-        tmp<Field<Type>> sampleField
-        (
-            const GeometricField<Type, fvPatchField, volMesh>& vField
-        ) const;
-
-        //- Interpolate field to vertices
-        template<class Type>
-        tmp<Field<Type>>
-        interpolateField(const interpolation<Type>&) const;
-
-
-public:
-
-    //- Runtime type information
-    TypeName("isoSurface");
-
-
-    // Constructors
-
-        //- Construct from dictionary
-        isoSurface
-        (
-            const word& name,
-            const polyMesh& mesh,
-            const dictionary& dict
-        );
-
-
-    //- Destructor
-    virtual ~isoSurface();
-
-
-    // Member Functions
-
-        //- Return the list of fields required
-        virtual wordList fields() const;
-
-        //- Does the surface need an update?
-        virtual bool needsUpdate() const;
-
-        //- Mark the surface as needing an update.
-        //  May also free up unneeded data.
-        //  Return false if surface was already marked as expired.
-        virtual bool expire();
-
-        //- Update the surface as required.
-        //  Do nothing (and return false) if no update was needed
-        virtual bool update();
-
-        //- Points of surface
-        virtual const pointField& points() const
-        {
-            return MeshedSurface<face>::points();
-        }
-
-        //- Faces of surface
-        virtual const faceList& faces() const
-        {
-            return *this;
-        }
-
-        //- Sample field on surface
-        virtual tmp<scalarField> sample
-        (
-            const volScalarField&
-        ) const;
-
-        //- Sample field on surface
-        virtual tmp<vectorField> sample
-        (
-            const volVectorField&
-        ) const;
-
-        //- Sample field on surface
-        virtual tmp<sphericalTensorField> sample
-        (
-            const volSphericalTensorField&
-        ) const;
-
-        //- Sample field on surface
-        virtual tmp<symmTensorField> sample
-        (
-            const volSymmTensorField&
-        ) const;
-
-        //- Sample field on surface
-        virtual tmp<tensorField> sample
-        (
-            const volTensorField&
-        ) const;
-
-        //- Interpolate field on surface
-        virtual tmp<scalarField> interpolate
-        (
-            const interpolation<scalar>&
-        ) const;
-
-        //- Interpolate field on surface
-        virtual tmp<vectorField> interpolate
-        (
-            const interpolation<vector>&
-        ) const;
-
-        //- Interpolate field on surface
-        virtual tmp<sphericalTensorField> interpolate
-        (
-            const interpolation<sphericalTensor>&
-        ) const;
-
-        //- Interpolate field on surface
-        virtual tmp<symmTensorField> interpolate
-        (
-            const interpolation<symmTensor>&
-        ) const;
-
-        //- Interpolate field on surface
-        virtual tmp<tensorField> interpolate
-        (
-            const interpolation<tensor>&
-        ) const;
-
-        //- Write
-        virtual void print(Ostream&) const;
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace sampledSurfaces
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#ifdef NoRepository
-    #include "sampledIsoSurfaceTemplates.C"
-#endif
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //

src/sampling/sampledSurface/isoSurface/sampledIsoSurfaceTemplates.C

@@ -1,85 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2019 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "sampledIsoSurface.H"
-#include "volPointInterpolation.H"
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-template<class Type>
-Foam::tmp<Foam::Field<Type>>
-Foam::sampledSurfaces::isoSurface::sampleField
-(
-    const GeometricField<Type, fvPatchField, volMesh>& vField
-) const
-{
-    // Recreate geometry if time has changed
-    updateGeometry();
-
-    return tmp<Field<Type>>(new Field<Type>(vField, meshCells_));
-}
-
-
-template<class Type>
-Foam::tmp<Foam::Field<Type>>
-Foam::sampledSurfaces::isoSurface::interpolateField
-(
-    const interpolation<Type>& interpolator
-) const
-{
-    // Recreate geometry if time has changed
-    updateGeometry();
-
-    // One value per point
-    tmp<Field<Type>> tvalues(new Field<Type>(points().size()));
-    Field<Type>& values = tvalues.ref();
-
-    boolList pointDone(points().size(), false);
-
-    forAll(faces(), cutFacei)
-    {
-        const face& f = faces()[cutFacei];
-
-        forAll(f, faceVertI)
-        {
-            label pointi = f[faceVertI];
-
-            if (!pointDone[pointi])
-            {
-                values[pointi] = interpolator.interpolate
-                (
-                    points()[pointi],
-                    meshCells_[cutFacei]
-                );
-                pointDone[pointi] = true;
-            }
-        }
-    }
-
-    return tvalues;
-}
-
-
-// ************************************************************************* //

src/sampling/sampledSurface/sampledCutPlane/sampledCutPlane.C

@@ -0,0 +1,111 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "sampledCutPlane.H"
+#include "emptyFvPatchFields.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace sampledSurfaces
+{
+    defineTypeNameAndDebug(cutPlane, 0);
+    addToRunTimeSelectionTable(sampledSurface, cutPlane, word);
+
+    // Backwards compatible lookup as "cuttingPlane"
+    addNamedToRunTimeSelectionTable
+    (
+        sampledSurface,
+        cutPlane,
+        word,
+        cuttingPlane
+    );
+}
+}
+
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+Foam::autoPtr<Foam::cutPolyIsoSurface>
+Foam::sampledSurfaces::cutPlane::calcIsoSurf() const
+{
+    // Compute the distance from the mesh points to the plane
+    scalarField pointDistance(mesh().nPoints());
+    {
+        forAll(mesh().points(), pointi)
+        {
+            pointDistance[pointi] =
+                plane_.signedDistance(mesh().points()[pointi]);
+        }
+    }
+
+    // Construct an iso-surface at the given distance
+    return autoPtr<cutPolyIsoSurface>
+    (
+        new cutPolyIsoSurface(mesh(), pointDistance, 0, zoneIDs())
+    );
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::sampledSurfaces::cutPlane::cutPlane
+(
+    const word& name,
+    const polyMesh& mesh,
+    const dictionary& dict
+)
+:
+    sampledIsoSurfaceSurface(name, mesh, dict),
+    plane_(dict)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::sampledSurfaces::cutPlane::~cutPlane()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+bool Foam::sampledSurfaces::cutPlane::needsUpdate() const
+{
+    return timeIndex() == -1;
+}
+
+
+void Foam::sampledSurfaces::cutPlane::print(Ostream& os) const
+{
+    os  << "cutPlane: " << name() << " :"
+        << "  plane:" << plane_
+        << "  faces:" << faces().size()
+        << "  points:" << points().size();
+}
+
+
+// ************************************************************************* //

src/sampling/sampledSurface/sampledCutPlane/sampledCutPlane.H

@@ -0,0 +1,128 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::sampledSurfaces::cutPlane
+
+Description
+    A sampledSurface defined by a plane.
+
+    Example:
+    \verbatim
+    {
+        type        cutPlane;
+        planeType   pointAndNormal;
+        point       (0 0 0);
+        normal      (0 0 1);
+        interpolate yes;
+    }
+    \endverbatim
+
+Usage
+    \table
+        Property     | Description             | Required    | Default value
+        planeType    | the method of specification of the plane | yes |
+        interpolate  | interpolate values to the surface points | no | no
+    \endtable
+
+See also
+    Foam::plane
+
+SourceFiles
+    sampledCutPlane.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef sampledCutPlane_H
+#define sampledCutPlane_H
+
+#include "sampledIsoSurfaceSurface.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace sampledSurfaces
+{
+
+/*---------------------------------------------------------------------------*\
+                        Class cutPlane Declaration
+\*---------------------------------------------------------------------------*/
+
+class cutPlane
+:
+    public sampledIsoSurfaceSurface
+{
+    // Private Data
+
+        //- Plane
+        const plane plane_;
+
+
+    // Private Member Functions
+
+        //- Generate the iso surface
+        virtual autoPtr<cutPolyIsoSurface> calcIsoSurf() const;
+
+
+public:
+
+    //- Runtime type information
+    TypeName("cutPlane");
+
+
+    // Constructors
+
+        //- Construct from dictionary
+        cutPlane
+        (
+            const word& name,
+            const polyMesh& mesh,
+            const dictionary& dict
+        );
+
+
+    //- Destructor
+    virtual ~cutPlane();
+
+
+    // Member Functions
+
+        //- Does the surface need an update?
+        virtual bool needsUpdate() const;
+
+        //- Write
+        virtual void print(Ostream&) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace sampledSurfaces
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/sampling/sampledSurface/sampledCuttingPlane/sampledCuttingPlane.C

@@ -1,458 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2021 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "sampledCuttingPlane.H"
-#include "emptyFvPatchFields.H"
-#include "addToRunTimeSelectionTable.H"
-
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-namespace Foam
-{
-namespace sampledSurfaces
-{
-    defineTypeNameAndDebug(cuttingPlane, 0);
-    addToRunTimeSelectionTable(sampledSurface, cuttingPlane, word);
-}
-}
-
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-void Foam::sampledSurfaces::cuttingPlane::createGeometry()
-{
-    if (debug)
-    {
-        Pout<< "cuttingPlane::createGeometry :updating geometry."
-            << endl;
-    }
-
-    // Clear any stored topologies
-    facesPtr_.clear();
-    isoSurfPtr_.clear();
-    pointDistance_.clear();
-    cellDistancePtr_.clear();
-
-    // Clear derived data
-    clearGeom();
-
-    // Get any subMesh
-    if (zoneKey_.size() && !subMeshPtr_.valid())
-    {
-        const polyBoundaryMesh& patches = mesh().boundaryMesh();
-
-        // Patch to put exposed internal faces into
-        const label exposedPatchi = patches.findPatchID(exposedPatchName_);
-
-        subMeshPtr_.reset
-        (
-            new fvMeshSubset(static_cast<const fvMesh&>(mesh()))
-        );
-        subMeshPtr_().setLargeCellSubset
-        (
-            labelHashSet(mesh().cellZones().findMatching(zoneKey_).used()),
-            exposedPatchi
-        );
-
-        DebugInfo
-            << "Allocating subset of size " << subMeshPtr_().subMesh().nCells()
-            << " with exposed faces into patch "
-            << patches[exposedPatchi].name() << endl;
-    }
-
-    // Select either the submesh or the underlying mesh
-    const fvMesh& mesh =
-    (
-        subMeshPtr_.valid()
-      ? subMeshPtr_().subMesh()
-      : static_cast<const fvMesh&>(this->mesh())
-    );
-
-
-    // Distance to cell centres
-    // ~~~~~~~~~~~~~~~~~~~~~~~~
-
-    cellDistancePtr_.reset
-    (
-        new volScalarField
-        (
-            IOobject
-            (
-                "cellDistance",
-                mesh.time().timeName(),
-                mesh.time(),
-                IOobject::NO_READ,
-                IOobject::NO_WRITE,
-                false
-            ),
-            mesh,
-            dimensionedScalar(dimLength, 0)
-        )
-    );
-    volScalarField& cellDistance = cellDistancePtr_();
-
-    // Internal field
-    {
-        const pointField& cc = mesh.cellCentres();
-        scalarField& fld = cellDistance.primitiveFieldRef();
-
-        forAll(cc, i)
-        {
-            // Signed distance
-            fld[i] = (cc[i] - plane_.refPoint()) & plane_.normal();
-        }
-    }
-
-    volScalarField::Boundary& cellDistanceBf =
-        cellDistance.boundaryFieldRef();
-
-    // Patch fields
-    {
-        forAll(cellDistanceBf, patchi)
-        {
-            if
-            (
-                isA<emptyFvPatchScalarField>
-                (
-                    cellDistanceBf[patchi]
-                )
-            )
-            {
-                cellDistanceBf.set
-                (
-                    patchi,
-                    new calculatedFvPatchScalarField
-                    (
-                        mesh.boundary()[patchi],
-                        cellDistance
-                    )
-                );
-
-                const polyPatch& pp = mesh.boundary()[patchi].patch();
-                pointField::subField cc = pp.patchSlice(mesh.faceCentres());
-
-                fvPatchScalarField& fld = cellDistanceBf[patchi];
-                fld.setSize(pp.size());
-                forAll(fld, i)
-                {
-                    fld[i] = (cc[i] - plane_.refPoint()) & plane_.normal();
-                }
-            }
-            else
-            {
-                const pointField& cc = mesh.C().boundaryField()[patchi];
-                fvPatchScalarField& fld = cellDistanceBf[patchi];
-
-                forAll(fld, i)
-                {
-                    fld[i] = (cc[i] - plane_.refPoint()) & plane_.normal();
-                }
-            }
-        }
-    }
-
-
-    // On processor patches the mesh.C() will already be the cell centre
-    // on the opposite side so no need to swap cellDistance.
-
-
-    // Distance to points
-    pointDistance_.setSize(mesh.nPoints());
-    {
-        const pointField& pts = mesh.points();
-
-        forAll(pointDistance_, i)
-        {
-            pointDistance_[i] = (pts[i] - plane_.refPoint()) & plane_.normal();
-        }
-    }
-
-
-    if (debug)
-    {
-        Pout<< "Writing cell distance:" << cellDistance.objectPath() << endl;
-        cellDistance.write();
-        pointScalarField pDist
-        (
-            IOobject
-            (
-                "pointDistance",
-                mesh.time().timeName(),
-                mesh.time(),
-                IOobject::NO_READ,
-                IOobject::NO_WRITE,
-                false
-            ),
-            pointMesh::New(mesh),
-            dimensionedScalar(dimLength, 0)
-        );
-        pDist.primitiveFieldRef() = pointDistance_;
-
-        Pout<< "Writing point distance:" << pDist.objectPath() << endl;
-        pDist.write();
-    }
-
-
-    //- Direct from cell field and point field.
-    isoSurfPtr_.reset
-    (
-        new isoSurface
-        (
-            mesh,
-            cellDistance,
-            pointDistance_,
-            0,
-            filter_
-        )
-    );
-
-    if (debug)
-    {
-        print(Pout);
-        Pout<< endl;
-    }
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::sampledSurfaces::cuttingPlane::cuttingPlane
-(
-    const word& name,
-    const polyMesh& mesh,
-    const dictionary& dict
-)
-:
-    sampledSurface(name, mesh, dict),
-    plane_(dict),
-    filter_
-    (
-        dict.found("filtering")
-      ? isoSurface::filterTypeNames_.read(dict.lookup("filtering"))
-      : isoSurface::filterType::full
-    ),
-    average_(dict.lookupOrDefault("average", false)),
-    zoneKey_(dict.lookupOrDefault("zone", wordRe::null)),
-    exposedPatchName_(word::null),
-    needsUpdate_(true),
-    subMeshPtr_(nullptr),
-    cellDistancePtr_(nullptr),
-    isoSurfPtr_(nullptr)
-{
-    if (zoneKey_.size())
-    {
-        dict.lookup("exposedPatchName") >> exposedPatchName_;
-
-        if (mesh.boundaryMesh().findPatchID(exposedPatchName_) == -1)
-        {
-            FatalErrorInFunction
-                << "Cannot find patch " << exposedPatchName_
-                << " in which to put exposed faces." << endl
-                << "Valid patches are " << mesh.boundaryMesh().names()
-                << exit(FatalError);
-        }
-
-        if (debug)
-        {
-            Info<< "Restricting to cellZone " << zoneKey_
-                << " with exposed internal faces into patch "
-                << exposedPatchName_ << endl;
-        }
-
-        if (mesh.cellZones().findIndex(zoneKey_) < 0)
-        {
-            WarningInFunction
-                << "cellZone " << zoneKey_
-                << " not found - using entire mesh" << endl;
-        }
-    }
-}
-
-
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-Foam::sampledSurfaces::cuttingPlane::~cuttingPlane()
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-bool Foam::sampledSurfaces::cuttingPlane::needsUpdate() const
-{
-    return needsUpdate_;
-}
-
-
-bool Foam::sampledSurfaces::cuttingPlane::expire()
-{
-    if (debug)
-    {
-        Pout<< "cuttingPlane::expire :"
-            << " needsUpdate_:" << needsUpdate_ << endl;
-    }
-
-    // Clear derived data
-    clearGeom();
-
-    // Already marked as expired
-    if (needsUpdate_)
-    {
-        return false;
-    }
-
-    needsUpdate_ = true;
-    return true;
-}
-
-
-bool Foam::sampledSurfaces::cuttingPlane::update()
-{
-    if (debug)
-    {
-        Pout<< "cuttingPlane::update :"
-            << " needsUpdate_:" << needsUpdate_ << endl;
-    }
-
-    if (!needsUpdate_)
-    {
-        return false;
-    }
-
-    createGeometry();
-
-    needsUpdate_ = false;
-    return true;
-}
-
-
-Foam::tmp<Foam::scalarField>
-Foam::sampledSurfaces::cuttingPlane::sample
-(
-    const volScalarField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::vectorField>
-Foam::sampledSurfaces::cuttingPlane::sample
-(
-    const volVectorField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::sphericalTensorField>
-Foam::sampledSurfaces::cuttingPlane::sample
-(
-    const volSphericalTensorField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::symmTensorField>
-Foam::sampledSurfaces::cuttingPlane::sample
-(
-    const volSymmTensorField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::tensorField>
-Foam::sampledSurfaces::cuttingPlane::sample
-(
-    const volTensorField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::scalarField>
-Foam::sampledSurfaces::cuttingPlane::interpolate
-(
-    const interpolation<scalar>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-
-Foam::tmp<Foam::vectorField>
-Foam::sampledSurfaces::cuttingPlane::interpolate
-(
-    const interpolation<vector>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-Foam::tmp<Foam::sphericalTensorField>
-Foam::sampledSurfaces::cuttingPlane::interpolate
-(
-    const interpolation<sphericalTensor>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-
-Foam::tmp<Foam::symmTensorField>
-Foam::sampledSurfaces::cuttingPlane::interpolate
-(
-    const interpolation<symmTensor>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-
-Foam::tmp<Foam::tensorField>
-Foam::sampledSurfaces::cuttingPlane::interpolate
-(
-    const interpolation<tensor>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-
-void Foam::sampledSurfaces::cuttingPlane::print(Ostream& os) const
-{
-    os  << "cuttingPlane: " << name() << " :"
-        << "  plane:" << plane_
-        << "  faces:" << faces().size()
-        << "  points:" << points().size();
-}
-
-
-// ************************************************************************* //

src/sampling/sampledSurface/sampledCuttingPlane/sampledCuttingPlane.H

@@ -1,270 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2021 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-Class
-    Foam::sampledSurfaces::cuttingPlane
-
-Description
-    A sampledSurface defined by a plane using the iso-surface algorithm
-    to 'cut' the mesh.
-
-    Example:
-    \verbatim
-    {
-        type        cuttingPlane;
-        planeType   pointAndNormal;
-        point       (0 0 0);
-        normal      (0 0 1);
-        filtering   full;
-        interpolate yes;
-    }
-    \endverbatim
-
-Usage
-    \table
-        Property     | Description             | Required    | Default value
-        planeType    | the method of specification of the plane | yes |
-        filtering    | the level of filtering to perform on the iso-surface \\
-                                               | no          | full
-        interpolate  | interpolate values to the surface points | no | no
-    \endtable
-
-See also
-    Foam::plane
-    Foam::isoSurface
-
-SourceFiles
-    sampledCuttingPlane.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef sampledCuttingPlane_H
-#define sampledCuttingPlane_H
-
-#include "sampledSurface.H"
-#include "isoSurface.H"
-#include "fvMeshSubset.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-namespace sampledSurfaces
-{
-
-/*---------------------------------------------------------------------------*\
-                        Class cuttingPlane Declaration
-\*---------------------------------------------------------------------------*/
-
-class cuttingPlane
-:
-    public sampledSurface
-{
-    // Private Data
-
-        //- Plane
-        const plane plane_;
-
-        //- Whether to coarsen
-        const isoSurface::filterType filter_;
-
-        //- Whether to recalculate cell values as average of point values
-        const Switch average_;
-
-        //- If restricted to zones, name of this zone or a regular expression
-        wordRe zoneKey_;
-
-        //- For zones: patch to put exposed faces into
-        mutable word exposedPatchName_;
-
-        //- Track if the surface needs an update
-        mutable bool needsUpdate_;
-
-        //- Optional subsetted mesh
-        autoPtr<fvMeshSubset> subMeshPtr_;
-
-        //- Distance to cell centres
-        autoPtr<volScalarField> cellDistancePtr_;
-
-        //- Distance to points
-        scalarField pointDistance_;
-
-        //- Constructed iso surface
-        autoPtr<isoSurface> isoSurfPtr_;
-
-        //- Triangles converted to faceList
-        mutable autoPtr<faceList> facesPtr_;
-
-
-    // Private Member Functions
-
-        //- Create iso surface
-        void createGeometry();
-
-        //- Sample field on faces
-        template<class Type>
-        tmp<Field<Type>> sampleField
-        (
-            const GeometricField<Type, fvPatchField, volMesh>& vField
-        ) const;
-
-
-        template<class Type>
-        tmp<Field<Type>>
-        interpolateField(const interpolation<Type>&) const;
-
-
-public:
-
-    //- Runtime type information
-    TypeName("cuttingPlane");
-
-
-    // Constructors
-
-        //- Construct from dictionary
-        cuttingPlane
-        (
-            const word& name,
-            const polyMesh& mesh,
-            const dictionary& dict
-        );
-
-
-    //- Destructor
-    virtual ~cuttingPlane();
-
-
-    // Member Functions
-
-        //- Does the surface need an update?
-        virtual bool needsUpdate() const;
-
-        //- Mark the surface as needing an update.
-        //  May also free up unneeded data.
-        //  Return false if surface was already marked as expired.
-        virtual bool expire();
-
-        //- Update the surface as required.
-        //  Do nothing (and return false) if no update was needed
-        virtual bool update();
-
-        //- Points of surface
-        virtual const pointField& points() const
-        {
-            return surface().points();
-        }
-
-        //- Faces of surface
-        virtual const faceList& faces() const
-        {
-            return surface().faces();
-        }
-
-        const isoSurface& surface() const
-        {
-            return isoSurfPtr_();
-        }
-
-        //- Sample field on surface
-        virtual tmp<scalarField> sample
-        (
-            const volScalarField&
-        ) const;
-
-        //- Sample field on surface
-        virtual tmp<vectorField> sample
-        (
-            const volVectorField&
-        ) const;
-
-        //- Sample field on surface
-        virtual tmp<sphericalTensorField> sample
-        (
-            const volSphericalTensorField&
-        ) const;
-
-        //- Sample field on surface
-        virtual tmp<symmTensorField> sample
-        (
-            const volSymmTensorField&
-        ) const;
-
-        //- Sample field on surface
-        virtual tmp<tensorField> sample
-        (
-            const volTensorField&
-        ) const;
-
-
-        //- Interpolate field on surface
-        virtual tmp<scalarField> interpolate
-        (
-            const interpolation<scalar>&
-        ) const;
-
-        //- Interpolate field on surface
-        virtual tmp<vectorField> interpolate
-        (
-            const interpolation<vector>&
-        ) const;
-
-        //- Interpolate field on surface
-        virtual tmp<sphericalTensorField> interpolate
-        (
-            const interpolation<sphericalTensor>&
-        ) const;
-
-        //- Interpolate field on surface
-        virtual tmp<symmTensorField> interpolate
-        (
-            const interpolation<symmTensor>&
-        ) const;
-
-        //- Interpolate field on surface
-        virtual tmp<tensorField> interpolate
-        (
-            const interpolation<tensor>&
-        ) const;
-
-        //- Write
-        virtual void print(Ostream&) const;
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace sampledSurfaces
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#ifdef NoRepository
-    #include "sampledCuttingPlaneTemplates.C"
-#endif
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //

src/sampling/sampledSurface/sampledCuttingPlane/sampledCuttingPlaneTemplates.C

@@ -1,94 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2020 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "sampledCuttingPlane.H"
-#include "volPointInterpolation.H"
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-template<class Type>
-Foam::tmp<Foam::Field<Type>>
-Foam::sampledSurfaces::cuttingPlane::sampleField
-(
-    const GeometricField<Type, fvPatchField, volMesh>& vField
-) const
-{
-    return tmp<Field<Type>>(new Field<Type>(vField, surface().meshCells()));
-}
-
-
-template<class Type>
-Foam::tmp<Foam::Field<Type>>
-Foam::sampledSurfaces::cuttingPlane::interpolateField
-(
-    const interpolation<Type>& interpolator
-) const
-{
-    // Get fields to sample. Assume volPointInterpolation!
-    const GeometricField<Type, fvPatchField, volMesh>& volFld =
-        interpolator.psi();
-
-    if (subMeshPtr_.valid())
-    {
-        tmp<GeometricField<Type, fvPatchField, volMesh>> tvolSubFld =
-            subMeshPtr_().interpolate(volFld);
-
-        const GeometricField<Type, fvPatchField, volMesh>& volSubFld =
-            tvolSubFld();
-
-        tmp<GeometricField<Type, pointPatchField, pointMesh>> tpointSubFld =
-            volPointInterpolation::New(volSubFld.mesh()).interpolate(volSubFld);
-
-        return surface().interpolate
-        (
-            (
-                average_
-              ? pointAverage(tpointSubFld())()
-              : volSubFld
-            ),
-            tpointSubFld()
-        );
-    }
-    else
-    {
-        tmp<GeometricField<Type, pointPatchField, pointMesh>> tpointFld
-        (
-            volPointInterpolation::New(volFld.mesh()).interpolate(volFld)
-        );
-
-        return surface().interpolate
-        (
-            (
-                average_
-              ? pointAverage(tpointFld())()
-              : volFld
-            ),
-            tpointFld()
-        );
-    }
-}
-
-
-// ************************************************************************* //

src/sampling/sampledSurface/sampledDistanceSurface/sampledDistanceSurface.C

@@ -0,0 +1,154 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "sampledDistanceSurface.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace sampledSurfaces
+{
+    defineTypeNameAndDebug(distanceSurface, 0);
+    addToRunTimeSelectionTable(sampledSurface, distanceSurface, word);
+}
+}
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+Foam::autoPtr<Foam::cutPolyIsoSurface>
+Foam::sampledSurfaces::distanceSurface::calcIsoSurf() const
+{
+    // Compute the distance from the mesh points to the surface
+    scalarField pointDistance(mesh().nPoints());
+    {
+        List<pointIndexHit> nearest;
+        surfPtr_().findNearest
+        (
+            mesh().points(),
+            scalarField(mesh().points().size(), great),
+            nearest
+        );
+
+        if (signed_)
+        {
+            List<volumeType> volType;
+            surfPtr_().getVolumeType(mesh().points(), volType);
+
+            forAll(nearest, i)
+            {
+                if
+                (
+                    volType[i] != volumeType::outside
+                 && volType[i] != volumeType::inside
+                )
+                {
+                    FatalErrorInFunction
+                        << "Point " << mesh().points()[i] << " could not be "
+                        << "classified as either inside or outside the surface "
+                        << surfPtr_->name() << exit(FatalError);
+                }
+
+                pointDistance[i] =
+                    (volType[i] == volumeType::outside ? +1 : -1)
+                   *mag(mesh().points()[i] - nearest[i].hitPoint());
+            }
+        }
+        else
+        {
+            forAll(nearest, i)
+            {
+                pointDistance[i] =
+                    mag(mesh().points()[i] - nearest[i].hitPoint());
+            }
+        }
+    }
+
+    // Construct an iso-surface at the given distance
+    return autoPtr<cutPolyIsoSurface>
+    (
+        new cutPolyIsoSurface(mesh(), pointDistance, distance_, zoneIDs())
+    );
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::sampledSurfaces::distanceSurface::distanceSurface
+(
+    const word& name,
+    const polyMesh& mesh,
+    const dictionary& dict
+)
+:
+    sampledIsoSurfaceSurface(name, mesh, dict),
+    surfPtr_
+    (
+        searchableSurface::New
+        (
+            dict.lookup("surfaceType"),
+            IOobject
+            (
+                dict.lookupOrDefault("surfaceName", name),
+                mesh.time().constant(),
+                searchableSurface::geometryDir(mesh.time()),
+                mesh.time(),
+                IOobject::MUST_READ,
+                IOobject::NO_WRITE
+            ),
+            dict
+        )
+    ),
+    distance_(dict.lookup<scalar>("distance")),
+    signed_(readBool(dict.lookup("signed")))
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::sampledSurfaces::distanceSurface::~distanceSurface()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+bool Foam::sampledSurfaces::distanceSurface::needsUpdate() const
+{
+    return timeIndex() == -1;
+}
+
+
+void Foam::sampledSurfaces::distanceSurface::print(Ostream& os) const
+{
+    os  << "distanceSurface: " << name() << " :"
+        << "  surface:" << surfPtr_().name()
+        << "  distance:" << distance_
+        << "  faces:" << faces().size()
+        << "  points:" << points().size();
+}
+
+
+// ************************************************************************* //

src/sampling/sampledSurface/sampledDistanceSurface/sampledDistanceSurface.H

@@ -0,0 +1,137 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+     Foam::sampledSurfaces::distanceSurface
+
+Description
+    A sampledSurface defined by a distance to a surface.
+
+    Example:
+    \verbatim
+    {
+        type        distanceSurface;
+        surfaceType searchableBox;
+        min         (-1 -1 -1);
+        max         (1 1 1);
+        distance    0.1;
+        signed      yes;
+        interpolate yes;
+    }
+    \endverbatim
+
+Usage
+    \table
+        Property     | Description             | Required    | Default value
+        surfaceType  | the type of surface to sample from | yes |
+        distance     | the distance from which to sample the surface | yes |
+        signed       | sample only on one side of the surface, specified by \\
+                       the sign of the distance | yes |
+        interpolate  | interpolate values to the surface points | no | no
+    \endtable
+
+SourceFiles
+    sampledDistanceSurface.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef sampledDistanceSurface_H
+#define sampledDistanceSurface_H
+
+#include "sampledIsoSurfaceSurface.H"
+#include "searchableSurface.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace sampledSurfaces
+{
+
+/*---------------------------------------------------------------------------*\
+                       Class distanceSurface Declaration
+\*---------------------------------------------------------------------------*/
+
+class distanceSurface
+:
+    public sampledIsoSurfaceSurface
+{
+    // Private Data
+
+        //- Surface
+        const autoPtr<searchableSurface> surfPtr_;
+
+        //- Distance value
+        const scalar distance_;
+
+        //- Is the distance signed? The surface must be closed if so.
+        const bool signed_;
+
+
+    // Private Member Functions
+
+        //- Generate the iso surface
+        virtual autoPtr<cutPolyIsoSurface> calcIsoSurf() const;
+
+
+public:
+
+    //- Runtime type information
+    TypeName("distanceSurface");
+
+
+    // Constructors
+
+        //- Construct from dictionary
+        distanceSurface
+        (
+            const word& name,
+            const polyMesh& mesh,
+            const dictionary& dict
+        );
+
+
+    //- Destructor
+    virtual ~distanceSurface();
+
+
+    // Member Functions
+
+        //- Does the surface need an update?
+        virtual bool needsUpdate() const;
+
+        //- Print
+        virtual void print(Ostream&) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace sampledSurfaces
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/sampling/sampledSurface/sampledIsoSurface/sampledIsoSurface.C

@@ -0,0 +1,136 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "sampledIsoSurface.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace sampledSurfaces
+{
+    defineTypeNameAndDebug(isoSurface, 0);
+    addToRunTimeSelectionTable(sampledSurface, isoSurface, word);
+}
+}
+
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+Foam::autoPtr<Foam::cutPolyIsoSurface>
+Foam::sampledSurfaces::isoSurface::calcIsoSurf() const
+{
+    // Lookup the field
+    const volScalarField& vField =
+        mesh().lookupObject<volScalarField>(isoFieldName_);
+
+    // Interpolate the field to the points
+    tmp<pointScalarField> pField
+    (
+        volPointInterpolation::New(vField.mesh()).interpolate(vField)
+    );
+
+    // Construct iso-surfaces for the given iso-values
+    PtrList<cutPolyIsoSurface> isoSurfs(isoValues_.size());
+    forAll(isoValues_, i)
+    {
+        isoSurfs.set
+        (
+            i,
+            new cutPolyIsoSurface(mesh(), pField, isoValues_[i], zoneIDs())
+        );
+    }
+
+    // Return or combine into a single surface
+    if (isoValues_.size() == 1)
+    {
+        return autoPtr<cutPolyIsoSurface>(isoSurfs.set(0, nullptr));
+    }
+    else
+    {
+        return autoPtr<cutPolyIsoSurface>(new cutPolyIsoSurface(isoSurfs));
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::sampledSurfaces::isoSurface::isoSurface
+(
+    const word& name,
+    const polyMesh& mesh,
+    const dictionary& dict
+)
+:
+    sampledIsoSurfaceSurface(name, mesh, dict),
+    isoFieldName_(dict.lookup("isoField")),
+    isoValues_
+    (
+        dict.found("isoValues")
+      ? scalarField(dict.lookup("isoValues"))
+      : scalarField(1, dict.lookup<scalar>("isoValue"))
+    )
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::sampledSurfaces::isoSurface::~isoSurface()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+Foam::wordList Foam::sampledSurfaces::isoSurface::fields() const
+{
+    return wordList(isoFieldName_);
+}
+
+
+bool Foam::sampledSurfaces::isoSurface::needsUpdate() const
+{
+    return timeIndex() != mesh().time().timeIndex();
+}
+
+
+void Foam::sampledSurfaces::isoSurface::print(Ostream& os) const
+{
+    os  << "isoSurface: " << name() << " :"
+        << "  field:" << isoFieldName_;
+    if (isoValues_.size() == 1)
+    {
+        os << "  value:" << isoValues_[0];
+    }
+    else
+    {
+        os << "  values:" << isoValues_;
+    }
+    os  << "  faces:" << faces().size()
+        << "  points:" << points().size();
+}
+
+
+// ************************************************************************* //

src/sampling/sampledSurface/sampledIsoSurface/sampledIsoSurface.H

@@ -0,0 +1,131 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::sampledSurfaces::isoSurface
+
+Description
+    A sampledSurface defined by an iso-value of a field.
+
+    Example:
+    \verbatim
+    {
+        type        isoSurface;
+        isoField    p;
+        isoValue    1e5;
+        interpolate yes;
+    }
+    \endverbatim
+
+Usage
+    \table
+        Property     | Description             | Required    | Default value
+        isoField     | the field to get an iso-surface of | yes |
+        isoValue     | the iso-surface value   | yes         |
+        interpolate  | interpolate values to the surface points | no | no
+    \endtable
+
+SourceFiles
+    sampledIsoSurface.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef sampledIsoSurface_H
+#define sampledIsoSurface_H
+
+#include "sampledIsoSurfaceSurface.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace sampledSurfaces
+{
+
+/*---------------------------------------------------------------------------*\
+                         Class isoSurface Declaration
+\*---------------------------------------------------------------------------*/
+
+class isoSurface
+:
+    public sampledIsoSurfaceSurface
+{
+    // Private Data
+
+        //- Field to get iso surface of
+        const word isoFieldName_;
+
+        //- Iso values
+        const scalarField isoValues_;
+
+
+    // Private Member Functions
+
+        //- Generate the iso surface
+        virtual autoPtr<cutPolyIsoSurface> calcIsoSurf() const;
+
+
+public:
+
+    //- Runtime type information
+    TypeName("isoSurface");
+
+
+    // Constructors
+
+        //- Construct from dictionary
+        isoSurface
+        (
+            const word& name,
+            const polyMesh& mesh,
+            const dictionary& dict
+        );
+
+
+    //- Destructor
+    virtual ~isoSurface();
+
+
+    // Member Functions
+
+        //- Return the list of fields required
+        virtual wordList fields() const;
+
+        //- Does the surface need an update?
+        virtual bool needsUpdate() const;
+
+        //- Write
+        virtual void print(Ostream&) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace sampledSurfaces
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/sampling/sampledSurface/sampledIsoSurfaceSurface/sampledIsoSurfaceSurface.C

@@ -0,0 +1,141 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "sampledIsoSurfaceSurface.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace sampledSurfaces
+{
+    defineTypeNameAndDebug(sampledIsoSurfaceSurface, 0);
+}
+}
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::sampledSurfaces::sampledIsoSurfaceSurface::sampledIsoSurfaceSurface
+(
+    const word& name,
+    const polyMesh& mesh,
+    const dictionary& dict
+)
+:
+    sampledSurface(name, mesh, dict),
+    zoneName_(dict.lookupOrDefault("zone", wordRe::null)),
+    zoneIDs_(mesh.cellZones().findIndices(zoneName_)),
+    isoSurfPtr_(nullptr),
+    isoSurfTimeIndex_(-1)
+{
+    if (zoneName_ != wordRe::null && zoneIDs_.empty())
+    {
+        WarningInFunction
+            << "Cell zone " << zoneName_
+            << " not found. Using the entire mesh" << endl;
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::sampledSurfaces::sampledIsoSurfaceSurface::~sampledIsoSurfaceSurface()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+bool Foam::sampledSurfaces::sampledIsoSurfaceSurface::expire()
+{
+    // Clear data
+    sampledSurface::clearGeom();
+    isoSurfPtr_.clear();
+
+    // Already marked as expired
+    if (isoSurfTimeIndex_ == -1)
+    {
+        return false;
+    }
+
+    // Force update
+    isoSurfTimeIndex_ = -1;
+    return true;
+}
+
+
+bool Foam::sampledSurfaces::sampledIsoSurfaceSurface::update() const
+{
+    // Quick return if no update needed
+    if (!needsUpdate())
+    {
+        return false;
+    }
+
+    // Clear any information in the base class
+    sampledSurface::clearGeom();
+
+    // Update the iso surface
+    isoSurfPtr_.reset(calcIsoSurf().ptr());
+
+    // Set the time index
+    isoSurfTimeIndex_ = mesh().time().timeIndex();
+
+    return true;
+}
+
+
+bool Foam::sampledSurfaces::sampledIsoSurfaceSurface::update()
+{
+    return static_cast<const sampledIsoSurfaceSurface&>(*this).update();
+}
+
+
+#define IMPLEMENT_SAMPLE(Type, nullArg)                                        \
+    Foam::tmp<Foam::Field<Foam::Type>>                                         \
+    Foam::sampledSurfaces::sampledIsoSurfaceSurface::sample                    \
+    (                                                                          \
+        const VolField<Type>& vField                                           \
+    ) const                                                                    \
+    {                                                                          \
+        return sampleField(vField);                                            \
+    }
+FOR_ALL_FIELD_TYPES(IMPLEMENT_SAMPLE);
+#undef IMPLEMENT_SAMPLE
+
+
+#define IMPLEMENT_INTERPOLATE(Type, nullArg)                                   \
+    Foam::tmp<Foam::Field<Foam::Type>>                                         \
+    Foam::sampledSurfaces::sampledIsoSurfaceSurface::interpolate               \
+    (                                                                          \
+        const interpolation<Type>& interpolator                                \
+    ) const                                                                    \
+    {                                                                          \
+        return interpolateField(interpolator);                                 \
+    }
+FOR_ALL_FIELD_TYPES(IMPLEMENT_INTERPOLATE);
+#undef IMPLEMENT_INTERPOLATE
+
+
+// ************************************************************************* //

src/sampling/sampledSurface/sampledIsoSurfaceSurface/sampledIsoSurfaceSurface.H

@@ -0,0 +1,190 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+     Foam::sampledSurfaces::sampledIsoSurfaceSurface
+
+Description
+    A base class for sampled surfaces constructed from iso-surfaces
+
+SourceFiles
+    sampledIsoSurfaceSurface.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef sampledIsoSurfaceSurface_H
+#define sampledIsoSurfaceSurface_H
+
+#include "sampledSurface.H"
+#include "cutPolyIsoSurface.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace sampledSurfaces
+{
+
+/*---------------------------------------------------------------------------*\
+                       Class sampledIsoSurfaceSurface Declaration
+\*---------------------------------------------------------------------------*/
+
+class sampledIsoSurfaceSurface
+:
+    public sampledSurface
+{
+private:
+
+    // Private Data
+
+        //- If restricted to zones, name of this zone or a regular expression
+        const wordRe zoneName_;
+
+        //- If restricted to zones, the indices of the zones
+        const labelList zoneIDs_;
+
+        //- Constructed iso surface
+        mutable autoPtr<cutPolyIsoSurface> isoSurfPtr_;
+
+        //- Time index at last surface generation
+        mutable label isoSurfTimeIndex_;
+
+
+    // Private Member Functions
+
+        //- Generate the iso surface
+        virtual autoPtr<cutPolyIsoSurface> calcIsoSurf() const = 0;
+
+        //- Sample field on the surface's faces
+        template<class Type>
+        tmp<Field<Type>> sampleField
+        (
+            const GeometricField<Type, fvPatchField, volMesh>& vField
+        ) const;
+
+        //- Interpolate field to the surface's points
+        template<class Type>
+        tmp<Field<Type>> interpolateField(const interpolation<Type>&) const;
+
+
+protected:
+
+    // Protected Member Functions
+
+        //- Access the zone indices
+        inline const labelList& zoneIDs() const
+        {
+            return zoneIDs_.empty() ? NullObjectRef<labelList>() : zoneIDs_;
+        }
+
+        //- Access the time index
+        inline const label& timeIndex() const
+        {
+            return isoSurfTimeIndex_;
+        }
+
+
+public:
+
+    //- Runtime type information
+    TypeName("sampledIsoSurfaceSurface");
+
+
+    // Constructors
+
+        //- Construct from dictionary
+        sampledIsoSurfaceSurface
+        (
+            const word& name,
+            const polyMesh& mesh,
+            const dictionary& dict
+        );
+
+
+    //- Destructor
+    virtual ~sampledIsoSurfaceSurface();
+
+
+    // Member Functions
+
+        //- Mark the surface as needing an update.
+        //  May also free up unneeded data.
+        //  Return false if surface was already marked as expired.
+        virtual bool expire();
+
+        //- Update the surface as required.
+        //  Do nothing (and return false) if no update was needed
+        virtual bool update();
+
+        //- Update the surface as required.
+        //  Do nothing (and return false) if no update was needed
+        virtual bool update() const;
+
+        //- Points of surface
+        virtual const pointField& points() const
+        {
+            return isoSurfPtr_->points();
+        }
+
+        //- Faces of surface
+        virtual const faceList& faces() const
+        {
+            return isoSurfPtr_->faces();
+        }
+
+        //- Sample field on the surface's faces
+        #define DEFINE_SAMPLE(Type, nullArg)                                   \
+            virtual tmp<Field<Type>> sample                                    \
+            (                                                                  \
+                const VolField<Type>&                                          \
+            ) const;
+        FOR_ALL_FIELD_TYPES(DEFINE_SAMPLE);
+        #undef DEFINE_SAMPLE
+
+        //- Interpolate field to the surface's points
+        #define DEFINE_INTERPOLATE(Type, nullArg)                              \
+            virtual tmp<Field<Type>> interpolate                               \
+            (                                                                  \
+                const interpolation<Type>&                                     \
+            ) const;
+        FOR_ALL_FIELD_TYPES(DEFINE_INTERPOLATE);
+        #undef DEFINE_INTERPOLATE
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace sampledSurfaces
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+    #include "sampledIsoSurfaceSurfaceTemplates.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/sampling/sampledSurface/sampledIsoSurfaceSurface/sampledIsoSurfaceSurfaceTemplates.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2019 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -23,55 +23,58 @@ License
 
 \*---------------------------------------------------------------------------*/
 
-#include "sampledPlane.H"
+#include "sampledIsoSurfaceSurface.H"
+#include "interpolationVolPointInterpolation.H"
 
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
 
 template<class Type>
 Foam::tmp<Foam::Field<Type>>
-Foam::sampledSurfaces::plane::sampleField
+Foam::sampledSurfaces::sampledIsoSurfaceSurface::sampleField
 (
     const GeometricField<Type, fvPatchField, volMesh>& vField
 ) const
 {
-    return tmp<Field<Type>>(new Field<Type>(vField, meshCells()));
+    update();
+
+    return isoSurfPtr_->sample(vField.primitiveField());
 }
 
 
 template<class Type>
 Foam::tmp<Foam::Field<Type>>
-Foam::sampledSurfaces::plane::interpolateField
+Foam::sampledSurfaces::sampledIsoSurfaceSurface::interpolateField
 (
     const interpolation<Type>& interpolator
 ) const
 {
-    // One value per point
-    tmp<Field<Type>> tvalues(new Field<Type>(points().size()));
-    Field<Type>& values = tvalues.ref();
+    update();
 
-    boolList pointDone(points().size(), false);
-
-    forAll(faces(), cutFacei)
+    if (isA<interpolationVolPointInterpolation<Type>>(interpolator))
     {
-        const face& f = faces()[cutFacei];
+        const GeometricField<Type, pointPatchField, pointMesh>& pField =
+            refCast<const interpolationVolPointInterpolation<Type>>
+            (interpolator).psip();
 
-        forAll(f, faceVertI)
+        return isoSurfPtr_().interpolate(pField.primitiveField());
+    }
+    else
+    {
+        const pointField& points = isoSurfPtr_->points();
+        const faceList& faces = isoSurfPtr_->faces();
+        const labelList& faceCells = isoSurfPtr_->faceCells();
+
+        labelList pointCells(points.size());
+        forAll(faces, facei)
         {
-            label pointi = f[faceVertI];
-
-            if (!pointDone[pointi])
+            forAll(faces[facei], facePointi)
             {
-                values[pointi] = interpolator.interpolate
-                (
-                    points()[pointi],
-                    meshCells()[cutFacei]
-                );
-                pointDone[pointi] = true;
+                pointCells[faces[facei][facePointi]] = faceCells[facei];
             }
         }
-    }
 
-    return tvalues;
+        return interpolator.interpolate(points, pointCells);
+    }
 }
 
 

src/sampling/sampledSurface/sampledPlane/sampledPlane.C

@@ -1,243 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2020 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "sampledPlane.H"
-#include "addToRunTimeSelectionTable.H"
-
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-namespace Foam
-{
-namespace sampledSurfaces
-{
-    defineTypeNameAndDebug(plane, 0);
-    addToRunTimeSelectionTable(sampledSurface, plane, word);
-}
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::sampledSurfaces::plane::plane
-(
-    const word& name,
-    const polyMesh& mesh,
-    const dictionary& dict
-)
-:
-    sampledSurface(name, mesh, dict),
-    cuttingPlane(Foam::plane(dict)),
-    zoneKey_(dict.lookupOrDefault("zone", wordRe::null)),
-    triangulate_(dict.lookupOrDefault("triangulate", true)),
-    needsUpdate_(true)
-{
-    // Make plane relative to the coordinateSystem (Cartesian)
-    // allow lookup from global coordinate systems
-    if (dict.found("coordinateSystem"))
-    {
-        autoPtr<coordinateSystem> cs
-        (
-            coordinateSystem::New(mesh, dict.subDict("coordinateSystem"))
-        );
-
-        point  base = cs->globalPosition(planeDesc().refPoint());
-        vector norm = cs->globalVector(planeDesc().normal());
-
-        // Assign the plane description
-        static_cast<Foam::plane&>(*this) = Foam::plane(base, norm);
-    }
-
-    if (debug && zoneKey_.size() && mesh.cellZones().findIndex(zoneKey_) < 0)
-    {
-        Info<< "cellZone " << zoneKey_
-            << " not found - using entire mesh" << endl;
-    }
-}
-
-
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-Foam::sampledSurfaces::plane::~plane()
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-bool Foam::sampledSurfaces::plane::needsUpdate() const
-{
-    return needsUpdate_;
-}
-
-
-bool Foam::sampledSurfaces::plane::expire()
-{
-    // Already marked as expired
-    if (needsUpdate_)
-    {
-        return false;
-    }
-
-    sampledSurface::clearGeom();
-
-    needsUpdate_ = true;
-    return true;
-}
-
-
-bool Foam::sampledSurfaces::plane::update()
-{
-    if (!needsUpdate_)
-    {
-        return false;
-    }
-
-    sampledSurface::clearGeom();
-
-    const labelList selectedCells
-    (
-        mesh().cellZones().findMatching(zoneKey_).used()
-    );
-
-    if (selectedCells.empty())
-    {
-        reCut(mesh(), triangulate_);
-    }
-    else
-    {
-        reCut(mesh(), triangulate_, selectedCells);
-    }
-
-    if (debug)
-    {
-        print(Pout);
-        Pout<< endl;
-    }
-
-    needsUpdate_ = false;
-    return true;
-}
-
-
-Foam::tmp<Foam::scalarField> Foam::sampledSurfaces::plane::sample
-(
-    const volScalarField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::vectorField> Foam::sampledSurfaces::plane::sample
-(
-    const volVectorField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::sphericalTensorField> Foam::sampledSurfaces::plane::sample
-(
-    const volSphericalTensorField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::symmTensorField> Foam::sampledSurfaces::plane::sample
-(
-    const volSymmTensorField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::tensorField> Foam::sampledSurfaces::plane::sample
-(
-    const volTensorField& vField
-) const
-{
-    return sampleField(vField);
-}
-
-
-Foam::tmp<Foam::scalarField> Foam::sampledSurfaces::plane::interpolate
-(
-    const interpolation<scalar>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-
-Foam::tmp<Foam::vectorField> Foam::sampledSurfaces::plane::interpolate
-(
-    const interpolation<vector>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-Foam::tmp<Foam::sphericalTensorField> Foam::sampledSurfaces::plane::interpolate
-(
-    const interpolation<sphericalTensor>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-
-Foam::tmp<Foam::symmTensorField> Foam::sampledSurfaces::plane::interpolate
-(
-    const interpolation<symmTensor>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-
-Foam::tmp<Foam::tensorField> Foam::sampledSurfaces::plane::interpolate
-(
-    const interpolation<tensor>& interpolator
-) const
-{
-    return interpolateField(interpolator);
-}
-
-
-void Foam::sampledSurfaces::plane::print(Ostream& os) const
-{
-    os  << "plane: " << name() << " :"
-        << "  base:" << refPoint()
-        << "  normal:" << normal()
-        << "  triangulate:" << triangulate_
-        << "  faces:" << faces().size()
-        << "  points:" << points().size();
-}
-
-
-// ************************************************************************* //

src/sampling/sampledSurface/sampledPlane/sampledPlane.H

@@ -1,243 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2021 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-Class
-    Foam::sampledSurfaces::plane
-
-Description
-    A sampledSurface defined by a plane which 'cuts' the mesh using the
-    cuttingPlane algorithm. The plane is triangulated by default.
-
-    Example:
-    \verbatim
-    {
-        type        plane;
-        planeType   pointAndNormal;
-        point       (0 0 0);
-        normal      (0 0 1);
-        interpolate yes;
-    }
-    \endverbatim
-
-Usage
-    \table
-        Property     | Description             | Required    | Default value
-        planeType    | the method of specification of the plane | yes |
-        triangulate  | triangulate the output  | no          | yes
-        interpolate  | interpolate values to the surface points | no | no
-    \endtable
-
-See also
-    Foam::plane
-
-SourceFiles
-    sampledPlane.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef sampledPlane_H
-#define sampledPlane_H
-
-#include "sampledSurface.H"
-#include "cuttingPlane.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-namespace sampledSurfaces
-{
-
-/*---------------------------------------------------------------------------*\
-                            Class plane Declaration
-\*---------------------------------------------------------------------------*/
-
-class plane
-:
-    public sampledSurface,
-    public cuttingPlane
-{
-    // Private Data
-
-        //- If restricted to zones, name of this zone or a regular expression
-        wordRe zoneKey_;
-
-        //- Triangulated faces or keep faces as is
-        const bool triangulate_;
-
-        //- Track if the surface needs an update
-        mutable bool needsUpdate_;
-
-    // Private Member Functions
-
-        //- Sample field on faces
-        template<class Type>
-        tmp<Field<Type>> sampleField
-        (
-            const GeometricField<Type, fvPatchField, volMesh>& vField
-        ) const;
-
-
-        template<class Type>
-        tmp<Field<Type>>
-        interpolateField(const interpolation<Type>&) const;
-
-
-public:
-
-    //- Runtime type information
-    TypeName("plane");
-
-
-    // Constructors
-
-        //- Construct from dictionary
-        plane
-        (
-            const word& name,
-            const polyMesh& mesh,
-            const dictionary& dict
-        );
-
-
-    //- Destructor
-    virtual ~plane();
-
-
-    // Member Functions
-
-        //- Does the surface need an update?
-        virtual bool needsUpdate() const;
-
-        //- Mark the surface as needing an update.
-        //  May also free up unneeded data.
-        //  Return false if surface was already marked as expired.
-        virtual bool expire();
-
-        //- Update the surface as required.
-        //  Do nothing (and return false) if no update was needed
-        virtual bool update();
-
-
-        //- Points of surface
-        virtual const pointField& points() const
-        {
-            return cuttingPlane::points();
-        }
-
-        //- Faces of surface
-        virtual const faceList& faces() const
-        {
-            return cuttingPlane::faces();
-        }
-
-        //- For every face original cell in mesh
-        const labelList& meshCells() const
-        {
-            return cuttingPlane::cutCells();
-        }
-
-        //- Sample field on surface
-        virtual tmp<scalarField> sample
-        (
-            const volScalarField&
-        ) const;
-
-
-        //- Sample field on surface
-        virtual tmp<vectorField> sample
-        (
-            const volVectorField&
-        ) const;
-
-        //- Sample field on surface
-        virtual tmp<sphericalTensorField> sample
-        (
-            const volSphericalTensorField&
-        ) const;
-
-        //- Sample field on surface
-        virtual tmp<symmTensorField> sample
-        (
-            const volSymmTensorField&
-        ) const;
-
-        //- Sample field on surface
-        virtual tmp<tensorField> sample
-        (
-            const volTensorField&
-        ) const;
-
-
-        //- Interpolate field on surface
-        virtual tmp<scalarField> interpolate
-        (
-            const interpolation<scalar>&
-        ) const;
-
-
-        //- Interpolate field on surface
-        virtual tmp<vectorField> interpolate
-        (
-            const interpolation<vector>&
-        ) const;
-
-        //- Interpolate field on surface
-        virtual tmp<sphericalTensorField> interpolate
-        (
-            const interpolation<sphericalTensor>&
-        ) const;
-
-        //- Interpolate field on surface
-        virtual tmp<symmTensorField> interpolate
-        (
-            const interpolation<symmTensor>&
-        ) const;
-
-        //- Interpolate field on surface
-        virtual tmp<tensorField> interpolate
-        (
-            const interpolation<tensor>&
-        ) const;
-
-        //- Write
-        virtual void print(Ostream&) const;
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace sampledSurfaces
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#ifdef NoRepository
-    #include "sampledPlaneTemplates.C"
-#endif
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //

src/sampling/sampledSurface/sampledSurface/sampledSurface.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2021 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -330,13 +330,6 @@ public:
         //- Project field onto surface
         tmp<Field<vector>> project(const Field<tensor>&) const;
 
-        //- Interpolate from points to cell centre
-        template<class Type>
-        tmp<GeometricField<Type, fvPatchField, volMesh>> pointAverage
-        (
-            const GeometricField<Type, pointPatchField, pointMesh>& pfld
-        ) const;
-
         //- Sample field on surface
         virtual tmp<scalarField> sample
         (

src/sampling/sampledSurface/sampledSurface/sampledSurfaceTemplates.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2019 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -153,50 +153,4 @@ Foam::sampledSurface::project
 }
 
 
-template<class Type>
-Foam::tmp<Foam::GeometricField<Type, Foam::fvPatchField, Foam::volMesh>>
-Foam::sampledSurface::pointAverage
-(
-    const GeometricField<Type, pointPatchField, pointMesh>& pfld
-) const
-{
-    const fvMesh& mesh = dynamic_cast<const fvMesh&>(pfld.mesh()());
-
-    tmp<GeometricField<Type, fvPatchField, volMesh>> tcellAvg
-    (
-        GeometricField<Type, fvPatchField, volMesh>::New
-        (
-            "cellAvg",
-            mesh,
-            dimensioned<Type>("zero", dimless, Zero)
-        )
-    );
-    GeometricField<Type, fvPatchField, volMesh>& cellAvg = tcellAvg.ref();
-
-    labelField nPointCells(mesh.nCells(), 0);
-    {
-        for (label pointi = 0; pointi < mesh.nPoints(); pointi++)
-        {
-            const labelList& pCells = mesh.pointCells(pointi);
-
-            forAll(pCells, i)
-            {
-                label celli = pCells[i];
-
-                cellAvg[celli] += pfld[pointi];
-                nPointCells[celli]++;
-            }
-        }
-    }
-    forAll(cellAvg, celli)
-    {
-        cellAvg[celli] /= nPointCells[celli];
-    }
-    // Give value to calculatedFvPatchFields
-    cellAvg.correctBoundaryConditions();
-
-    return tcellAvg;
-}
-
-
 // ************************************************************************* //

src/sampling/sampledSurface/sampledThresholdCellFaces/sampledThresholdCellFaces.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2021 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

src/sampling/sampledSurface/sampledThresholdCellFaces/sampledThresholdCellFaces.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2021 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

src/sampling/sampledSurface/sampledThresholdCellFaces/sampledThresholdCellFacesTemplates.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2019 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

src/sampling/sampledSurface/sampledThresholdCellFaces/thresholdCellFaces.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Copyright (C) 2011-2020 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2022 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License

tutorials/compressible/rhoCentralFoam/movingCone/system/cuttingPlane

@@ -1,35 +0,0 @@
-/*--------------------------------*- C++ -*----------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     | Website:  https://openfoam.org
-    \\  /    A nd           | Version:  dev
-     \\/     M anipulation  |
-\*---------------------------------------------------------------------------*/
-
-cuttingPlane
-{
-    type            surfaces;
-    libs            ("libsampling.so");
-
-    writeControl    writeTime;
-
-    surfaceFormat   vtk;
-    fields          ( p U );
-
-    interpolationScheme cellPoint;
-
-    surfaces
-    (
-        zNormal
-        {
-            type            cuttingPlane;
-            planeType       pointAndNormal;
-            point           (0 0 0);
-            normal          (0 0 1);
-            interpolate     true;
-        }
-    );
-}
-
-
-// ************************************************************************* //

tutorials/modules/incompressibleFluid/motorBike/lesFiles/controlDict

@@ -66,7 +66,7 @@ functions
         );
     }
 
-    #include "cuttingPlane"
+    #include "cutPlane"
     #include "streamlines"
     #include "forceCoeffs"
 }

tutorials/modules/incompressibleFluid/motorBike/motorBike/system/controlDict

@@ -48,7 +48,7 @@ runTimeModifiable true;
 
 functions
 {
-    #include "cuttingPlane"
+    #include "cutPlane"
     #include "streamlines"
     #include "forceCoeffs"
 }

tutorials/modules/incompressibleFluid/motorBike/motorBike/system/cutPlane

@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 \*---------------------------------------------------------------------------*/
 
-cuttingPlane
+cutPlane
 {
     type            surfaces;
     libs            ("libsampling.so");
@@ -22,7 +22,7 @@ cuttingPlane
     (
         yNormal
         {
-            type            cuttingPlane;
+            type            cutPlane;
             planeType       pointAndNormal;
             point           (0 0 0);
             normal          (0 1 0);

tutorials/modules/incompressibleFluid/motorBikeSteady/system/controlDict

@@ -47,8 +47,8 @@ runTimeModifiable true;
 
 functions
 {
+    #include "cutPlane"
     #include "streamlines"
-    #include "cuttingPlane"
     #include "forceCoeffs"
 }
 

tutorials/modules/incompressibleFluid/motorBikeSteady/system/cutPlane

@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 \*---------------------------------------------------------------------------*/
 
-cuttingPlane
+cutPlane
 {
     type            surfaces;
     libs            ("libsampling.so");
@@ -21,7 +21,7 @@ cuttingPlane
     (
         yNormal
         {
-            type            cuttingPlane;
+            type            cutPlane;
             planeType       pointAndNormal;
             point           (0 0 0);
             normal          (0 1 0);

tutorials/modules/incompressibleFluid/movingCone/system/controlDict

@@ -52,7 +52,7 @@ maxCo           0.5;
 
 functions
 {
-    #include "cuttingPlane"
+    #include "cutPlane"
 }
 
 // ************************************************************************* //

tutorials/modules/incompressibleFluid/movingCone/system/cutPlane

@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 \*---------------------------------------------------------------------------*/
 
-cuttingPlane
+cutPlane
 {
     type            surfaces;
     libs            ("libsampling.so");
@@ -23,7 +23,7 @@ cuttingPlane
     (
         zNormal
         {
-            type            cuttingPlane;
+            type            cutPlane;
             planeType       pointAndNormal;
             point           (0 0 0);
             normal          (0 0 1);

tutorials/modules/incompressibleFluid/propeller/system/surfaces

@@ -22,7 +22,7 @@ surfaces
     (
         zNormal
         {
-            type            cuttingPlane;
+            type            cutPlane;
             planeType       pointAndNormal;
             point           (0 0 0);
             normal          (0 0 1);