From 0ba5f5b8a90a8d303df48300230f0aa307d1f752 Mon Sep 17 00:00:00 2001
From: Will Bainbridge <http://cfd.direct>
Date: Thu, 17 Jun 2021 10:05:37 +0100
Subject: [PATCH] fileFormats: Added generic write functions for VTK poly data

The new write functions are currently being utilised by setSet and the
vtkSurfaceWriter, but it should eventually be possible for more examples
of VTK poly data writing to be converted to use these functions.
---
 .../mesh/manipulation/setSet/Make/files       |   2 -
 .../mesh/manipulation/setSet/setSet.C         |  47 +--
 .../mesh/manipulation/setSet/writePatch.C     | 119 -------
 .../mesh/manipulation/setSet/writePatch.H     |  68 ----
 .../mesh/manipulation/setSet/writePointSet.C  | 109 -------
 .../mesh/manipulation/setSet/writePointSet.H  |  67 ----
 src/fileFormats/vtk/vtkWritePolyData.H        | 222 +++++++++++++
 .../vtk/vtkWritePolyDataTemplates.C           | 297 ++++++++++++++++++
 .../writers/vtk/vtkSurfaceWriter.C            | 131 ++------
 .../writers/vtk/vtkSurfaceWriter.H            |  10 +-
 10 files changed, 574 insertions(+), 498 deletions(-)
 delete mode 100644 applications/utilities/mesh/manipulation/setSet/writePatch.C
 delete mode 100644 applications/utilities/mesh/manipulation/setSet/writePatch.H
 delete mode 100644 applications/utilities/mesh/manipulation/setSet/writePointSet.C
 delete mode 100644 applications/utilities/mesh/manipulation/setSet/writePointSet.H
 create mode 100644 src/fileFormats/vtk/vtkWritePolyData.H
 create mode 100644 src/fileFormats/vtk/vtkWritePolyDataTemplates.C

diff --git a/applications/utilities/mesh/manipulation/setSet/Make/files b/applications/utilities/mesh/manipulation/setSet/Make/files
index a14c37f5ee..75c578a1e5 100644
--- a/applications/utilities/mesh/manipulation/setSet/Make/files
+++ b/applications/utilities/mesh/manipulation/setSet/Make/files
@@ -1,5 +1,3 @@
-writePointSet.C
-writePatch.C
 setSet.C
 
 EXE = $(FOAM_APPBIN)/setSet
diff --git a/applications/utilities/mesh/manipulation/setSet/setSet.C b/applications/utilities/mesh/manipulation/setSet/setSet.C
index e63a8b03e3..861e754e7b 100644
--- a/applications/utilities/mesh/manipulation/setSet/setSet.C
+++ b/applications/utilities/mesh/manipulation/setSet/setSet.C
@@ -41,8 +41,7 @@ Description
 #include "OFstream.H"
 #include "IFstream.H"
 #include "demandDrivenData.H"
-#include "writePatch.H"
-#include "writePointSet.H"
+#include "vtkWritePolyData.H"
 #include "IOobjectList.H"
 #include "cellZoneSet.H"
 #include "faceZoneSet.H"
@@ -80,7 +79,7 @@ void writeVTK
         // Faces of set with OpenFOAM faceID as value
 
         faceList setFaces(currentSet.size());
-        labelList faceValues(currentSet.size());
+        labelField faceValues(currentSet.size());
         label setFacei = 0;
 
         forAllConstIter(topoSet, currentSet, iter)
@@ -92,14 +91,18 @@ void writeVTK
 
         primitiveFacePatch fp(setFaces, mesh.points());
 
-        vtkWriteOps::writePatch
+        vtkWritePolyData::write
         (
-            true,
+            mesh.time().path()/vtkName,
             currentSet.name(),
-            fp,
+            true,
+            fp.localPoints(),
+            labelList(),
+            edgeList(),
+            fp.localFaces(),
             "faceID",
-            faceValues,
-            mesh.time().path()/vtkName
+            false,
+            faceValues
         );
     }
     else if (isA<cellSet>(currentSet))
@@ -140,7 +143,7 @@ void writeVTK
         }
 
         faceList setFaces(cellFaces.size());
-        labelList faceValues(cellFaces.size());
+        labelField faceValues(cellFaces.size());
         label setFacei = 0;
 
         forAllConstIter(Map<label>, cellFaces, iter)
@@ -152,24 +155,32 @@ void writeVTK
 
         primitiveFacePatch fp(setFaces, mesh.points());
 
-        vtkWriteOps::writePatch
+        vtkWritePolyData::write
         (
-            true,
+            mesh.time().path()/vtkName,
             currentSet.name(),
-            fp,
+            true,
+            fp.localPoints(),
+            labelList(),
+            edgeList(),
+            fp.localFaces(),
             "cellID",
-            faceValues,
-            mesh.time().path()/vtkName
+            false,
+            faceValues
         );
     }
     else if (isA<pointSet>(currentSet))
     {
-        vtkWriteOps::writePointSet
+        std::ofstream os(mesh.time().path()/vtkName);
+        vtkWritePolyData::write
         (
+            mesh.time().path()/vtkName,
+            currentSet.name(),
             true,
-            mesh,
-            currentSet,
-            mesh.time().path()/vtkName
+            pointField(mesh.points(), currentSet.toc()),
+            identity(currentSet.size()),
+            edgeList(),
+            faceList()
         );
     }
     else
diff --git a/applications/utilities/mesh/manipulation/setSet/writePatch.C b/applications/utilities/mesh/manipulation/setSet/writePatch.C
deleted file mode 100644
index 0aea556367..0000000000
--- a/applications/utilities/mesh/manipulation/setSet/writePatch.C
+++ /dev/null
@@ -1,119 +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 "writePatch.H"
-#include "OFstream.H"
-#include "vtkWriteOps.H"
-#include "primitiveFacePatch.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-void Foam::vtkWriteOps::writePatch
-(
-    const bool binary,
-    const word& setName,
-    const primitiveFacePatch& fp,
-    const word& fieldName,
-    labelList& fieldValues,
-    const fileName& fileName
-)
-{
-    std::ofstream pStream(fileName.c_str());
-
-    pStream
-        << "# vtk DataFile Version 2.0" << std::endl
-        << setName << std::endl;
-    if (binary)
-    {
-        pStream << "BINARY" << std::endl;
-    }
-    else
-    {
-        pStream << "ASCII" << std::endl;
-    }
-    pStream << "DATASET POLYDATA" << std::endl;
-
-
-    //------------------------------------------------------------------
-    //
-    // Write topology
-    //
-    //------------------------------------------------------------------
-
-
-    // Write points and faces as polygons
-
-    pStream << "POINTS " << fp.nPoints() << " float" << std::endl;
-
-    DynamicList<floatScalar> ptField(3*fp.nPoints());
-
-    vtkWriteOps::insert(fp.localPoints(), ptField);
-
-    vtkWriteOps::write(pStream, binary, ptField);
-
-
-    label nFaceVerts = 0;
-
-    forAll(fp.localFaces(), facei)
-    {
-        nFaceVerts += fp.localFaces()[facei].size() + 1;
-    }
-    pStream << "POLYGONS " << fp.size() << ' ' << nFaceVerts
-        << std::endl;
-
-
-    DynamicList<label> vertLabels(nFaceVerts);
-
-    forAll(fp.localFaces(), facei)
-    {
-        const face& f = fp.localFaces()[facei];
-
-        vertLabels.append(f.size());
-
-        vtkWriteOps::insert(f, vertLabels);
-    }
-    vtkWriteOps::write(pStream, binary, vertLabels);
-
-
-    //-----------------------------------------------------------------
-    //
-    // Write data
-    //
-    //-----------------------------------------------------------------
-
-    // Write faceID
-
-    pStream
-        << "CELL_DATA " << fp.size() << std::endl
-        << "FIELD attributes 1" << std::endl;
-
-    // Cell ids first
-    pStream << fieldName << " 1 " << fp.size() << " int" << std::endl;
-
-    vtkWriteOps::write(pStream, binary, fieldValues);
-}
-
-
-// ************************************************************************* //
diff --git a/applications/utilities/mesh/manipulation/setSet/writePatch.H b/applications/utilities/mesh/manipulation/setSet/writePatch.H
deleted file mode 100644
index 9bf616e3fc..0000000000
--- a/applications/utilities/mesh/manipulation/setSet/writePatch.H
+++ /dev/null
@@ -1,68 +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/>.
-
-InClass
-    Foam::writePatch
-
-Description
-    Write faceSet to vtk polydata file. Only one data which is original
-    faceID.
-
-SourceFiles
-    writePatch.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef writePatch_H
-#define writePatch_H
-
-#include "primitiveMesh.H"
-#include "primitiveFacePatch.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-namespace vtkWriteOps
-{
-
-//- Write patch
-void writePatch
-(
-    const bool binary,
-    const word& setName,
-    const primitiveFacePatch& fp,
-    const word& fieldName,
-    labelList& fieldValues,
-    const fileName& fileName
-);
-
-} // End namespace vtkWriteOps
-} // End namespace Foam
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/applications/utilities/mesh/manipulation/setSet/writePointSet.C b/applications/utilities/mesh/manipulation/setSet/writePointSet.C
deleted file mode 100644
index 0b3e52b08f..0000000000
--- a/applications/utilities/mesh/manipulation/setSet/writePointSet.C
+++ /dev/null
@@ -1,109 +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 "writePointSet.H"
-#include "OFstream.H"
-#include "vtkWriteOps.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-void Foam::vtkWriteOps::writePointSet
-(
-    const bool binary,
-    const primitiveMesh& mesh,
-    const topoSet& set,
-    const fileName& fileName
-)
-{
-    std::ofstream pStream(fileName.c_str());
-
-    pStream
-        << "# vtk DataFile Version 2.0" << std::endl
-        << set.name() << std::endl;
-    if (binary)
-    {
-        pStream << "BINARY" << std::endl;
-    }
-    else
-    {
-        pStream << "ASCII" << std::endl;
-    }
-    pStream << "DATASET POLYDATA" << std::endl;
-
-
-    //------------------------------------------------------------------
-    //
-    // Write topology
-    //
-    //------------------------------------------------------------------
-
-
-    labelList pointLabels(set.toc());
-
-    pointField setPoints(mesh.points(), pointLabels);
-
-    // Write points
-
-    pStream << "POINTS " << pointLabels.size() << " float" << std::endl;
-
-    DynamicList<floatScalar> ptField(3*pointLabels.size());
-
-    vtkWriteOps::insert(setPoints, ptField);
-
-    vtkWriteOps::write(pStream, binary, ptField);
-
-    // Write vertices
-
-    pStream << "VERTICES " << pointLabels.size() << ' '
-        << 2*pointLabels.size() << std::endl;
-
-    DynamicList<label> vertexPoints(2*pointLabels.size());
-    forAll(pointLabels, pointi)
-    {
-        vertexPoints.append(1);
-        vertexPoints.append(pointi);
-    }
-    vtkWriteOps::write(pStream, binary, vertexPoints);
-
-    //-----------------------------------------------------------------
-    //
-    // Write data
-    //
-    //-----------------------------------------------------------------
-
-    // Write pointID
-
-    pStream
-        << "POINT_DATA " << pointLabels.size() << std::endl
-        << "FIELD attributes 1" << std::endl;
-
-    // Cell ids first
-    pStream << "pointID 1 " << pointLabels.size() << " int" << std::endl;
-
-    vtkWriteOps::write(pStream, binary, pointLabels);
-}
-
-
-// ************************************************************************* //
diff --git a/applications/utilities/mesh/manipulation/setSet/writePointSet.H b/applications/utilities/mesh/manipulation/setSet/writePointSet.H
deleted file mode 100644
index 6309461971..0000000000
--- a/applications/utilities/mesh/manipulation/setSet/writePointSet.H
+++ /dev/null
@@ -1,67 +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/>.
-
-InNamespace
-    Foam
-
-Description
-    Write pointSet to vtk polydata file. Only one data which is original
-    pointID.
-
-SourceFiles
-    writePointSet.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef writePointSet_H
-#define writePointSet_H
-
-#include "primitiveMesh.H"
-#include "pointSet.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-namespace vtkWriteOps
-{
-
-//- Write pointSet to vtk polydata file.
-//  Only one data which is original pointID.
-void writePointSet
-(
-    const bool binary,
-    const primitiveMesh& mesh,
-    const topoSet& set,
-    const fileName& fileName
-);
-
-} // End namespace vtkWriteOps
-} // End namespace Foam
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
diff --git a/src/fileFormats/vtk/vtkWritePolyData.H b/src/fileFormats/vtk/vtkWritePolyData.H
new file mode 100644
index 0000000000..685f163428
--- /dev/null
+++ b/src/fileFormats/vtk/vtkWritePolyData.H
@@ -0,0 +1,222 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 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::vtkWritePolyData
+
+Description
+    General write functions for vtk polygonal data files
+
+SourceFiles
+    vtkWritePolyData.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef vtkWritePolyData_H
+#define vtkWritePolyData_H
+
+#include "vtkWriteOps.H"
+#include "fileName.H"
+#include "Field.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                     Namespace vtkWritePolyData Declaration
+\*---------------------------------------------------------------------------*/
+
+namespace vtkWritePolyData
+{
+
+//- Helper for templated write
+template<class Types, class Type>
+inline void setFieldTypeValue
+(
+    UPtrList<const Field<Types>>& fieldTypeValues,
+    const label fieldi,
+    const Field<Type>& fieldTypeValue
+)
+{}
+
+
+//- Helper for templated write
+template<class Type>
+inline void setFieldTypeValue
+(
+    UPtrList<const Field<Type>>& fieldTypeValues,
+    const label fieldi,
+    const Field<Type>& fieldTypeValue
+)
+{
+    fieldTypeValues.set(fieldi, &fieldTypeValue);
+}
+
+
+//- Helper for templated write
+inline void unpackFieldTypeValues
+(
+    wordList& fieldNames,
+    boolList& fieldIsPointValues,
+    UPtrList<const Field<label>>& fieldLabelValues
+    #define FieldTypeValuesNonConstArg(Type, nullArg) \
+    , UPtrList<const Field<Type>>& field##Type##Values
+    FOR_ALL_FIELD_TYPES(FieldTypeValuesNonConstArg)
+    #undef FieldTypeValuesNonConstArg
+)
+{}
+
+
+//- Helper for templated write
+template<class Type, class ... Args>
+inline void unpackFieldTypeValues
+(
+    wordList& fieldNames,
+    boolList& fieldIsPointValues,
+    UPtrList<const Field<label>>& fieldLabelValues
+    #define FieldTypeValuesNonConstArg(Type, nullArg) \
+    , UPtrList<const Field<Type>>& field##Type##Values
+    FOR_ALL_FIELD_TYPES(FieldTypeValuesNonConstArg),
+    #undef FieldTypeValuesNonConstArg
+    const word& fieldName,
+    const bool fieldIsPointValue,
+    const Field<Type>& fieldTypeValue,
+    Args& ... args
+)
+{
+    const label fieldi = fieldNames.size() - 1 - sizeof...(Args)/3;
+
+    fieldNames[fieldi] = fieldName;
+    fieldIsPointValues[fieldi] = fieldIsPointValue;
+    setFieldTypeValue(fieldLabelValues, fieldi, fieldTypeValue);
+    #define SetFieldTypeValue(Type, nullArg) \
+        setFieldTypeValue(field##Type##Values, fieldi, fieldTypeValue);
+    FOR_ALL_FIELD_TYPES(SetFieldTypeValue);
+    #undef SetFieldTypeValue
+
+    unpackFieldTypeValues
+    (
+        fieldNames,
+        fieldIsPointValues,
+        fieldLabelValues
+        #define FieldTypeValuesParameter(Type, nullArg) , field##Type##Values
+        FOR_ALL_FIELD_TYPES(FieldTypeValuesParameter),
+        #undef FieldTypeValuesParameter
+        args ...
+    );
+}
+
+
+//- Write the field values out for a type
+template<class Type, class DataType>
+void writeFieldTypeValues
+(
+    std::ostream& os,
+    const bool binary,
+    const wordList& fieldNames,
+    const boolList& fieldIsPointValues,
+    const UPtrList<const Field<Type>>& fieldTypeValues,
+    const bool writePointValues
+);
+
+
+//- Write VTK polygonal data to a file. Takes a PtrList of fields of labels and
+//  of every primitive type. Each PtrList should be the length of the total
+//  number of fields and only one PtrList should be non-null for each field
+//  index.
+template<class PointField, class VertexList, class LineList, class FaceList>
+void write
+(
+    const fileName& file,
+    const word& title,
+    const bool binary,
+    const PointField& points,
+    const VertexList& vertices,
+    const LineList& lines,
+    const FaceList& faces,
+    const wordList& fieldNames,
+    const boolList& fieldIsPointValues,
+    const UPtrList<const Field<label>>& fieldLabelValues
+    #define FieldTypeValuesConstArg(Type, nullArg) \
+        , const UPtrList<const Field<Type>>& field##Type##Values
+    FOR_ALL_FIELD_TYPES(FieldTypeValuesConstArg)
+    #undef FieldTypeValuesConstArg
+);
+
+
+//- Write VTK polygonal data to a file. Takes any number of name,
+//  isPointValues, values arguments at the end. E.g.;
+//
+//  write
+//  (
+//      // Output options
+//      "myPolyData.vtk", "myPolyData", false,
+//
+//      // Geometry
+//      pp.localPoints(), labelList(), labelListList(), pp.localFaces(),
+//
+//      // Fields
+//      "faceIDs", true, pp.addressing(),
+//      "facePressures", false, Field<scalar>(pp.size(), ...),
+//      "pointVelocities", true, Field<vector>(pp.nPoints(), ...)
+//  );
+//
+template
+<
+    class PointField,
+    class VertexList,
+    class LineList,
+    class FaceList,
+    class ... Args
+>
+inline void write
+(
+    const fileName& file,
+    const word& title,
+    const bool binary,
+    const PointField& points,
+    const VertexList& vertices,
+    const LineList& lines,
+    const FaceList& faces,
+    const Args& ... args
+);
+
+};
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+    #include "vtkWritePolyDataTemplates.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/fileFormats/vtk/vtkWritePolyDataTemplates.C b/src/fileFormats/vtk/vtkWritePolyDataTemplates.C
new file mode 100644
index 0000000000..94596eebeb
--- /dev/null
+++ b/src/fileFormats/vtk/vtkWritePolyDataTemplates.C
@@ -0,0 +1,297 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Copyright (C) 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 "vtkWritePolyData.H"
+#include "OFstream.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+template<class Type, class DataType>
+void Foam::vtkWritePolyData::writeFieldTypeValues
+(
+    std::ostream& os,
+    const bool binary,
+    const wordList& fieldNames,
+    const boolList& fieldIsPointValues,
+    const UPtrList<const Field<Type>>& fieldTypeValues,
+    const bool writePointValues
+)
+{
+    forAll(fieldNames, fieldi)
+    {
+        if
+        (
+            fieldIsPointValues[fieldi] == writePointValues
+         && fieldTypeValues.set(fieldi)
+        )
+        {
+            const label nCmpt = pTraits<Type>::nComponents;
+
+            os  << fieldNames[fieldi] << ' ' << pTraits<Type>::nComponents
+                << ' ' << fieldTypeValues[fieldi].size() << ' '
+                << (std::is_integral<DataType>::value ? "int" : "float") << nl;
+
+            List<DataType> data(nCmpt*fieldTypeValues[fieldi].size());
+            label i = 0;
+            forAll(fieldTypeValues[fieldi], fieldValuei)
+            {
+                for (direction cmpt = 0; cmpt < nCmpt; ++ cmpt)
+                {
+                    data[i ++] =
+                        component
+                        (
+                            fieldTypeValues[fieldi][fieldValuei],
+                            cmpt
+                        );
+                }
+            }
+
+            vtkWriteOps::write(os, binary, data);
+        }
+    }
+}
+
+
+template<class PointField, class VertexList, class LineList, class FaceList>
+void Foam::vtkWritePolyData::write
+(
+    const fileName& file,
+    const word& title,
+    const bool binary,
+    const PointField& points,
+    const VertexList& vertices,
+    const LineList& lines,
+    const FaceList& faces,
+    const wordList& fieldNames,
+    const boolList& fieldIsPointValues,
+    const UPtrList<const Field<label>>& fieldLabelValues
+    #define FieldTypeValuesConstArg(Type, nullArg) \
+        , const UPtrList<const Field<Type>>& field##Type##Values
+    FOR_ALL_FIELD_TYPES(FieldTypeValuesConstArg)
+    #undef FieldTypeValuesConstArg
+)
+{
+    // Open the file
+    std::ofstream os(file, std::ios::binary);
+
+    // Write the header
+    vtkWriteOps::writeHeader(os, binary, title);
+    os << "DATASET POLYDATA" << nl;
+
+    // Write the points
+    {
+        os  << "POINTS " << points.size() << " float" << nl;
+        List<floatScalar> coordinates(points.size()*3);
+        forAll(points, pointi)
+        {
+            const point& p = points[pointi];
+            forAll(p, i)
+            {
+                coordinates[3*pointi + i] = float(p[i]);
+            }
+        }
+        vtkWriteOps::write(os, binary, coordinates);
+    }
+
+    // Write the vertices
+    if (vertices.size())
+    {
+        os  << "VERTICES " << vertices.size() << ' '
+            << 2*vertices.size() << nl;
+        labelList data(2*vertices.size());
+        forAll(vertices, vertexi)
+        {
+            data[2*vertexi] = 1;
+            data[2*vertexi + 1] = vertices[vertexi];
+        }
+        vtkWriteOps::write(os, binary, data);
+    }
+
+    // Write the lines
+    if (lines.size())
+    {
+        label nLineNodes = 0;
+        forAll(lines, facei)
+        {
+            nLineNodes += lines[facei].size();
+        }
+        os  << "LINES " << lines.size() << ' '
+            << lines.size() + nLineNodes << nl;
+        labelList data(lines.size() + nLineNodes);
+        label i = 0;
+        forAll(lines, linei)
+        {
+            data[i ++] = lines[linei].size();
+            forAll(lines[linei], linePointi)
+            {
+                data[i ++] = lines[linei][linePointi];
+            }
+        }
+        vtkWriteOps::write(os, binary, data);
+    }
+
+    // Write the faces
+    if (faces.size())
+    {
+        label nFaceNodes = 0;
+        forAll(faces, facei)
+        {
+            nFaceNodes += faces[facei].size();
+        }
+        os  << "POLYGONS " << faces.size() << ' '
+            << faces.size() + nFaceNodes << nl;
+        labelList data(faces.size() + nFaceNodes);
+        label i = 0;
+        forAll(faces, facei)
+        {
+            data[i ++] = faces[facei].size();
+            forAll(faces[facei], facePointi)
+            {
+                data[i ++] = faces[facei][facePointi];
+            }
+        }
+        vtkWriteOps::write(os, binary, data);
+    }
+
+    // Write the fields
+    const label nPointFields = count(fieldIsPointValues, true);
+    const label nFaceFields = count(fieldIsPointValues, false);
+    if (nPointFields > 0)
+    {
+        os  << "POINT_DATA " << points.size() << nl
+            << "FIELD attributes " << nPointFields << nl;
+        writeFieldTypeValues<label, label>
+        (
+            os,
+            binary,
+            fieldNames,
+            fieldIsPointValues,
+            fieldLabelValues,
+            true
+        );
+        #define WriteFieldTypeValues(Type, nullArg)                            \
+            writeFieldTypeValues<Type, floatScalar>                            \
+            (                                                                  \
+                os,                                                            \
+                binary,                                                        \
+                fieldNames,                                                    \
+                fieldIsPointValues,                                            \
+                field##Type##Values,                                           \
+                true                                                           \
+            );
+        FOR_ALL_FIELD_TYPES(WriteFieldTypeValues)
+        #undef WriteFieldTypeValues
+    }
+    if (nFaceFields > 0)
+    {
+        os  << "CELL_DATA "
+            << vertices.size() + lines.size() + faces.size() << nl
+            << "FIELD attributes " << nFaceFields << nl;
+        writeFieldTypeValues<label, label>
+        (
+            os,
+            binary,
+            fieldNames,
+            fieldIsPointValues,
+            fieldLabelValues,
+            false
+        );
+        #define WriteFieldTypeValues(Type, nullArg)                            \
+            writeFieldTypeValues<Type, floatScalar>                            \
+            (                                                                  \
+                os,                                                            \
+                binary,                                                        \
+                fieldNames,                                                    \
+                fieldIsPointValues,                                            \
+                field##Type##Values,                                           \
+                false                                                          \
+            );
+        FOR_ALL_FIELD_TYPES(WriteFieldTypeValues)
+        #undef WriteFieldTypeValues
+    }
+}
+
+
+template
+<
+    class PointField,
+    class VertexList,
+    class LineList,
+    class FaceList,
+    class ... Args
+>
+inline void Foam::vtkWritePolyData::write
+(
+    const fileName& file,
+    const word& title,
+    const bool binary,
+    const PointField& points,
+    const VertexList& vertices,
+    const LineList& lines,
+    const FaceList& faces,
+    const Args& ... args
+)
+{
+    const label nFields = sizeof...(Args)/3;
+
+    wordList fieldNames(nFields);
+    boolList fieldIsPointValues(nFields);
+    UPtrList<const Field<label>> fieldLabelValues(nFields);
+    #define DeclareFieldTypeValues(Type, nullArg) \
+        UPtrList<const Field<Type>> field##Type##Values(nFields);
+    FOR_ALL_FIELD_TYPES(DeclareFieldTypeValues);
+    #undef DeclareFieldTypeValues
+
+    unpackFieldTypeValues
+    (
+        fieldNames,
+        fieldIsPointValues,
+        fieldLabelValues
+        #define FieldTypeValuesParameter(Type, nullArg) , field##Type##Values
+        FOR_ALL_FIELD_TYPES(FieldTypeValuesParameter),
+        #undef FieldTypeValuesParameter
+        args ...
+    );
+
+    write
+    (
+        file,
+        title,
+        binary,
+        points,
+        vertices,
+        lines,
+        faces,
+        fieldNames,
+        fieldIsPointValues,
+        fieldLabelValues
+        #define FieldTypeValuesParameter(Type, nullArg) , field##Type##Values
+        FOR_ALL_FIELD_TYPES(FieldTypeValuesParameter)
+        #undef FieldTypeValuesParameter
+    );
+}
+
+
+// ************************************************************************* //
diff --git a/src/sampling/sampledSurface/writers/vtk/vtkSurfaceWriter.C b/src/sampling/sampledSurface/writers/vtk/vtkSurfaceWriter.C
index 8f02d97f8b..a1d7d08c8e 100644
--- a/src/sampling/sampledSurface/writers/vtk/vtkSurfaceWriter.C
+++ b/src/sampling/sampledSurface/writers/vtk/vtkSurfaceWriter.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-2021 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -27,7 +27,7 @@ License
 #include "OFstream.H"
 #include "OSspecific.H"
 #include "makeSurfaceWriterMethods.H"
-#include "vtkWriteOps.H"
+#include "vtkWritePolyData.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -39,59 +39,6 @@ namespace Foam
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
-void Foam::vtkSurfaceWriter::writeGeometry
-(
-    std::ostream& os,
-    const pointField& points,
-    const faceList& faces
-) const
-{
-    const bool binary = (writeFormat_ == IOstream::BINARY);
-
-    // VTK header
-    vtkWriteOps::writeHeader(os, binary, "sampleSurface");
-    os << "DATASET POLYDATA" << nl;
-
-    // Write vertex coords
-    os  << "POINTS " << points.size() << " float" << nl;
-
-    List<floatScalar> po(points.size()*3);
-    label ind = 0;
-    forAll(points, pointi)
-    {
-        const point& pt = points[pointi];
-        forAll(pt, cmpt)
-        {
-            po[ind++] = float(pt[cmpt]);
-        }
-    }
-    vtkWriteOps::write(os, binary, po);
-
-    // Write faces
-    label nNodes = 0;
-    forAll(faces, facei)
-    {
-        nNodes += faces[facei].size();
-    }
-
-    os  << "POLYGONS " << faces.size() << ' '
-        << faces.size() + nNodes << nl;
-
-    labelList polygons(faces.size() + nNodes);
-    ind = 0;
-    forAll(faces, facei)
-    {
-        const face& f = faces[facei];
-        polygons[ind++] = f.size();
-        forAll(f, fp)
-        {
-            polygons[ind++] = f[fp];
-        }
-    }
-    vtkWriteOps::write(os, binary, polygons);
-}
-
-
 template<class Type>
 void Foam::vtkSurfaceWriter::Write
 (
@@ -104,53 +51,24 @@ void Foam::vtkSurfaceWriter::Write
     const bool isNodeValues
 ) const
 {
-    const bool binary = (writeFormat_ == IOstream::BINARY);
-
     if (!isDir(outputDir))
     {
         mkDir(outputDir);
     }
 
-    const word filePath = outputDir/fieldName + '_' + surfaceName + ".vtk";
-
-    ofstream os(filePath, std::ios::binary);
-
-    if (debug)
-    {
-        Info<< "Writing field " << fieldName << " to " << filePath << endl;
-    }
-
-    writeGeometry(os, points, faces);
-
-    // Write data
-    if (isNodeValues)
-    {
-        os  << "POINT_DATA ";
-    }
-    else
-    {
-        os  << "CELL_DATA ";
-    }
-
-    os  << values.size() << nl
-        << "FIELD attributes 1" << nl
-        << fieldName << " ";
-
-    const label nComp = pTraits<Type>::nComponents;
-
-    os  << nComp << " " << values.size() << " float" << nl;
-
-    List<floatScalar> vals(values.size()*nComp);
-    label ind = 0;
-    forAll(values, elemI)
-    {
-        for (direction cmpt=0; cmpt < nComp; ++cmpt)
-        {
-            vals[ind++] = component(values[elemI], cmpt);
-        }
-    }
-
-    vtkWriteOps::write(os, binary, vals);
+    vtkWritePolyData::write
+    (
+        outputDir/fieldName + '_' + surfaceName + ".vtk",
+        "sampleSurface",
+        writeFormat_ == IOstream::BINARY,
+        points,
+        labelList(),
+        edgeList(),
+        faces,
+        fieldName,
+        isNodeValues,
+        values
+    );
 }
 
 
@@ -186,15 +104,16 @@ void Foam::vtkSurfaceWriter::write
         mkDir(outputDir);
     }
 
-    word filePath =  outputDir/surfaceName + ".vtk";
-    ofstream os(filePath, std::ios::binary);
-
-    if (debug)
-    {
-        Info<< "Writing geometry to " << filePath << endl;
-    }
-
-    writeGeometry(os, points, faces);
+    vtkWritePolyData::write
+    (
+        outputDir/surfaceName + ".vtk",
+        "sampleSurface",
+        writeFormat_ == IOstream::BINARY,
+        points,
+        labelList(),
+        edgeList(),
+        faces
+    );
 }
 
 
diff --git a/src/sampling/sampledSurface/writers/vtk/vtkSurfaceWriter.H b/src/sampling/sampledSurface/writers/vtk/vtkSurfaceWriter.H
index a44fed6d49..514ca31075 100644
--- a/src/sampling/sampledSurface/writers/vtk/vtkSurfaceWriter.H
+++ b/src/sampling/sampledSurface/writers/vtk/vtkSurfaceWriter.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-2021 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -53,13 +53,6 @@ class vtkSurfaceWriter
 {
     // Private Member Functions
 
-        void writeGeometry
-        (
-            std::ostream&,
-            const pointField&,
-            const faceList&
-        ) const;
-
         //- Templated write operation
         template<class Type>
         void Write
@@ -101,7 +94,6 @@ public:
             const faceList& faces
         ) const;
 
-
         //- Write scalarField for a single surface to file.
         //  One value per face or vertex (isNodeValues = true)
         virtual void write