diff --git a/dune/vtk/CMakeLists.txt b/dune/vtk/CMakeLists.txt
index 4779d6df3f8bbc8c4e17129bb24eeec1dd8b8f95..e7aec1ea0c585cd84ee4ed1b9169d4398dd4a780 100644
--- a/dune/vtk/CMakeLists.txt
+++ b/dune/vtk/CMakeLists.txt
@@ -1,5 +1,5 @@
 dune_add_library("vtktypes" OBJECT
-  vtktypes.cc)
+  types.cc)
 
 #install headers
 install(FILES
@@ -8,20 +8,22 @@ install(FILES
   defaultvtkfunction.hh
   filereader.hh
   filewriter.hh
+  forward.hh
+  gridcreatorinterface.hh
   legacyvtkfunction.hh
   pvdwriter.hh
   pvdwriter.impl.hh
-  vtkfunction.hh
-  vtklocalfunction.hh
-  vtklocalfunctioninterface.hh
-  vtkreader.hh
-  vtkreader.impl.hh
-  vtktimeserieswriter.hh
-  vtktimeserieswriter.impl.hh
-  vtktypes.hh
-  vtkwriter.hh
-  vtkwriterinterface.hh
-  vtkwriterinterface.impl.hh
+  function.hh
+  localfunction.hh
+  localfunctioninterface.hh
+  reader.hh
+  reader.impl.hh
+  timeserieswriter.hh
+  timeserieswriter.impl.hh
+  types.hh
+  writer.hh
+  writerinterface.hh
+  writerinterface.impl.hh
   DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dune/vtkwriter)
 
 add_subdirectory(datacollectors)
diff --git a/dune/vtk/datacollectorinterface.hh b/dune/vtk/datacollectorinterface.hh
index d37ca636ba8cb36c6e24b4c29d562df54b5ded54..6afcba8377767ab6aee12a45b84299ada4d38572 100644
--- a/dune/vtk/datacollectorinterface.hh
+++ b/dune/vtk/datacollectorinterface.hh
@@ -5,130 +5,135 @@
 
 namespace Dune
 {
-  /// Base class for data collectors in a CRTP style.
-  /**
-   * \tparam GridViewType  Model of Dune::GridView
-   * \tparam Derived       Implementation of a concrete DataCollector.
-   * \tparam Partition     Dune::PartitionType [Partitions::InteriorBorder]
-   **/
-  template <class GridViewType, class Derived, class Partition>
-  class DataCollectorInterface
-  {
-  public:
-    /// The partitionset to collect data from
-    static constexpr auto partition = Partition{};
-
-    using GridView = GridViewType;
-
-    /// The dimension of the grid
-    enum { dim = GridView::dimension };
-
-    /// The dimension of the world
-    enum { dow = GridView::dimensionworld };
-
-  public:
-    /// Store a copy of the GridView
-    DataCollectorInterface (GridView const& gridView)
-      : gridView_(gridView)
-    {}
-
-    /// Update the DataCollector on the current GridView
-    void update ()
-    {
-      asDerived().updateImpl();
-    }
-
-    /// Return the number of ghost elements
-    int ghostLevel () const
-    {
-      return asDerived().ghostLevelImpl();
-    }
 
-    /// Return the number of cells in (this partition of the) grid
-    std::uint64_t numCells () const
-    {
-      return asDerived().numCellsImpl();
-    }
-
-    /// Return the number of points in (this partition of the) grid
-    std::uint64_t numPoints () const
-    {
-      return asDerived().numPointsImpl();
-    }
-
-    /// \brief Return a flat vector of point coordinates
-    /**
-    * All coordinates are extended to 3 components and concatenated.
-    * [p0_x, p0_y, p0_z, p1_x, p1_y, p1_z, ...]
-    * If the GridView::dimensionworld < 3, the remaining components are
-    * set to 0
-    **/
-    template <class T>
-    std::vector<T> points () const
-    {
-      return asDerived().template pointsImpl<T>();
-    }
+  namespace Vtk
+  {
 
-    /// \brief Return a flat vector of function values evaluated at the points.
+    /// Base class for data collectors in a CRTP style.
     /**
-    * In case of a vector valued function, flat the vector entries:
-    * [fct(p0)_0, fct(p0)_1, fct(p0)_2, fct(p1)_0, ...]
-    * where the vector dimension must be 3 (possible extended by 0s)
-    * In case of tensor valued function, flat the tensor row-wise:
-    * [fct(p0)_00, fct(p0)_01, fct(p0)_02, fct(p0)_10, fct(p0)_11, fct(p0)_12, fct(p0)_20...]
-    * where the tensor dimension must be 3x3 (possible extended by 0s)
+    * \tparam GridViewType  Model of Dune::GridView
+    * \tparam Derived       Implementation of a concrete DataCollector.
+    * \tparam Partition     Dune::PartitionType [Partitions::InteriorBorder]
     **/
-    template <class T, class VtkFunction>
-    std::vector<T> pointData (VtkFunction const& fct) const
+    template <class GridViewType, class Derived, class Partition>
+    class DataCollectorInterface
     {
-      return asDerived().template pointDataImpl<T>(fct);
-    }
-
-    /// \brief Return a flat vector of function values evaluated at the cells in the order of traversal.
-    /**
-    * \see pointData.
-    * Note: Cells might be described explicitly by connectivity, offsets, and types, e.g. in
-    * an UnstructuredGrid, or might be described implicitly by the grid type, e.g. in
-    * StructuredGrid.
-    */
-    template <class T, class VtkFunction>
-    std::vector<T> cellData (VtkFunction const& fct) const
-    {
-      return asDerived().template cellDataImpl<T>(fct);
-    }
-
-  protected: // cast to derived type
-
-    Derived& asDerived ()
-    {
-      return static_cast<Derived&>(*this);
-    }
-
-    const Derived& asDerived () const
-    {
-      return static_cast<const Derived&>(*this);
-    }
-
-  public: // default implementations
-
-    void updateImpl ()
-    {
-      /* do nothing */
-    }
-
-    int ghostLevelImpl () const
-    {
-      return gridView_.overlapSize(0);
-    }
-
-    // Evaluate `fct` in center of cell.
-    template <class T, class VtkFunction>
-    std::vector<T> cellDataImpl (VtkFunction const& fct) const;
-
-  protected:
-    GridView gridView_;
-  };
-
+    public:
+      /// The partitionset to collect data from
+      static constexpr auto partition = Partition{};
+
+      using GridView = GridViewType;
+
+      /// The dimension of the grid
+      enum { dim = GridView::dimension };
+
+      /// The dimension of the world
+      enum { dow = GridView::dimensionworld };
+
+    public:
+      /// Store a copy of the GridView
+      DataCollectorInterface (GridView const& gridView)
+        : gridView_(gridView)
+      {}
+
+      /// Update the DataCollector on the current GridView
+      void update ()
+      {
+        asDerived().updateImpl();
+      }
+
+      /// Return the number of ghost elements
+      int ghostLevel () const
+      {
+        return asDerived().ghostLevelImpl();
+      }
+
+      /// Return the number of cells in (this partition of the) grid
+      std::uint64_t numCells () const
+      {
+        return asDerived().numCellsImpl();
+      }
+
+      /// Return the number of points in (this partition of the) grid
+      std::uint64_t numPoints () const
+      {
+        return asDerived().numPointsImpl();
+      }
+
+      /// \brief Return a flat vector of point coordinates
+      /**
+      * All coordinates are extended to 3 components and concatenated.
+      * [p0_x, p0_y, p0_z, p1_x, p1_y, p1_z, ...]
+      * If the GridView::dimensionworld < 3, the remaining components are
+      * set to 0
+      **/
+      template <class T>
+      std::vector<T> points () const
+      {
+        return asDerived().template pointsImpl<T>();
+      }
+
+      /// \brief Return a flat vector of function values evaluated at the points.
+      /**
+      * In case of a vector valued function, flat the vector entries:
+      * [fct(p0)_0, fct(p0)_1, fct(p0)_2, fct(p1)_0, ...]
+      * where the vector dimension must be 3 (possible extended by 0s)
+      * In case of tensor valued function, flat the tensor row-wise:
+      * [fct(p0)_00, fct(p0)_01, fct(p0)_02, fct(p0)_10, fct(p0)_11, fct(p0)_12, fct(p0)_20...]
+      * where the tensor dimension must be 3x3 (possible extended by 0s)
+      **/
+      template <class T, class VtkFunction>
+      std::vector<T> pointData (VtkFunction const& fct) const
+      {
+        return asDerived().template pointDataImpl<T>(fct);
+      }
+
+      /// \brief Return a flat vector of function values evaluated at the cells in the order of traversal.
+      /**
+      * \see pointData.
+      * Note: Cells might be described explicitly by connectivity, offsets, and types, e.g. in
+      * an UnstructuredGrid, or might be described implicitly by the grid type, e.g. in
+      * StructuredGrid.
+      */
+      template <class T, class VtkFunction>
+      std::vector<T> cellData (VtkFunction const& fct) const
+      {
+        return asDerived().template cellDataImpl<T>(fct);
+      }
+
+    protected: // cast to derived type
+
+      Derived& asDerived ()
+      {
+        return static_cast<Derived&>(*this);
+      }
+
+      const Derived& asDerived () const
+      {
+        return static_cast<const Derived&>(*this);
+      }
+
+    public: // default implementations
+
+      void updateImpl ()
+      {
+        /* do nothing */
+      }
+
+      int ghostLevelImpl () const
+      {
+        return gridView_.overlapSize(0);
+      }
+
+      // Evaluate `fct` in center of cell.
+      template <class T, class VtkFunction>
+      std::vector<T> cellDataImpl (VtkFunction const& fct) const;
+
+    protected:
+      GridView gridView_;
+    };
+
+  } // end namespace Vtk
 } // end namespace Dune
 
 #include "datacollectorinterface.impl.hh"
diff --git a/dune/vtk/datacollectorinterface.impl.hh b/dune/vtk/datacollectorinterface.impl.hh
index 957a6c0d4714e752edafa7885ebc6986de5ee05a..10f558fdbb3548640876d7abef4c12f241dc990c 100644
--- a/dune/vtk/datacollectorinterface.impl.hh
+++ b/dune/vtk/datacollectorinterface.impl.hh
@@ -2,25 +2,30 @@
 
 #include <dune/geometry/referenceelements.hh>
 
-namespace Dune {
-
-template <class GV, class D, class P>
-  template <class T, class VtkFunction>
-std::vector<T> DataCollectorInterface<GV,D,P>
-  ::cellDataImpl (VtkFunction const& fct) const
+namespace Dune
 {
-  std::vector<T> data;
-  data.reserve(this->numCells() * fct.ncomps());
 
-  auto localFct = localFunction(fct);
-  for (auto const& e : elements(gridView_, partition)) {
-    localFct.bind(e);
-    auto refElem = referenceElement<T,dim>(e.type());
-    for (int comp = 0; comp < fct.ncomps(); ++comp)
-      data.emplace_back(localFct.evaluate(comp, refElem.position(0,0)));
-    localFct.unbind();
-  }
-  return data;
-}
+  namespace Vtk
+  {
+
+    template <class GV, class D, class P>
+      template <class T, class VtkFunction>
+    std::vector<T> DataCollectorInterface<GV,D,P>
+      ::cellDataImpl (VtkFunction const& fct) const
+    {
+      std::vector<T> data;
+      data.reserve(this->numCells() * fct.ncomps());
+
+      auto localFct = localFunction(fct);
+      for (auto const& e : elements(gridView_, partition)) {
+        localFct.bind(e);
+        auto refElem = referenceElement<T,dim>(e.type());
+        for (int comp = 0; comp < fct.ncomps(); ++comp)
+          data.emplace_back(localFct.evaluate(comp, refElem.position(0,0)));
+        localFct.unbind();
+      }
+      return data;
+    }
 
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/datacollectors/continuousdatacollector.hh b/dune/vtk/datacollectors/continuousdatacollector.hh
index e6460d1accd6b41d2a44198954c0627a669f6d41..ce2e0f68356d73b79b9a12363b305e514805b6cf 100644
--- a/dune/vtk/datacollectors/continuousdatacollector.hh
+++ b/dune/vtk/datacollectors/continuousdatacollector.hh
@@ -7,140 +7,144 @@
 #include <dune/grid/utility/globalindexset.hh>
 
 #include <dune/vtk/forward.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
 
 #include "unstructureddatacollector.hh"
 
 namespace Dune
 {
 
-/// Implementation of \ref DataCollector for linear cells, with continuous data.
-template <class GridView, class Partition>
-class ContinuousDataCollector
-    : public UnstructuredDataCollectorInterface<GridView, ContinuousDataCollector<GridView,Partition>, Partition>
-{
-  using Self = ContinuousDataCollector;
-  using Super = UnstructuredDataCollectorInterface<GridView, Self, Partition>;
+  namespace Vtk
+  {
 
-public:
-  using Super::dim;
-  using Super::partition;
+    /// Implementation of \ref DataCollector for linear cells, with continuous data.
+    template <class GridView, class Partition>
+    class ContinuousDataCollector
+        : public UnstructuredDataCollectorInterface<GridView, ContinuousDataCollector<GridView,Partition>, Partition>
+    {
+      using Self = ContinuousDataCollector;
+      using Super = UnstructuredDataCollectorInterface<GridView, Self, Partition>;
+
+    public:
+      using Super::dim;
+      using Super::partition;
+
+    public:
+      ContinuousDataCollector (GridView const& gridView)
+        : Super(gridView)
+      {}
+
+      /// Collect the vertex indices
+      void updateImpl ()
+      {
+        numPoints_ = 0;
+        indexMap_.resize(gridView_.size(dim));
+        auto const& indexSet = gridView_.indexSet();
+        for (auto const& vertex : vertices(gridView_, partition))
+          indexMap_[indexSet.index(vertex)] = std::int64_t(numPoints_++);
+
+        if (gridView_.comm().size() > 1) {
+          auto&& e = elements(gridView_, partition);
+          numCells_ = std::distance(std::begin(e), std::end(e));
+        } else {
+          numCells_ = gridView_.size(0);
+        }
+      }
 
-public:
-  ContinuousDataCollector (GridView const& gridView)
-    : Super(gridView)
-  {}
+      /// Return number of grid vertices
+      std::uint64_t numPointsImpl () const
+      {
+        return numPoints_;
+      }
 
-  /// Collect the vertex indices
-  void updateImpl ()
-  {
-    numPoints_ = 0;
-    indexMap_.resize(gridView_.size(dim));
-    auto const& indexSet = gridView_.indexSet();
-    for (auto const& vertex : vertices(gridView_, partition))
-      indexMap_[indexSet.index(vertex)] = std::int64_t(numPoints_++);
-
-    if (gridView_.comm().size() > 1) {
-      auto&& e = elements(gridView_, partition);
-      numCells_ = std::distance(std::begin(e), std::end(e));
-    } else {
-      numCells_ = gridView_.size(0);
-    }
-  }
-
-  /// Return number of grid vertices
-  std::uint64_t numPointsImpl () const
-  {
-    return numPoints_;
-  }
+      /// Return the coordinates of all grid vertices in the order given by the indexSet
+      template <class T>
+      std::vector<T> pointsImpl () const
+      {
+        std::vector<T> data;
+        data.reserve(numPoints_ * 3);
+        for (auto const& vertex : vertices(gridView_, partition)) {
+          auto v = vertex.geometry().center();
+          for (std::size_t j = 0; j < v.size(); ++j)
+            data.emplace_back(v[j]);
+          for (std::size_t j = v.size(); j < 3u; ++j)
+            data.emplace_back(0);
+        }
+        return data;
+      }
 
-  /// Return the coordinates of all grid vertices in the order given by the indexSet
-  template <class T>
-  std::vector<T> pointsImpl () const
-  {
-    std::vector<T> data;
-    data.reserve(numPoints_ * 3);
-    for (auto const& vertex : vertices(gridView_, partition)) {
-      auto v = vertex.geometry().center();
-      for (std::size_t j = 0; j < v.size(); ++j)
-        data.emplace_back(v[j]);
-      for (std::size_t j = v.size(); j < 3u; ++j)
-        data.emplace_back(0);
-    }
-    return data;
-  }
-
-  /// Return a vector of global unique ids of the points
-  std::vector<std::uint64_t> pointIdsImpl () const
-  {
-    std::vector<std::uint64_t> data;
-    data.reserve(numPoints_);
-    GlobalIndexSet<GridView> globalIndexSet(gridView_, dim);
-    for (auto const& vertex : vertices(gridView_, partition)) {
-      data.emplace_back(globalIndexSet.index(vertex));
-    }
-    return data;
-  }
-
-  /// Return number of grid cells
-  std::uint64_t numCellsImpl () const
-  {
-    return numCells_;
-  }
+      /// Return a vector of global unique ids of the points
+      std::vector<std::uint64_t> pointIdsImpl () const
+      {
+        std::vector<std::uint64_t> data;
+        data.reserve(numPoints_);
+        GlobalIndexSet<GridView> globalIndexSet(gridView_, dim);
+        for (auto const& vertex : vertices(gridView_, partition)) {
+          data.emplace_back(globalIndexSet.index(vertex));
+        }
+        return data;
+      }
 
-  /// Return the types, offsets and connectivity of the cells, using the same connectivity as
-  /// given by the grid.
-  Cells cellsImpl () const
-  {
-    auto const& indexSet = gridView_.indexSet();
-    auto types = indexSet.types(0);
-    int maxVertices = std::accumulate(types.begin(), types.end(), 1, [](int m, GeometryType t) {
-      auto refElem = referenceElement<double,dim>(t);
-      return std::max(m, refElem.size(dim));
-    });
-
-    Cells cells;
-    cells.connectivity.reserve(numCells_ * maxVertices);
-    cells.offsets.reserve(numCells_);
-    cells.types.reserve(numCells_);
-
-    std::int64_t old_o = 0;
-    for (auto const& c : elements(gridView_, partition)) {
-      Vtk::CellType cellType(c.type());
-      for (unsigned int j = 0; j < c.subEntities(dim); ++j)
-        cells.connectivity.emplace_back(indexMap_[indexSet.subIndex(c,cellType.permutation(j),dim)]);
-      cells.offsets.push_back(old_o += c.subEntities(dim));
-      cells.types.push_back(cellType.type());
-    }
-    return cells;
-  }
-
-  /// Evaluate the `fct` at the corners of the elements
-  template <class T, class GlobalFunction>
-  std::vector<T> pointDataImpl (GlobalFunction const& fct) const
-  {
-    std::vector<T> data(numPoints_ * fct.ncomps());
-    auto const& indexSet = gridView_.indexSet();
-    auto localFct = localFunction(fct);
-    for (auto const& e : elements(gridView_, partition)) {
-      localFct.bind(e);
-      Vtk::CellType cellType{e.type()};
-      auto refElem = referenceElement(e.geometry());
-      for (unsigned int j = 0; j < e.subEntities(dim); ++j) {
-        std::size_t idx = fct.ncomps() * indexMap_[indexSet.subIndex(e,cellType.permutation(j),dim)];
-        for (int comp = 0; comp < fct.ncomps(); ++comp)
-          data[idx + comp] = T(localFct.evaluate(comp, refElem.position(cellType.permutation(j),dim)));
+      /// Return number of grid cells
+      std::uint64_t numCellsImpl () const
+      {
+        return numCells_;
       }
-      localFct.unbind();
-    }
-    return data;
-  }
-
-protected:
-  using Super::gridView_;
-  std::uint64_t numPoints_ = 0;
-  std::uint64_t numCells_ = 0;
-  std::vector<std::int64_t> indexMap_;
-};
 
+      /// Return the types, offsets and connectivity of the cells, using the same connectivity as
+      /// given by the grid.
+      Cells cellsImpl () const
+      {
+        auto const& indexSet = gridView_.indexSet();
+        auto types = indexSet.types(0);
+        int maxVertices = std::accumulate(types.begin(), types.end(), 1, [](int m, GeometryType t) {
+          auto refElem = referenceElement<double,dim>(t);
+          return std::max(m, refElem.size(dim));
+        });
+
+        Cells cells;
+        cells.connectivity.reserve(numCells_ * maxVertices);
+        cells.offsets.reserve(numCells_);
+        cells.types.reserve(numCells_);
+
+        std::int64_t old_o = 0;
+        for (auto const& c : elements(gridView_, partition)) {
+          Vtk::CellType cellType(c.type());
+          for (unsigned int j = 0; j < c.subEntities(dim); ++j)
+            cells.connectivity.emplace_back(indexMap_[indexSet.subIndex(c,cellType.permutation(j),dim)]);
+          cells.offsets.push_back(old_o += c.subEntities(dim));
+          cells.types.push_back(cellType.type());
+        }
+        return cells;
+      }
+
+      /// Evaluate the `fct` at the corners of the elements
+      template <class T, class GlobalFunction>
+      std::vector<T> pointDataImpl (GlobalFunction const& fct) const
+      {
+        std::vector<T> data(numPoints_ * fct.ncomps());
+        auto const& indexSet = gridView_.indexSet();
+        auto localFct = localFunction(fct);
+        for (auto const& e : elements(gridView_, partition)) {
+          localFct.bind(e);
+          Vtk::CellType cellType{e.type()};
+          auto refElem = referenceElement(e.geometry());
+          for (unsigned int j = 0; j < e.subEntities(dim); ++j) {
+            std::size_t idx = fct.ncomps() * indexMap_[indexSet.subIndex(e,cellType.permutation(j),dim)];
+            for (int comp = 0; comp < fct.ncomps(); ++comp)
+              data[idx + comp] = T(localFct.evaluate(comp, refElem.position(cellType.permutation(j),dim)));
+          }
+          localFct.unbind();
+        }
+        return data;
+      }
+
+    protected:
+      using Super::gridView_;
+      std::uint64_t numPoints_ = 0;
+      std::uint64_t numCells_ = 0;
+      std::vector<std::int64_t> indexMap_;
+    };
+
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/datacollectors/discontinuousdatacollector.hh b/dune/vtk/datacollectors/discontinuousdatacollector.hh
index eab3952372315889a5a226ad892a65c3a9c96c8a..adafba58b25f3033be46907eb4f50cba19ebfc03 100644
--- a/dune/vtk/datacollectors/discontinuousdatacollector.hh
+++ b/dune/vtk/datacollectors/discontinuousdatacollector.hh
@@ -3,126 +3,130 @@
 #include <vector>
 
 #include <dune/vtk/forward.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
 
 #include "unstructureddatacollector.hh"
 
 namespace Dune
 {
 
-/// Implementation of \ref DataCollector for linear cells, with discontinuous data.
-template <class GridView, class Partition>
-class DiscontinuousDataCollector
-    : public UnstructuredDataCollectorInterface<GridView, DiscontinuousDataCollector<GridView,Partition>, Partition>
-{
-  using Self = DiscontinuousDataCollector;
-  using Super = UnstructuredDataCollectorInterface<GridView, Self, Partition>;
-
-public:
-  using Super::dim;
-  using Super::partition;
-
-public:
-  DiscontinuousDataCollector (GridView const& gridView)
-    : Super(gridView)
-  {}
-
-  /// Create an index map the uniquely assigns an index to each pair (element,corner)
-  void updateImpl ()
+  namespace Vtk
   {
-    numPoints_ = 0;
-    numCells_ = 0;
-    indexMap_.resize(gridView_.size(dim));
-    std::int64_t vertex_idx = 0;
-    auto const& indexSet = gridView_.indexSet();
-    for (auto const& c : elements(gridView_, partition)) {
-      numCells_++;
-      numPoints_ += c.subEntities(dim);
-      for (unsigned int i = 0; i < c.subEntities(dim); ++i)
-        indexMap_[indexSet.subIndex(c, i, dim)] = vertex_idx++;
-    }
-  }
-
-  /// The number of points approx. #cell * #corners-per-cell
-  std::uint64_t numPointsImpl () const
-  {
-    return numPoints_;
-  }
 
-  /// Return the coordinates of the corners of all cells
-  template <class T>
-  std::vector<T> pointsImpl () const
-  {
-    std::vector<T> data(numPoints_ * 3);
-    auto const& indexSet = gridView_.indexSet();
-    for (auto const& element : elements(gridView_, partition)) {
-      for (unsigned int i = 0; i < element.subEntities(dim); ++i) {
-        std::size_t idx = 3 * indexMap_[indexSet.subIndex(element, i, dim)];
-        auto v = element.geometry().corner(i);
-        for (std::size_t j = 0; j < v.size(); ++j)
-          data[idx + j] = T(v[j]);
-        for (std::size_t j = v.size(); j < 3u; ++j)
-          data[idx + j] = T(0);
+    /// Implementation of \ref DataCollector for linear cells, with discontinuous data.
+    template <class GridView, class Partition>
+    class DiscontinuousDataCollector
+        : public UnstructuredDataCollectorInterface<GridView, DiscontinuousDataCollector<GridView,Partition>, Partition>
+    {
+      using Self = DiscontinuousDataCollector;
+      using Super = UnstructuredDataCollectorInterface<GridView, Self, Partition>;
+
+    public:
+      using Super::dim;
+      using Super::partition;
+
+    public:
+      DiscontinuousDataCollector (GridView const& gridView)
+        : Super(gridView)
+      {}
+
+      /// Create an index map the uniquely assigns an index to each pair (element,corner)
+      void updateImpl ()
+      {
+        numPoints_ = 0;
+        numCells_ = 0;
+        indexMap_.resize(gridView_.size(dim));
+        std::int64_t vertex_idx = 0;
+        auto const& indexSet = gridView_.indexSet();
+        for (auto const& c : elements(gridView_, partition)) {
+          numCells_++;
+          numPoints_ += c.subEntities(dim);
+          for (unsigned int i = 0; i < c.subEntities(dim); ++i)
+            indexMap_[indexSet.subIndex(c, i, dim)] = vertex_idx++;
+        }
       }
-    }
-    return data;
-  }
 
-  /// Return number of grid cells
-  std::uint64_t numCellsImpl () const
-  {
-    return numCells_;
-  }
+      /// The number of points approx. #cell * #corners-per-cell
+      std::uint64_t numPointsImpl () const
+      {
+        return numPoints_;
+      }
 
-  /// Connect the corners of each cell. The leads to a global discontinuous grid
-  Cells cellsImpl () const
-  {
-    Cells cells;
-    cells.connectivity.reserve(numPoints_);
-    cells.offsets.reserve(numCells_);
-    cells.types.reserve(numCells_);
-
-    std::int64_t old_o = 0;
-    auto const& indexSet = gridView_.indexSet();
-    for (auto const& c : elements(gridView_, partition)) {
-      Vtk::CellType cellType(c.type());
-      for (unsigned int j = 0; j < c.subEntities(dim); ++j) {
-        std::int64_t vertex_idx = indexMap_[indexSet.subIndex(c,cellType.permutation(j),dim)];
-        cells.connectivity.push_back(vertex_idx);
+      /// Return the coordinates of the corners of all cells
+      template <class T>
+      std::vector<T> pointsImpl () const
+      {
+        std::vector<T> data(numPoints_ * 3);
+        auto const& indexSet = gridView_.indexSet();
+        for (auto const& element : elements(gridView_, partition)) {
+          for (unsigned int i = 0; i < element.subEntities(dim); ++i) {
+            std::size_t idx = 3 * indexMap_[indexSet.subIndex(element, i, dim)];
+            auto v = element.geometry().corner(i);
+            for (std::size_t j = 0; j < v.size(); ++j)
+              data[idx + j] = T(v[j]);
+            for (std::size_t j = v.size(); j < 3u; ++j)
+              data[idx + j] = T(0);
+          }
+        }
+        return data;
       }
-      cells.offsets.push_back(old_o += c.subEntities(dim));
-      cells.types.push_back(cellType.type());
-    }
 
-    return cells;
-  }
+      /// Return number of grid cells
+      std::uint64_t numCellsImpl () const
+      {
+        return numCells_;
+      }
 
-  /// Evaluate the `fct` in the corners of each cell
-  template <class T, class GlobalFunction>
-  std::vector<T> pointDataImpl (GlobalFunction const& fct) const
-  {
-    std::vector<T> data(numPoints_ * fct.ncomps());
-    auto const& indexSet = gridView_.indexSet();
-    auto localFct = localFunction(fct);
-    for (auto const& e : elements(gridView_, partition)) {
-      localFct.bind(e);
-      Vtk::CellType cellType{e.type()};
-      auto refElem = referenceElement(e.geometry());
-      for (unsigned int j = 0; j < e.subEntities(dim); ++j) {
-        std::size_t idx = fct.ncomps() * indexMap_[indexSet.subIndex(e, cellType.permutation(j), dim)];
-        for (int comp = 0; comp < fct.ncomps(); ++comp)
-          data[idx + comp] = T(localFct.evaluate(comp, refElem.position(cellType.permutation(j),dim)));
+      /// Connect the corners of each cell. The leads to a global discontinuous grid
+      Cells cellsImpl () const
+      {
+        Cells cells;
+        cells.connectivity.reserve(numPoints_);
+        cells.offsets.reserve(numCells_);
+        cells.types.reserve(numCells_);
+
+        std::int64_t old_o = 0;
+        auto const& indexSet = gridView_.indexSet();
+        for (auto const& c : elements(gridView_, partition)) {
+          Vtk::CellType cellType(c.type());
+          for (unsigned int j = 0; j < c.subEntities(dim); ++j) {
+            std::int64_t vertex_idx = indexMap_[indexSet.subIndex(c,cellType.permutation(j),dim)];
+            cells.connectivity.push_back(vertex_idx);
+          }
+          cells.offsets.push_back(old_o += c.subEntities(dim));
+          cells.types.push_back(cellType.type());
+        }
+
+        return cells;
       }
-      localFct.unbind();
-    }
-    return data;
-  }
-
-protected:
-  using Super::gridView_;
-  std::uint64_t numCells_ = 0;
-  std::uint64_t numPoints_ = 0;
-  std::vector<std::int64_t> indexMap_;
-};
 
+      /// Evaluate the `fct` in the corners of each cell
+      template <class T, class GlobalFunction>
+      std::vector<T> pointDataImpl (GlobalFunction const& fct) const
+      {
+        std::vector<T> data(numPoints_ * fct.ncomps());
+        auto const& indexSet = gridView_.indexSet();
+        auto localFct = localFunction(fct);
+        for (auto const& e : elements(gridView_, partition)) {
+          localFct.bind(e);
+          Vtk::CellType cellType{e.type()};
+          auto refElem = referenceElement(e.geometry());
+          for (unsigned int j = 0; j < e.subEntities(dim); ++j) {
+            std::size_t idx = fct.ncomps() * indexMap_[indexSet.subIndex(e, cellType.permutation(j), dim)];
+            for (int comp = 0; comp < fct.ncomps(); ++comp)
+              data[idx + comp] = T(localFct.evaluate(comp, refElem.position(cellType.permutation(j),dim)));
+          }
+          localFct.unbind();
+        }
+        return data;
+      }
+
+    protected:
+      using Super::gridView_;
+      std::uint64_t numCells_ = 0;
+      std::uint64_t numPoints_ = 0;
+      std::vector<std::int64_t> indexMap_;
+    };
+
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/datacollectors/lagrangedatacollector.hh b/dune/vtk/datacollectors/lagrangedatacollector.hh
index d38f3411026548a241761bf6924c7cda87154a57..091df3de90572ac137f65af50bbb17e974ef12e6 100644
--- a/dune/vtk/datacollectors/lagrangedatacollector.hh
+++ b/dune/vtk/datacollectors/lagrangedatacollector.hh
@@ -8,189 +8,194 @@
 #include <dune/grid/common/partitionset.hh>
 
 #include <dune/vtk/forward.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
 #include <dune/vtk/utility/lagrangepoints.hh>
 
 #include "unstructureddatacollector.hh"
 
-namespace Dune {
-
-/// Implementation of \ref DataCollector for lagrange cells
-template <class GridView, int ORDER = -1>
-class LagrangeDataCollector
-    : public UnstructuredDataCollectorInterface<GridView, LagrangeDataCollector<GridView,ORDER>, Partitions::All>
+namespace Dune
 {
-  using Self = LagrangeDataCollector;
-  using Super = UnstructuredDataCollectorInterface<GridView, Self, Partitions::All>;
-
-public:
-  static_assert(ORDER != 0, "Order 0 not supported");
-  using Super::dim;
-  using Super::partition; // NOTE: lagrange data-collector currently implemented for the All partition only
-
-public:
-  LagrangeDataCollector (GridView const& gridView, int order = (ORDER < 0 ? 2 : ORDER))
-    : Super(gridView)
-    , order_(order)
-  {
-    assert(order > 0 && "Order 0 not supported");
-    assert(ORDER < 0 || order == ORDER);
-  }
 
-  /// Construct the point sets
-  void updateImpl ()
+  namespace Vtk
   {
-    auto const& indexSet = gridView_.indexSet();
 
-    pointSets_.clear();
-    for (auto gt : indexSet.types(0))
-      pointSets_.emplace(gt, order_);
+    /// Implementation of \ref DataCollector for lagrange cells
+    template <class GridView, int ORDER = -1>
+    class LagrangeDataCollector
+        : public UnstructuredDataCollectorInterface<GridView, LagrangeDataCollector<GridView,ORDER>, Partitions::All>
+    {
+      using Self = LagrangeDataCollector;
+      using Super = UnstructuredDataCollectorInterface<GridView, Self, Partitions::All>;
+
+    public:
+      static_assert(ORDER != 0, "Order 0 not supported");
+      using Super::dim;
+      using Super::partition; // NOTE: lagrange data-collector currently implemented for the All partition only
+
+    public:
+      LagrangeDataCollector (GridView const& gridView, int order = (ORDER < 0 ? 2 : ORDER))
+        : Super(gridView)
+        , order_(order)
+      {
+        assert(order > 0 && "Order 0 not supported");
+        assert(ORDER < 0 || order == ORDER);
+      }
 
-    for (auto& pointSet : pointSets_)
-      pointSet.second.build(pointSet.first);
+      /// Construct the point sets
+      void updateImpl ()
+      {
+        auto const& indexSet = gridView_.indexSet();
 
-    numPoints_ = indexSet.size(dim);
-    for (auto const& pointSet : pointSets_) {
-      auto gt = pointSet.first;
-      auto refElem = referenceElement<double,dim>(gt);
-      numPoints_ += (pointSet.second.size() - refElem.size(dim)) * indexSet.size(gt);
-    }
-  }
+        pointSets_.clear();
+        for (auto gt : indexSet.types(0))
+          pointSets_.emplace(gt, order_);
 
-  /// Return number of lagrange nodes
-  std::uint64_t numPointsImpl () const
-  {
-    return numPoints_;
-  }
-
-  /// Return a vector of point coordinates.
-  /**
-   * The vector of point coordinates is composed of vertex coordinates first and second
-   * edge center coordinates.
-   **/
-  template <class T>
-  std::vector<T> pointsImpl () const
-  {
-    std::vector<T> data(this->numPoints() * 3);
-    auto const& indexSet = gridView_.indexSet();
-
-    std::size_t shift = indexSet.size(dim);
-
-    for (auto const& element : elements(gridView_, partition)) {
-      auto geometry = element.geometry();
-      auto refElem = referenceElement<T,dim>(element.type());
-
-      auto const& pointSet = pointSets_.at(element.type());
-      unsigned int vertexDOFs = refElem.size(dim);
-      unsigned int innerDOFs = pointSet.size() - vertexDOFs;
-
-      for (std::size_t i = 0; i < pointSet.size(); ++i) {
-        auto const& p = pointSet[i];
-        if (i < vertexDOFs)
-          assert(p.localKey().codim() == dim);
-
-        auto const& localKey = p.localKey();
-        std::size_t idx = 3 * (localKey.codim() == dim
-            ? indexSet.subIndex(element, localKey.subEntity(), dim)
-            : innerDOFs*indexSet.index(element) + (i - vertexDOFs) + shift);
-
-        auto v = geometry.global(p.point());
-        for (std::size_t j = 0; j < v.size(); ++j)
-          data[idx + j] = T(v[j]);
-        for (std::size_t j = v.size(); j < 3u; ++j)
-          data[idx + j] = T(0);
+        for (auto& pointSet : pointSets_)
+          pointSet.second.build(pointSet.first);
+
+        numPoints_ = indexSet.size(dim);
+        for (auto const& pointSet : pointSets_) {
+          auto gt = pointSet.first;
+          auto refElem = referenceElement<double,dim>(gt);
+          numPoints_ += (pointSet.second.size() - refElem.size(dim)) * indexSet.size(gt);
+        }
       }
-    }
-    return data;
-  }
 
-  /// Return number of grid cells
-  std::uint64_t numCellsImpl () const
-  {
-    return gridView_.size(0);
-  }
-
-  /// \brief Return cell types, offsets, and connectivity. \see Cells
-  /**
-   * The cell connectivity is composed of cell vertices first and second cell edges,
-   * where the indices are grouped [vertex-indices..., (#vertices)+edge-indices...]
-   **/
-  Cells cellsImpl () const
-  {
-    Cells cells;
-    cells.connectivity.reserve(this->numPoints());
-    cells.offsets.reserve(this->numCells());
-    cells.types.reserve(this->numCells());
-
-    auto const& indexSet = gridView_.indexSet();
-    std::size_t shift = indexSet.size(dim);
-
-    std::int64_t old_o = 0;
-    for (auto const& element : elements(gridView_, partition)) {
-      auto refElem = referenceElement<double,dim>(element.type());
-      Vtk::CellType cellType(element.type(), Vtk::LAGRANGE);
-
-      auto const& pointSet = pointSets_.at(element.type());
-      unsigned int vertexDOFs = refElem.size(dim);
-      unsigned int innerDOFs = pointSet.size() - vertexDOFs;
-
-      for (std::size_t i = 0; i < pointSet.size(); ++i) {
-        auto const& p = pointSet[i];
-        auto const& localKey = p.localKey();
-        std::size_t idx = (localKey.codim() == dim
-            ? indexSet.subIndex(element, localKey.subEntity(), dim)
-            : innerDOFs*indexSet.index(element) + (i - vertexDOFs) + shift);
-        cells.connectivity.push_back(idx);
+      /// Return number of lagrange nodes
+      std::uint64_t numPointsImpl () const
+      {
+        return numPoints_;
       }
 
-      cells.offsets.push_back(old_o += pointSet.size());
-      cells.types.push_back(cellType.type());
-    }
-    return cells;
-  }
+      /// Return a vector of point coordinates.
+      /**
+      * The vector of point coordinates is composed of vertex coordinates first and second
+      * edge center coordinates.
+      **/
+      template <class T>
+      std::vector<T> pointsImpl () const
+      {
+        std::vector<T> data(this->numPoints() * 3);
+        auto const& indexSet = gridView_.indexSet();
+
+        std::size_t shift = indexSet.size(dim);
+
+        for (auto const& element : elements(gridView_, partition)) {
+          auto geometry = element.geometry();
+          auto refElem = referenceElement<T,dim>(element.type());
+
+          auto const& pointSet = pointSets_.at(element.type());
+          unsigned int vertexDOFs = refElem.size(dim);
+          unsigned int innerDOFs = pointSet.size() - vertexDOFs;
+
+          for (std::size_t i = 0; i < pointSet.size(); ++i) {
+            auto const& p = pointSet[i];
+            if (i < vertexDOFs)
+              assert(p.localKey().codim() == dim);
+
+            auto const& localKey = p.localKey();
+            std::size_t idx = 3 * (localKey.codim() == dim
+                ? indexSet.subIndex(element, localKey.subEntity(), dim)
+                : innerDOFs*indexSet.index(element) + (i - vertexDOFs) + shift);
+
+            auto v = geometry.global(p.point());
+            for (std::size_t j = 0; j < v.size(); ++j)
+              data[idx + j] = T(v[j]);
+            for (std::size_t j = v.size(); j < 3u; ++j)
+              data[idx + j] = T(0);
+          }
+        }
+        return data;
+      }
 
-  /// Evaluate the `fct` at element vertices and edge centers in the same order as the point coords.
-  template <class T, class GlobalFunction>
-  std::vector<T> pointDataImpl (GlobalFunction const& fct) const
-  {
-    int nComps = fct.ncomps();
-    std::vector<T> data(this->numPoints() * nComps);
-    auto const& indexSet = gridView_.indexSet();
-
-    std::size_t shift = indexSet.size(dim);
-
-    auto localFct = localFunction(fct);
-    for (auto const& element : elements(gridView_, partition)) {
-      localFct.bind(element);
-      auto refElem = referenceElement<T,dim>(element.type());
-
-      auto const& pointSet = pointSets_.at(element.type());
-      unsigned int vertexDOFs = refElem.size(dim);
-      unsigned int innerDOFs = pointSet.size() - vertexDOFs;
-
-      for (std::size_t i = 0; i < pointSet.size(); ++i) {
-        auto const& p = pointSet[i];
-        auto const& localKey = p.localKey();
-        std::size_t idx = nComps * (localKey.codim() == dim
-            ? indexSet.subIndex(element, localKey.subEntity(), dim)
-            : innerDOFs*indexSet.index(element) + (i - vertexDOFs) + shift);
-
-        for (int comp = 0; comp < nComps; ++comp)
-          data[idx + comp] = T(localFct.evaluate(comp, p.point()));
+      /// Return number of grid cells
+      std::uint64_t numCellsImpl () const
+      {
+        return gridView_.size(0);
+      }
+
+      /// \brief Return cell types, offsets, and connectivity. \see Cells
+      /**
+      * The cell connectivity is composed of cell vertices first and second cell edges,
+      * where the indices are grouped [vertex-indices..., (#vertices)+edge-indices...]
+      **/
+      Cells cellsImpl () const
+      {
+        Cells cells;
+        cells.connectivity.reserve(this->numPoints());
+        cells.offsets.reserve(this->numCells());
+        cells.types.reserve(this->numCells());
+
+        auto const& indexSet = gridView_.indexSet();
+        std::size_t shift = indexSet.size(dim);
+
+        std::int64_t old_o = 0;
+        for (auto const& element : elements(gridView_, partition)) {
+          auto refElem = referenceElement<double,dim>(element.type());
+          Vtk::CellType cellType(element.type(), Vtk::LAGRANGE);
+
+          auto const& pointSet = pointSets_.at(element.type());
+          unsigned int vertexDOFs = refElem.size(dim);
+          unsigned int innerDOFs = pointSet.size() - vertexDOFs;
+
+          for (std::size_t i = 0; i < pointSet.size(); ++i) {
+            auto const& p = pointSet[i];
+            auto const& localKey = p.localKey();
+            std::size_t idx = (localKey.codim() == dim
+                ? indexSet.subIndex(element, localKey.subEntity(), dim)
+                : innerDOFs*indexSet.index(element) + (i - vertexDOFs) + shift);
+            cells.connectivity.push_back(idx);
+          }
+
+          cells.offsets.push_back(old_o += pointSet.size());
+          cells.types.push_back(cellType.type());
+        }
+        return cells;
+      }
+
+      /// Evaluate the `fct` at element vertices and edge centers in the same order as the point coords.
+      template <class T, class GlobalFunction>
+      std::vector<T> pointDataImpl (GlobalFunction const& fct) const
+      {
+        int nComps = fct.ncomps();
+        std::vector<T> data(this->numPoints() * nComps);
+        auto const& indexSet = gridView_.indexSet();
+
+        std::size_t shift = indexSet.size(dim);
+
+        auto localFct = localFunction(fct);
+        for (auto const& element : elements(gridView_, partition)) {
+          localFct.bind(element);
+          auto refElem = referenceElement<T,dim>(element.type());
+
+          auto const& pointSet = pointSets_.at(element.type());
+          unsigned int vertexDOFs = refElem.size(dim);
+          unsigned int innerDOFs = pointSet.size() - vertexDOFs;
+
+          for (std::size_t i = 0; i < pointSet.size(); ++i) {
+            auto const& p = pointSet[i];
+            auto const& localKey = p.localKey();
+            std::size_t idx = nComps * (localKey.codim() == dim
+                ? indexSet.subIndex(element, localKey.subEntity(), dim)
+                : innerDOFs*indexSet.index(element) + (i - vertexDOFs) + shift);
+
+            for (int comp = 0; comp < nComps; ++comp)
+              data[idx + comp] = T(localFct.evaluate(comp, p.point()));
+          }
+          localFct.unbind();
+        }
+        return data;
       }
-      localFct.unbind();
-    }
-    return data;
-  }
 
-protected:
-  using Super::gridView_;
+    protected:
+      using Super::gridView_;
 
-  unsigned int order_;
-  std::uint64_t numPoints_ = 0;
+      unsigned int order_;
+      std::uint64_t numPoints_ = 0;
 
-  using PointSet = VtkLagrangePointSet<typename GridView::ctype, GridView::dimension>;
-  std::map<GeometryType, PointSet> pointSets_;
-};
+      using PointSet = LagrangePointSet<typename GridView::ctype, GridView::dimension>;
+      std::map<GeometryType, PointSet> pointSets_;
+    };
 
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/datacollectors/quadraticdatacollector.hh b/dune/vtk/datacollectors/quadraticdatacollector.hh
index 0ddee00b9028c156d017ca9f430681439b4bbb8b..90232a3b78a0baf20611e9e6ad2ea799b8a3ec55 100644
--- a/dune/vtk/datacollectors/quadraticdatacollector.hh
+++ b/dune/vtk/datacollectors/quadraticdatacollector.hh
@@ -6,140 +6,144 @@
 #include <dune/grid/common/partitionset.hh>
 
 #include <dune/vtk/forward.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
 
 #include "unstructureddatacollector.hh"
 
 namespace Dune
 {
 
-/// Implementation of \ref DataCollector for quadratic cells, with continuous data.
-template <class GridView>
-class QuadraticDataCollector
-    : public UnstructuredDataCollectorInterface<GridView, QuadraticDataCollector<GridView>, Partitions::All>
-{
-  using Self = QuadraticDataCollector;
-  using Super = UnstructuredDataCollectorInterface<GridView, Self, Partitions::All>;
-
-public:
-  using Super::dim;
-  using Super::partition; // NOTE: quadratic data-collector currently implemented for the All partition only
-
-public:
-  QuadraticDataCollector (GridView const& gridView)
-    : Super(gridView)
-  {}
-
-  /// Return number of vertices + number of edge
-  std::uint64_t numPointsImpl () const
+  namespace Vtk
   {
-    return gridView_.size(dim) + gridView_.size(dim-1);
-  }
-
-  /// Return a vector of point coordinates.
-  /**
-   * The vector of point coordinates is composed of vertex coordinates first and second
-   * edge center coordinates.
-   **/
-  template <class T>
-  std::vector<T> pointsImpl () const
-  {
-    std::vector<T> data(this->numPoints() * 3);
-    auto const& indexSet = gridView_.indexSet();
-    for (auto const& element : elements(gridView_, partition)) {
-      auto geometry = element.geometry();
-      auto refElem = referenceElement<T,dim>(element.type());
-
-      // vertices
-      for (unsigned int i = 0; i < element.subEntities(dim); ++i) {
-        std::size_t idx = 3 * indexSet.subIndex(element, i, dim);
-        auto v = geometry.global(refElem.position(i,dim));
-        for (std::size_t j = 0; j < v.size(); ++j)
-          data[idx + j] = T(v[j]);
-        for (std::size_t j = v.size(); j < 3u; ++j)
-          data[idx + j] = T(0);
-      }
-      // edge centers
-      for (unsigned int i = 0; i < element.subEntities(dim-1); ++i) {
-        std::size_t idx = 3 * (indexSet.subIndex(element, i, dim-1) + gridView_.size(dim));
-        auto v = geometry.global(refElem.position(i,dim-1));
-        for (std::size_t j = 0; j < v.size(); ++j)
-          data[idx + j] = T(v[j]);
-        for (std::size_t j = v.size(); j < 3u; ++j)
-          data[idx + j] = T(0);
+
+    /// Implementation of \ref DataCollector for quadratic cells, with continuous data.
+    template <class GridView>
+    class QuadraticDataCollector
+        : public UnstructuredDataCollectorInterface<GridView, QuadraticDataCollector<GridView>, Partitions::All>
+    {
+      using Self = QuadraticDataCollector;
+      using Super = UnstructuredDataCollectorInterface<GridView, Self, Partitions::All>;
+
+    public:
+      using Super::dim;
+      using Super::partition; // NOTE: quadratic data-collector currently implemented for the All partition only
+
+    public:
+      QuadraticDataCollector (GridView const& gridView)
+        : Super(gridView)
+      {}
+
+      /// Return number of vertices + number of edge
+      std::uint64_t numPointsImpl () const
+      {
+        return gridView_.size(dim) + gridView_.size(dim-1);
       }
-    }
-    return data;
-  }
 
-  /// Return number of grid cells
-  std::uint64_t numCellsImpl () const
-  {
-    return gridView_.size(0);
-  }
-
-  /// \brief Return cell types, offsets, and connectivity. \see Cells
-  /**
-   * The cell connectivity is composed of cell vertices first and second cell edges,
-   * where the indices are grouped [vertex-indices..., (#vertices)+edge-indices...]
-   **/
-  Cells cellsImpl () const
-  {
-    Cells cells;
-    cells.connectivity.reserve(this->numPoints());
-    cells.offsets.reserve(this->numCells());
-    cells.types.reserve(this->numCells());
-
-    std::int64_t old_o = 0;
-    auto const& indexSet = gridView_.indexSet();
-    for (auto const& c : elements(gridView_, partition)) {
-      Vtk::CellType cellType(c.type(), Vtk::QUADRATIC);
-      for (unsigned int j = 0; j < c.subEntities(dim); ++j) {
-        int k = cellType.permutation(j);
-        std::int64_t point_idx = indexSet.subIndex(c,k,dim);
-        cells.connectivity.push_back(point_idx);
+      /// Return a vector of point coordinates.
+      /**
+      * The vector of point coordinates is composed of vertex coordinates first and second
+      * edge center coordinates.
+      **/
+      template <class T>
+      std::vector<T> pointsImpl () const
+      {
+        std::vector<T> data(this->numPoints() * 3);
+        auto const& indexSet = gridView_.indexSet();
+        for (auto const& element : elements(gridView_, partition)) {
+          auto geometry = element.geometry();
+          auto refElem = referenceElement<T,dim>(element.type());
+
+          // vertices
+          for (unsigned int i = 0; i < element.subEntities(dim); ++i) {
+            std::size_t idx = 3 * indexSet.subIndex(element, i, dim);
+            auto v = geometry.global(refElem.position(i,dim));
+            for (std::size_t j = 0; j < v.size(); ++j)
+              data[idx + j] = T(v[j]);
+            for (std::size_t j = v.size(); j < 3u; ++j)
+              data[idx + j] = T(0);
+          }
+          // edge centers
+          for (unsigned int i = 0; i < element.subEntities(dim-1); ++i) {
+            std::size_t idx = 3 * (indexSet.subIndex(element, i, dim-1) + gridView_.size(dim));
+            auto v = geometry.global(refElem.position(i,dim-1));
+            for (std::size_t j = 0; j < v.size(); ++j)
+              data[idx + j] = T(v[j]);
+            for (std::size_t j = v.size(); j < 3u; ++j)
+              data[idx + j] = T(0);
+          }
+        }
+        return data;
       }
-      for (unsigned int j = 0; j < c.subEntities(dim-1); ++j) {
-        int k = cellType.permutation(c.subEntities(dim) + j);
-        std::int64_t point_idx = (indexSet.subIndex(c,k,dim-1) + gridView_.size(dim));
-        cells.connectivity.push_back(point_idx);
+
+      /// Return number of grid cells
+      std::uint64_t numCellsImpl () const
+      {
+        return gridView_.size(0);
       }
-      cells.offsets.push_back(old_o += c.subEntities(dim)+c.subEntities(dim-1));
-      cells.types.push_back(cellType.type());
-    }
-    return cells;
-  }
-
-  /// Evaluate the `fct` at element vertices and edge centers in the same order as the point coords.
-  template <class T, class GlobalFunction>
-  std::vector<T> pointDataImpl (GlobalFunction const& fct) const
-  {
-    std::vector<T> data(this->numPoints() * fct.ncomps());
-    auto const& indexSet = gridView_.indexSet();
-    auto localFct = localFunction(fct);
-    for (auto const& e : elements(gridView_, partition)) {
-      localFct.bind(e);
-      Vtk::CellType cellType{e.type(), Vtk::QUADRATIC};
-      auto refElem = referenceElement(e.geometry());
-      for (unsigned int j = 0; j < e.subEntities(dim); ++j) {
-        int k = cellType.permutation(j);
-        std::size_t idx = fct.ncomps() * indexSet.subIndex(e, k, dim);
-        for (int comp = 0; comp < fct.ncomps(); ++comp)
-          data[idx + comp] = T(localFct.evaluate(comp, refElem.position(k, dim)));
+
+      /// \brief Return cell types, offsets, and connectivity. \see Cells
+      /**
+      * The cell connectivity is composed of cell vertices first and second cell edges,
+      * where the indices are grouped [vertex-indices..., (#vertices)+edge-indices...]
+      **/
+      Cells cellsImpl () const
+      {
+        Cells cells;
+        cells.connectivity.reserve(this->numPoints());
+        cells.offsets.reserve(this->numCells());
+        cells.types.reserve(this->numCells());
+
+        std::int64_t old_o = 0;
+        auto const& indexSet = gridView_.indexSet();
+        for (auto const& c : elements(gridView_, partition)) {
+          Vtk::CellType cellType(c.type(), Vtk::QUADRATIC);
+          for (unsigned int j = 0; j < c.subEntities(dim); ++j) {
+            int k = cellType.permutation(j);
+            std::int64_t point_idx = indexSet.subIndex(c,k,dim);
+            cells.connectivity.push_back(point_idx);
+          }
+          for (unsigned int j = 0; j < c.subEntities(dim-1); ++j) {
+            int k = cellType.permutation(c.subEntities(dim) + j);
+            std::int64_t point_idx = (indexSet.subIndex(c,k,dim-1) + gridView_.size(dim));
+            cells.connectivity.push_back(point_idx);
+          }
+          cells.offsets.push_back(old_o += c.subEntities(dim)+c.subEntities(dim-1));
+          cells.types.push_back(cellType.type());
+        }
+        return cells;
       }
-      for (unsigned int j = 0; j < e.subEntities(dim-1); ++j) {
-        int k = cellType.permutation(e.subEntities(dim) + j);
-        std::size_t idx = fct.ncomps() * (indexSet.subIndex(e, k, dim-1) + gridView_.size(dim));
-        for (int comp = 0; comp < fct.ncomps(); ++comp)
-          data[idx + comp] = T(localFct.evaluate(comp, refElem.position(k, dim-1)));
+
+      /// Evaluate the `fct` at element vertices and edge centers in the same order as the point coords.
+      template <class T, class GlobalFunction>
+      std::vector<T> pointDataImpl (GlobalFunction const& fct) const
+      {
+        std::vector<T> data(this->numPoints() * fct.ncomps());
+        auto const& indexSet = gridView_.indexSet();
+        auto localFct = localFunction(fct);
+        for (auto const& e : elements(gridView_, partition)) {
+          localFct.bind(e);
+          Vtk::CellType cellType{e.type(), Vtk::QUADRATIC};
+          auto refElem = referenceElement(e.geometry());
+          for (unsigned int j = 0; j < e.subEntities(dim); ++j) {
+            int k = cellType.permutation(j);
+            std::size_t idx = fct.ncomps() * indexSet.subIndex(e, k, dim);
+            for (int comp = 0; comp < fct.ncomps(); ++comp)
+              data[idx + comp] = T(localFct.evaluate(comp, refElem.position(k, dim)));
+          }
+          for (unsigned int j = 0; j < e.subEntities(dim-1); ++j) {
+            int k = cellType.permutation(e.subEntities(dim) + j);
+            std::size_t idx = fct.ncomps() * (indexSet.subIndex(e, k, dim-1) + gridView_.size(dim));
+            for (int comp = 0; comp < fct.ncomps(); ++comp)
+              data[idx + comp] = T(localFct.evaluate(comp, refElem.position(k, dim-1)));
+          }
+          localFct.unbind();
+        }
+        return data;
       }
-      localFct.unbind();
-    }
-    return data;
-  }
 
-protected:
-  using Super::gridView_;
-};
+    protected:
+      using Super::gridView_;
+    };
 
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/datacollectors/spdatacollector.hh b/dune/vtk/datacollectors/spdatacollector.hh
index 25ab8be2ab1d222f8db359ac122e819ff712a010..462ddf50bf756dd8247eff2c6a0fa36f517cb175 100644
--- a/dune/vtk/datacollectors/spdatacollector.hh
+++ b/dune/vtk/datacollectors/spdatacollector.hh
@@ -16,86 +16,92 @@
 
 namespace Dune
 {
-#if HAVE_DUNE_SPGRID
-
-// Specialization for SPGrid
-template <class GridView>
-class SPDataCollector
-    : public StructuredDataCollectorInterface<GridView, SPDataCollector<GridView>>
-{
-  using Self = SPDataCollector;
-  using Super = StructuredDataCollectorInterface<GridView, Self>;
-  using ctype = typename GridView::ctype;
-
-public:
-  using Super::dim;
-  using Super::partition;
-
-public:
-  SPDataCollector (GridView const& gridView)
-    : Super(gridView)
-    , wholeExtent_(filledArray<6,int>(0))
-    , extent_(filledArray<6,int>(0))
-    , origin_(0.0)
-    , spacing_(0.0)
-  {}
-
-  std::array<int, 6> const& wholeExtentImpl () const
-  {
-    return wholeExtent_;
-  }
-
-  std::array<int, 6> const& extentImpl () const
-  {
-    return extent_;
-  }
 
-  auto const& originImpl () const
+  namespace Vtk
   {
-    return origin_;
-  }
 
-  auto const& spacingImpl () const
-  {
-    return spacing_;
-  }
-
-  void updateImpl ()
-  {
-    Super::updateImpl();
-
-    auto const& localMesh = gridView_.impl().gridLevel().localMesh();
-    auto const& begin = gridView_.impl().gridLevel().globalMesh().begin();
-    auto const& end = gridView_.impl().gridLevel().globalMesh().end();
-    auto const& cube = gridView_.grid().domain().cube();
-
-    for (int i = 0; i < dim; ++i) {
-      wholeExtent_[2*i] = begin[i];
-      wholeExtent_[2*i+1] = end[i];
-      extent_[2*i] = localMesh.begin()[i];
-      extent_[2*i+1] = localMesh.end()[i];
-      spacing_[i] = cube.width()[i] / (end[i] - begin[i]);
-      origin_[i] = cube.origin()[i];
+    #if HAVE_DUNE_SPGRID
+
+    // Specialization for SPGrid
+    template <class GridView>
+    class SPDataCollector
+        : public StructuredDataCollectorInterface<GridView, SPDataCollector<GridView>>
+    {
+      using Self = SPDataCollector;
+      using Super = StructuredDataCollectorInterface<GridView, Self>;
+      using ctype = typename GridView::ctype;
+
+    public:
+      using Super::dim;
+      using Super::partition;
+
+    public:
+      SPDataCollector (GridView const& gridView)
+        : Super(gridView)
+        , wholeExtent_(filledArray<6,int>(0))
+        , extent_(filledArray<6,int>(0))
+        , origin_(0.0)
+        , spacing_(0.0)
+      {}
+
+      std::array<int, 6> const& wholeExtentImpl () const
+      {
+        return wholeExtent_;
+      }
+
+      std::array<int, 6> const& extentImpl () const
+      {
+        return extent_;
+      }
+
+      auto const& originImpl () const
+      {
+        return origin_;
+      }
+
+      auto const& spacingImpl () const
+      {
+        return spacing_;
+      }
+
+      void updateImpl ()
+      {
+        Super::updateImpl();
+
+        auto const& localMesh = gridView_.impl().gridLevel().localMesh();
+        auto const& begin = gridView_.impl().gridLevel().globalMesh().begin();
+        auto const& end = gridView_.impl().gridLevel().globalMesh().end();
+        auto const& cube = gridView_.grid().domain().cube();
+
+        for (int i = 0; i < dim; ++i) {
+          wholeExtent_[2*i] = begin[i];
+          wholeExtent_[2*i+1] = end[i];
+          extent_[2*i] = localMesh.begin()[i];
+          extent_[2*i+1] = localMesh.end()[i];
+          spacing_[i] = cube.width()[i] / (end[i] - begin[i]);
+          origin_[i] = cube.origin()[i];
+        }
+      }
+
+    protected:
+      using Super::gridView_;
+      std::array<int, 6> wholeExtent_;
+      std::array<int, 6> extent_;
+      FieldVector<ctype,3> origin_;
+      FieldVector<ctype,3> spacing_;
+      std::vector<std::size_t> indexMap_;
+    };
+
+    namespace Impl
+    {
+      template <class GridView, class ct, int dim, template< int > class Ref, class Comm>
+      struct StructuredDataCollectorImpl<GridView, SPGrid<ct,dim,Ref,Comm>>
+      {
+        using type = SPDataCollector<GridView>;
+      };
     }
-  }
-
-protected:
-  using Super::gridView_;
-  std::array<int, 6> wholeExtent_;
-  std::array<int, 6> extent_;
-  FieldVector<ctype,3> origin_;
-  FieldVector<ctype,3> spacing_;
-  std::vector<std::size_t> indexMap_;
-};
-
-namespace Impl
-{
-  template <class GridView, class ct, int dim, template< int > class Ref, class Comm>
-  struct StructuredDataCollectorImpl<GridView, SPGrid<ct,dim,Ref,Comm>>
-  {
-    using type = SPDataCollector<GridView>;
-  };
-}
 
-#endif // HAVE_DUNE_SPGRID
+    #endif // HAVE_DUNE_SPGRID
+
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/datacollectors/structureddatacollector.hh b/dune/vtk/datacollectors/structureddatacollector.hh
index 0b1d5ab5d56f42f7a5aaa932a192c1eb0cb44e6d..6222710cca4eff53eee0729f8d0113b2ab966c34 100644
--- a/dune/vtk/datacollectors/structureddatacollector.hh
+++ b/dune/vtk/datacollectors/structureddatacollector.hh
@@ -12,225 +12,230 @@
 
 namespace Dune
 {
-/// The Interface for structured data-collectors
-template <class GridView, class Derived>
-class StructuredDataCollectorInterface
-    : public DataCollectorInterface<GridView, Derived>
-{
-protected:
-  using Super = DataCollectorInterface<GridView, Derived>;
-  using SubDataCollector = ContinuousDataCollector<GridView>;
-  using ctype = typename GridView::ctype;
-
-public:
-  using Super::dim;
-  using Super::partition;
-
-public:
-  StructuredDataCollectorInterface (GridView const& gridView)
-    : Super(gridView)
-    , subDataCollector_(gridView)
-  {}
-
-  /// Sequence of Index pairs [begin, end) for the cells in each direction
-  std::array<int, 6> wholeExtent () const
-  {
-    return this->asDerived().wholeExtentImpl();
-  }
 
-  /// Sequence of Index pairs [begin, end) for the cells in each direction of the local partition
-  std::array<int, 6> extent () const
-  {
-    return this->asDerived().extentImpl();
-  }
+  namespace Vtk
+  {
+
+    /// The Interface for structured data-collectors
+    template <class GridView, class Derived>
+    class StructuredDataCollectorInterface
+        : public DataCollectorInterface<GridView, Derived>
+    {
+    protected:
+      using Super = DataCollectorInterface<GridView, Derived>;
+      using SubDataCollector = ContinuousDataCollector<GridView>;
+      using ctype = typename GridView::ctype;
+
+    public:
+      using Super::dim;
+      using Super::partition;
+
+    public:
+      StructuredDataCollectorInterface (GridView const& gridView)
+        : Super(gridView)
+        , subDataCollector_(gridView)
+      {}
+
+      /// Sequence of Index pairs [begin, end) for the cells in each direction
+      std::array<int, 6> wholeExtent () const
+      {
+        return this->asDerived().wholeExtentImpl();
+      }
 
-  /// Call the `writer` with extent
-  template <class Writer>
-  void writeLocalPiece (Writer const& writer) const
-  {
-    this->asDerived().writeLocalPieceImpl(writer);
-  }
+      /// Sequence of Index pairs [begin, end) for the cells in each direction of the local partition
+      std::array<int, 6> extent () const
+      {
+        return this->asDerived().extentImpl();
+      }
 
-  /// Call the `writer` with piece number and piece extent
-  template <class Writer>
-  void writePieces (Writer const& writer) const
-  {
-    this->asDerived().writePiecesImpl(writer);
-  }
+      /// Call the `writer` with extent
+      template <class Writer>
+      void writeLocalPiece (Writer const& writer) const
+      {
+        this->asDerived().writeLocalPieceImpl(writer);
+      }
 
-  /// Interface for ImageData:
-  /// @{
+      /// Call the `writer` with piece number and piece extent
+      template <class Writer>
+      void writePieces (Writer const& writer) const
+      {
+        this->asDerived().writePiecesImpl(writer);
+      }
 
-  /// Lower left corner of the grid
-  FieldVector<ctype, 3> origin () const
-  {
-    return this->asDerived().originImpl();
-  }
+      /// Interface for ImageData:
+      /// @{
 
-  /// Constant grid spacing in each coordinate direction
-  FieldVector<ctype, 3> spacing () const
-  {
-    return this->asDerived().spacingImpl();
-  }
+      /// Lower left corner of the grid
+      FieldVector<ctype, 3> origin () const
+      {
+        return this->asDerived().originImpl();
+      }
 
-  /// @}
+      /// Constant grid spacing in each coordinate direction
+      FieldVector<ctype, 3> spacing () const
+      {
+        return this->asDerived().spacingImpl();
+      }
 
-  /// Interface for RectilinearGrid
-  /// @{
+      /// @}
 
-  /// The coordinates defines point coordinates for an extent by specifying the ordinate along each axis.
-  template <class T>
-  std::array<std::vector<T>, 3> coordinates () const
-  {
-    return this->asDerived().template coordinatesImpl<T>();
-  }
+      /// Interface for RectilinearGrid
+      /// @{
 
-  /// @}
+      /// The coordinates defines point coordinates for an extent by specifying the ordinate along each axis.
+      template <class T>
+      std::array<std::vector<T>, 3> coordinates () const
+      {
+        return this->asDerived().template coordinatesImpl<T>();
+      }
 
+      /// @}
 
-public: // default implementation:
 
-  /// \copyref DefaultDataCollector::update.
-  void updateImpl ()
-  {
-    subDataCollector_.update();
+    public: // default implementation:
 
-#if HAVE_MPI
-    int rank = gridView_.comm().rank();
-    int numRanks = gridView_.comm().size();
+      /// \copyref DefaultDataCollector::update.
+      void updateImpl ()
+      {
+        subDataCollector_.update();
 
-    if (rank == 0) {
-      extents_.resize(numRanks);
-      requests_.resize(numRanks, MPI_REQUEST_NULL);
-      for (int i = 1; i < numRanks; ++i)
-        MPI_Irecv(extents_[i].data(), extents_[i].size(), MPI_INT, i, /*tag=*/6, gridView_.comm(), &requests_[i]);
-    }
+    #if HAVE_MPI
+        int rank = gridView_.comm().rank();
+        int numRanks = gridView_.comm().size();
 
-    sendRequest_ = MPI_REQUEST_NULL;
-#endif
-  }
+        if (rank == 0) {
+          extents_.resize(numRanks);
+          requests_.resize(numRanks, MPI_REQUEST_NULL);
+          for (int i = 1; i < numRanks; ++i)
+            MPI_Irecv(extents_[i].data(), extents_[i].size(), MPI_INT, i, /*tag=*/6, gridView_.comm(), &requests_[i]);
+        }
 
-  /// Return number of grid cells
-  std::uint64_t numCellsImpl () const
-  {
-    auto extent = this->extent();
-    std::uint64_t num = 1;
-    for (int d = 0; d < dim; ++d)
-      num *= extent[2*d+1] - extent[2*d];
-    return num;
-  }
-
-  /// Return number of grid vertices
-  std::uint64_t numPointsImpl () const
-  {
-    return subDataCollector_.numPoints();
-  }
+        sendRequest_ = MPI_REQUEST_NULL;
+    #endif
+      }
 
-  /// \copydoc DefaultDataCollector::points.
-  template <class T>
-  std::vector<T> pointsImpl () const
-  {
-    return subDataCollector_.template points<T>();
-  }
+      /// Return number of grid cells
+      std::uint64_t numCellsImpl () const
+      {
+        auto extent = this->extent();
+        std::uint64_t num = 1;
+        for (int d = 0; d < dim; ++d)
+          num *= extent[2*d+1] - extent[2*d];
+        return num;
+      }
 
-  /// \copydoc DefaultDataCollector::pointData
-  template <class T, class GlobalFunction>
-  std::vector<T> pointDataImpl (GlobalFunction const& fct) const
-  {
-    return subDataCollector_.template pointData<T>(fct);
-  }
+      /// Return number of grid vertices
+      std::uint64_t numPointsImpl () const
+      {
+        return subDataCollector_.numPoints();
+      }
 
-  // Calculates the extent and communicates it to rank 0.
-  template <class Writer>
-  void writeLocalPieceImpl (Writer const& writer) const
-  {
-    auto&& extent = this->extent();
+      /// \copydoc DefaultDataCollector::points.
+      template <class T>
+      std::vector<T> pointsImpl () const
+      {
+        return subDataCollector_.template points<T>();
+      }
+
+      /// \copydoc DefaultDataCollector::pointData
+      template <class T, class GlobalFunction>
+      std::vector<T> pointDataImpl (GlobalFunction const& fct) const
+      {
+        return subDataCollector_.template pointData<T>(fct);
+      }
 
-#if HAVE_MPI
-    int sendFlag = 0;
-    MPI_Status sendStatus;
-    MPI_Test(&sendRequest_, &sendFlag, &sendStatus);
+      // Calculates the extent and communicates it to rank 0.
+      template <class Writer>
+      void writeLocalPieceImpl (Writer const& writer) const
+      {
+        auto&& extent = this->extent();
+
+    #if HAVE_MPI
+        int sendFlag = 0;
+        MPI_Status sendStatus;
+        MPI_Test(&sendRequest_, &sendFlag, &sendStatus);
+
+        if (sendFlag) {
+          int rank = gridView_.comm().rank();
+          if (rank != 0) {
+            MPI_Isend(extent.data(), extent.size(), MPI_INT, 0, /*tag=*/6, gridView_.comm(), &sendRequest_);
+          } else {
+            extents_[0] = extent;
+          }
+        }
+    #endif
+
+        writer(extent);
+      }
 
-    if (sendFlag) {
-      int rank = gridView_.comm().rank();
-      if (rank != 0) {
-        MPI_Isend(extent.data(), extent.size(), MPI_INT, 0, /*tag=*/6, gridView_.comm(), &sendRequest_);
-      } else {
-        extents_[0] = extent;
+      // Receive extent from all ranks and call the `writer` with the rank's extent vector
+      template <class Writer>
+      void writePiecesImpl (Writer const& writer) const
+      {
+    #if HAVE_MPI
+        writer(0, extents_[0], true);
+
+        int numRanks = gridView_.comm().size();
+        for (int p = 1; p < numRanks; ++p) {
+          int idx = -1;
+          MPI_Status status;
+          MPI_Waitany(numRanks, requests_.data(), &idx, &status);
+          if (idx != MPI_UNDEFINED) {
+            assert(idx == status.MPI_SOURCE && status.MPI_TAG == 6);
+            writer(idx, extents_[idx], true);
+          }
+        }
+    #else
+        writer(0, this->extent(), true);
+    #endif
       }
-    }
-#endif
 
-    writer(extent);
-  }
+      // Origin (0,0,0)
+      FieldVector<ctype, 3> originImpl () const
+      {
+        FieldVector<ctype, 3> vec; vec = ctype(0);
+        return vec;
+      }
 
-  // Receive extent from all ranks and call the `writer` with the rank's extent vector
-  template <class Writer>
-  void writePiecesImpl (Writer const& writer) const
-  {
-#if HAVE_MPI
-    writer(0, extents_[0], true);
-
-    int numRanks = gridView_.comm().size();
-    for (int p = 1; p < numRanks; ++p) {
-      int idx = -1;
-      MPI_Status status;
-      MPI_Waitany(numRanks, requests_.data(), &idx, &status);
-      if (idx != MPI_UNDEFINED) {
-        assert(idx == status.MPI_SOURCE && status.MPI_TAG == 6);
-        writer(idx, extents_[idx], true);
-      }
-    }
-#else
-    writer(0, this->extent(), true);
-#endif
-  }
-
-  // Origin (0,0,0)
-  FieldVector<ctype, 3> originImpl () const
-  {
-    FieldVector<ctype, 3> vec; vec = ctype(0);
-    return vec;
-  }
+      // Grid spacing (0,0,0)
+      FieldVector<ctype, 3> spacingImpl () const
+      {
+        FieldVector<ctype, 3> vec; vec = ctype(0);
+        return vec;
+      }
 
-  // Grid spacing (0,0,0)
-  FieldVector<ctype, 3> spacingImpl () const
-  {
-    FieldVector<ctype, 3> vec; vec = ctype(0);
-    return vec;
-  }
+      // Ordinate along each axis with constant \ref spacing from the \ref origin
+      template <class T>
+      std::array<std::vector<T>, 3> coordinatesImpl () const
+      {
+        auto origin = this->origin();
+        auto spacing = this->spacing();
+        auto extent = this->extent();
+
+        std::array<std::vector<T>, 3> ordinates{};
+        for (int d = 0; d < dim; ++d) {
+          auto s = extent[2*d+1] - extent[2*d] + 1;
+          ordinates[d].resize(s);
+          for (int i = 0; i < s; ++i)
+            ordinates[d][i] = origin[d] + (extent[2*d] + i)*spacing[d];
+        }
+
+        for (int d = dim; d < 3; ++d)
+          ordinates[d].resize(1, T(0));
+
+        return ordinates;
+      }
 
-  // Ordinate along each axis with constant \ref spacing from the \ref origin
-  template <class T>
-  std::array<std::vector<T>, 3> coordinatesImpl () const
-  {
-    auto origin = this->origin();
-    auto spacing = this->spacing();
-    auto extent = this->extent();
-
-    std::array<std::vector<T>, 3> ordinates{};
-    for (int d = 0; d < dim; ++d) {
-      auto s = extent[2*d+1] - extent[2*d] + 1;
-      ordinates[d].resize(s);
-      for (int i = 0; i < s; ++i)
-        ordinates[d][i] = origin[d] + (extent[2*d] + i)*spacing[d];
-    }
-
-    for (int d = dim; d < 3; ++d)
-      ordinates[d].resize(1, T(0));
-
-    return ordinates;
-  }
-
-protected:
-  using Super::gridView_;
-  SubDataCollector subDataCollector_;
-
-#if HAVE_MPI
-  mutable std::vector<std::array<int,6>> extents_;
-  mutable std::vector<MPI_Request> requests_;
-  mutable MPI_Request sendRequest_ = MPI_REQUEST_NULL;
-#endif
-};
+    protected:
+      using Super::gridView_;
+      SubDataCollector subDataCollector_;
+
+    #if HAVE_MPI
+      mutable std::vector<std::array<int,6>> extents_;
+      mutable std::vector<MPI_Request> requests_;
+      mutable MPI_Request sendRequest_ = MPI_REQUEST_NULL;
+    #endif
+    };
 
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/datacollectors/test/test-lagrangedatacollector.cc b/dune/vtk/datacollectors/test/test-lagrangedatacollector.cc
index 99d22c5fca646eebb7005f5d2f99bb847790ee6b..32d6870db952f757fd1aecd6cee108ba1cc3e625 100644
--- a/dune/vtk/datacollectors/test/test-lagrangedatacollector.cc
+++ b/dune/vtk/datacollectors/test/test-lagrangedatacollector.cc
@@ -18,7 +18,7 @@
 #include <dune/vtk/datacollectors/lagrangedatacollector.hh>
 
 /**
- * This test compares LagrangeDataCollector<1> against the ContinuousDataCollector, thus 
+ * This test compares LagrangeDataCollector<1> against the ContinuousDataCollector, thus
  * just compares the  linear order implementations
  **/
 
@@ -34,7 +34,7 @@ void testDataCollector (std::string prefix, Dune::TestSuite& testSuite, DataColl
   // check that number of points and cells are equal
   testSuite.check(dataCollector1.numPoints() == dataCollector2.numPoints(), prefix + "_numPoints");
   testSuite.check(dataCollector1.numCells() == dataCollector2.numCells(), prefix + "_numCells");
-  
+
   auto points1 = dataCollector1.template points<double>();
   auto points2 = dataCollector2.template points<double>();
 
@@ -71,14 +71,14 @@ void testGridView (std::string prefix, Dune::TestSuite& testSuite, GridView cons
 {
   // 1. test linear order lagrange data-collector
   {
-    Dune::ContinuousDataCollector<GridView> linearDataCollector(gridView);
-    
+    Dune::Vtk::ContinuousDataCollector<GridView> linearDataCollector(gridView);
+
     // data collector with template order
-    Dune::LagrangeDataCollector<GridView, 1> lagrangeDataCollector1a(gridView);
+    Dune::Vtk::LagrangeDataCollector<GridView, 1> lagrangeDataCollector1a(gridView);
     testDataCollector(prefix + "_linear_template", testSuite, lagrangeDataCollector1a, linearDataCollector);
 
     // data collector with runtime order
-    Dune::LagrangeDataCollector<GridView> lagrangeDataCollector1b(gridView, 1);
+    Dune::Vtk::LagrangeDataCollector<GridView> lagrangeDataCollector1b(gridView, 1);
     testDataCollector(prefix + "_linear_runtime", testSuite, lagrangeDataCollector1b, linearDataCollector);
   }
 }
diff --git a/dune/vtk/datacollectors/unstructureddatacollector.hh b/dune/vtk/datacollectors/unstructureddatacollector.hh
index b5ff27a5dc786649ecb719887e3d08025b369405..3e7b5a1235e626c61e3388649d1a5b12a32d844b 100644
--- a/dune/vtk/datacollectors/unstructureddatacollector.hh
+++ b/dune/vtk/datacollectors/unstructureddatacollector.hh
@@ -6,50 +6,55 @@
 #include <dune/vtk/forward.hh>
 #include <dune/vtk/datacollectorinterface.hh>
 
-namespace Dune {
-
-struct Cells
-{
-  std::vector<std::uint8_t> types;
-  std::vector<std::int64_t> offsets;
-  std::vector<std::int64_t> connectivity;
-};
-
-template <class GridView, class Derived, class Partition>
-class UnstructuredDataCollectorInterface
-    : public DataCollectorInterface<GridView, Derived, Partition>
+namespace Dune
 {
-  using Super = DataCollectorInterface<GridView, Derived, Partition>;
-
-public:
-  using Super::dim;
-  using Super::partition;
 
-public:
-  UnstructuredDataCollectorInterface (GridView const& gridView)
-    : Super(gridView)
-  {}
-
-  /// \brief Return cell types, offsets, and connectivity. \see Cells
-  Cells cells () const
+  namespace Vtk
   {
-    return this->asDerived().cellsImpl();
-  }
-
-  std::vector<std::uint64_t> pointIds () const
-  {
-    return this->asDerived().pointIdsImpl();
-  }
-
-protected:
-  // default implementation
-  std::vector<std::uint64_t> pointIdsImpl () const
-  {
-    return {};
-  }
-
-protected:
-  using Super::gridView_;
-};
 
+    struct Cells
+    {
+      std::vector<std::uint8_t> types;
+      std::vector<std::int64_t> offsets;
+      std::vector<std::int64_t> connectivity;
+    };
+
+    template <class GridView, class Derived, class Partition>
+    class UnstructuredDataCollectorInterface
+        : public DataCollectorInterface<GridView, Derived, Partition>
+    {
+      using Super = DataCollectorInterface<GridView, Derived, Partition>;
+
+    public:
+      using Super::dim;
+      using Super::partition;
+
+    public:
+      UnstructuredDataCollectorInterface (GridView const& gridView)
+        : Super(gridView)
+      {}
+
+      /// \brief Return cell types, offsets, and connectivity. \see Cells
+      Cells cells () const
+      {
+        return this->asDerived().cellsImpl();
+      }
+
+      std::vector<std::uint64_t> pointIds () const
+      {
+        return this->asDerived().pointIdsImpl();
+      }
+
+    protected:
+      // default implementation
+      std::vector<std::uint64_t> pointIdsImpl () const
+      {
+        return {};
+      }
+
+    protected:
+      using Super::gridView_;
+    };
+
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/datacollectors/yaspdatacollector.hh b/dune/vtk/datacollectors/yaspdatacollector.hh
index 818f200f1a5410295ac198d8ae941acee681a203..f246854594e14995d55061d457998b7b342ffda0 100644
--- a/dune/vtk/datacollectors/yaspdatacollector.hh
+++ b/dune/vtk/datacollectors/yaspdatacollector.hh
@@ -17,153 +17,158 @@
 
 namespace Dune
 {
-// Specialization for YaspGrid
-template <class GridView>
-class YaspDataCollector
-    : public StructuredDataCollectorInterface<GridView, YaspDataCollector<GridView>>
-{
-  using Self = YaspDataCollector;
-  using Super = StructuredDataCollectorInterface<GridView, Self>;
-  using ctype = typename GridView::ctype;
-
-public:
-  using Super::dim;
-  using Super::partition;
-
-public:
-  YaspDataCollector (GridView const& gridView)
-    : Super(gridView)
-    , wholeExtent_(filledArray<6,int>(0))
-    , extent_(filledArray<6,int>(0))
-    , origin_(0.0)
-    , spacing_(0.0)
-    , level_(0)
-  {}
-
-  std::array<int, 6> const& wholeExtentImpl () const
-  {
-    return wholeExtent_;
-  }
-
-  std::array<int, 6> const& extentImpl () const
-  {
-    return extent_;
-  }
-
-  auto const& originImpl () const
-  {
-    return origin_;
-  }
-
-  auto const& spacingImpl () const
-  {
-    return spacing_;
-  }
-
-  void updateImpl ()
-  {
-    Super::updateImpl();
-
-    level_ = gridView_.template begin<0,All_Partition>()->level();
-    for (int i = 0; i < dim; ++i) {
-      wholeExtent_[2*i] = 0;
-      wholeExtent_[2*i+1] = grid(gridView_.grid()).levelSize(level_,i);
-    }
-
-    auto const& gl = *grid(gridView_.grid()).begin(level_);
-    auto const& g = gl.interior[0];
-    auto const& gc = *g.dataBegin();
-    for (int i = 0; i < dim; ++i) {
-      extent_[2*i] = gc.min(i);
-      extent_[2*i+1] = gc.max(i)+1;
-    }
 
-    auto it = grid(gridView_.grid()).begin(level_);
-    initGeometry(it->coords);
-  }
-
-  void initGeometry (EquidistantCoordinates<ctype,dim> const& coords)
+  namespace Vtk
   {
-    for (int i = 0; i < dim; ++i) {
-      spacing_[i] = coords.meshsize(i,0);
-      origin_[i] = 0;
-    }
-  }
 
-  void initGeometry (EquidistantOffsetCoordinates<ctype,dim> const& coords)
-  {
-    for (int i = 0; i < dim; ++i) {
-      spacing_[i] = coords.meshsize(i,0);
-      origin_[i] = coords.origin(i);
+    // Specialization for YaspGrid
+    template <class GridView>
+    class YaspDataCollector
+        : public StructuredDataCollectorInterface<GridView, YaspDataCollector<GridView>>
+    {
+      using Self = YaspDataCollector;
+      using Super = StructuredDataCollectorInterface<GridView, Self>;
+      using ctype = typename GridView::ctype;
+
+    public:
+      using Super::dim;
+      using Super::partition;
+
+    public:
+      YaspDataCollector (GridView const& gridView)
+        : Super(gridView)
+        , wholeExtent_(filledArray<6,int>(0))
+        , extent_(filledArray<6,int>(0))
+        , origin_(0.0)
+        , spacing_(0.0)
+        , level_(0)
+      {}
+
+      std::array<int, 6> const& wholeExtentImpl () const
+      {
+        return wholeExtent_;
+      }
+
+      std::array<int, 6> const& extentImpl () const
+      {
+        return extent_;
+      }
+
+      auto const& originImpl () const
+      {
+        return origin_;
+      }
+
+      auto const& spacingImpl () const
+      {
+        return spacing_;
+      }
+
+      void updateImpl ()
+      {
+        Super::updateImpl();
+
+        level_ = gridView_.template begin<0,All_Partition>()->level();
+        for (int i = 0; i < dim; ++i) {
+          wholeExtent_[2*i] = 0;
+          wholeExtent_[2*i+1] = grid(gridView_.grid()).levelSize(level_,i);
+        }
+
+        auto const& gl = *grid(gridView_.grid()).begin(level_);
+        auto const& g = gl.interior[0];
+        auto const& gc = *g.dataBegin();
+        for (int i = 0; i < dim; ++i) {
+          extent_[2*i] = gc.min(i);
+          extent_[2*i+1] = gc.max(i)+1;
+        }
+
+        auto it = grid(gridView_.grid()).begin(level_);
+        initGeometry(it->coords);
+      }
+
+      void initGeometry (EquidistantCoordinates<ctype,dim> const& coords)
+      {
+        for (int i = 0; i < dim; ++i) {
+          spacing_[i] = coords.meshsize(i,0);
+          origin_[i] = 0;
+        }
+      }
+
+      void initGeometry (EquidistantOffsetCoordinates<ctype,dim> const& coords)
+      {
+        for (int i = 0; i < dim; ++i) {
+          spacing_[i] = coords.meshsize(i,0);
+          origin_[i] = coords.origin(i);
+        }
+      }
+
+      void initGeometry (TensorProductCoordinates<ctype,dim> const& coords)
+      {
+        for (int i = 0; i < dim; ++i) {
+          spacing_[i] = coords.meshsize(i,0); // is not constant, but also not used.
+          origin_[i] = coords.coordinate(i,0);
+        }
+      }
+
+
+      /// Extract the ordinates from the coordinates object of the current level
+      template <class T>
+      std::array<std::vector<T>, 3> coordinatesImpl () const
+      {
+        auto it = grid(gridView_.grid()).begin(level_);
+        auto const& coords = it->coords;
+
+        std::array<std::vector<T>, 3> ordinates{};
+        for (int d = 0; d < dim; ++d) {
+          auto s = extent_[2*d+1] - extent_[2*d] + 1;
+          ordinates[d].resize(s);
+          for (int i = 0; i < s; ++i)
+            ordinates[d][i] = coords.coordinate(d, extent_[2*d] + i);
+        }
+
+        for (int d = dim; d < 3; ++d)
+          ordinates[d].resize(1, T(0));
+
+        return ordinates;
+      }
+
+
+    private:
+
+      template <class G>
+      using HostGrid = decltype(std::declval<G>().hostGrid());
+
+      template <class G,
+        std::enable_if_t<not Std::is_detected<HostGrid, G>::value, int> = 0>
+      auto const& grid (G const& g) const
+      {
+        return g;
+      }
+
+      template <class G,
+        std::enable_if_t<Std::is_detected<HostGrid, G>::value, int> = 0>
+      auto const& grid (G const& g) const
+      {
+        return grid(g.hostGrid());
+      }
+
+    protected:
+      using Super::gridView_;
+      std::array<int, 6> wholeExtent_;
+      std::array<int, 6> extent_;
+      FieldVector<ctype,3> origin_;
+      FieldVector<ctype,3> spacing_;
+      int level_;
+    };
+
+    namespace Impl
+    {
+      template <class GridView, int dim, class Coordinates>
+      struct StructuredDataCollectorImpl<GridView, YaspGrid<dim,Coordinates>>
+      {
+        using type = YaspDataCollector<GridView>;
+      };
     }
-  }
-
-  void initGeometry (TensorProductCoordinates<ctype,dim> const& coords)
-  {
-    for (int i = 0; i < dim; ++i) {
-      spacing_[i] = coords.meshsize(i,0); // is not constant, but also not used.
-      origin_[i] = coords.coordinate(i,0);
-    }
-  }
-
-
-  /// Extract the ordinates from the coordinates object of the current level
-  template <class T>
-  std::array<std::vector<T>, 3> coordinatesImpl () const
-  {
-    auto it = grid(gridView_.grid()).begin(level_);
-    auto const& coords = it->coords;
-
-    std::array<std::vector<T>, 3> ordinates{};
-    for (int d = 0; d < dim; ++d) {
-      auto s = extent_[2*d+1] - extent_[2*d] + 1;
-      ordinates[d].resize(s);
-      for (int i = 0; i < s; ++i)
-        ordinates[d][i] = coords.coordinate(d, extent_[2*d] + i);
-    }
-
-    for (int d = dim; d < 3; ++d)
-      ordinates[d].resize(1, T(0));
-
-    return ordinates;
-  }
-
-
-private:
-
-  template <class G>
-  using HostGrid = decltype(std::declval<G>().hostGrid());
-
-  template <class G,
-    std::enable_if_t<not Std::is_detected<HostGrid, G>::value, int> = 0>
-  auto const& grid (G const& g) const
-  {
-    return g;
-  }
-
-  template <class G,
-    std::enable_if_t<Std::is_detected<HostGrid, G>::value, int> = 0>
-  auto const& grid (G const& g) const
-  {
-    return grid(g.hostGrid());
-  }
-
-protected:
-  using Super::gridView_;
-  std::array<int, 6> wholeExtent_;
-  std::array<int, 6> extent_;
-  FieldVector<ctype,3> origin_;
-  FieldVector<ctype,3> spacing_;
-  int level_;
-};
-
-namespace Impl
-{
-  template <class GridView, int dim, class Coordinates>
-  struct StructuredDataCollectorImpl<GridView, YaspGrid<dim,Coordinates>>
-  {
-    using type = YaspDataCollector<GridView>;
-  };
-}
 
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/defaultvtkfunction.hh b/dune/vtk/defaultvtkfunction.hh
index 32a02225e906ef974bdcd37789d4be0ed2fdb5ed..88f3afd1f849b90322e2a159535da14abd3afc31 100644
--- a/dune/vtk/defaultvtkfunction.hh
+++ b/dune/vtk/defaultvtkfunction.hh
@@ -5,76 +5,81 @@
 #include <dune/common/fmatrix.hh>
 #include <dune/common/fvector.hh>
 
-#include "vtklocalfunctioninterface.hh"
+#include "localfunctioninterface.hh"
 
 namespace Dune
 {
-  /// Type erasure for dune-functions LocalFunction interface
-  template <class GridView, class LocalFunction>
-  class LocalFunctionWrapper final
-      : public VtkLocalFunctionInterface<GridView>
+
+  namespace Vtk
   {
-    using Self = LocalFunctionWrapper;
-    using Interface = VtkLocalFunctionInterface<GridView>;
-    using Entity = typename Interface::Entity;
-    using LocalCoordinate = typename Interface::LocalCoordinate;
-
-    template <class F, class D>
-    using Range = std::decay_t<decltype(std::declval<F>()(std::declval<D>()))>;
-
-  public:
-    /// Constructor. Stores a copy of the passed `localFct` in a local variable.
-    template <class LocalFct,
-      disableCopyMove<Self, LocalFct> = 0>
-    LocalFunctionWrapper (LocalFct&& localFct)
-      : localFct_(std::forward<LocalFct>(localFct))
-    {}
-
-    /// Bind the LocalFunction to the Entity
-    virtual void bind (Entity const& entity) override
-    {
-      localFct_.bind(entity);
-    }
 
-    /// Unbind the LocalFunction from the Entity
-    virtual void unbind () override
+    /// Type erasure for dune-functions LocalFunction interface
+    template <class GridView, class LocalFunction>
+    class LocalFunctionWrapper final
+        : public LocalFunctionInterface<GridView>
     {
-      localFct_.unbind();
-    }
+      using Self = LocalFunctionWrapper;
+      using Interface = LocalFunctionInterface<GridView>;
+      using Entity = typename Interface::Entity;
+      using LocalCoordinate = typename Interface::LocalCoordinate;
 
-    /// Evaluate the LocalFunction in LocalCoordinates
-    virtual double evaluate (int comp, LocalCoordinate const& xi) const override
-    {
-      return evaluateImpl(comp, localFct_(xi));
-    }
+      template <class F, class D>
+      using Range = std::decay_t<decltype(std::declval<F>()(std::declval<D>()))>;
 
-  private:
-    // Evaluate a component of a vector valued data
-    template <class T, int N, int M>
-    double evaluateImpl (int comp, FieldMatrix<T,N,M> const& mat) const
-    {
-      int r = comp / 3;
-      int c = comp % 3;
-      return r < N && c < M ? mat[r][c] : 0.0;
-    }
-
-    // Evaluate a component of a vector valued data
-    template <class T, int N>
-    double evaluateImpl (int comp, FieldVector<T,N> const& vec) const
-    {
-      return comp < N ? vec[comp] : 0.0;
-    }
+    public:
+      /// Constructor. Stores a copy of the passed `localFct` in a local variable.
+      template <class LocalFct,
+        disableCopyMove<Self, LocalFct> = 0>
+      LocalFunctionWrapper (LocalFct&& localFct)
+        : localFct_(std::forward<LocalFct>(localFct))
+      {}
 
-    // Return the scalar values
-    template <class T>
-    double evaluateImpl (int comp, T const& value) const
-    {
-      assert(comp == 0);
-      return value;
-    }
+      /// Bind the LocalFunction to the Entity
+      virtual void bind (Entity const& entity) override
+      {
+        localFct_.bind(entity);
+      }
+
+      /// Unbind the LocalFunction from the Entity
+      virtual void unbind () override
+      {
+        localFct_.unbind();
+      }
+
+      /// Evaluate the LocalFunction in LocalCoordinates
+      virtual double evaluate (int comp, LocalCoordinate const& xi) const override
+      {
+        return evaluateImpl(comp, localFct_(xi));
+      }
+
+    private:
+      // Evaluate a component of a vector valued data
+      template <class T, int N, int M>
+      double evaluateImpl (int comp, FieldMatrix<T,N,M> const& mat) const
+      {
+        int r = comp / 3;
+        int c = comp % 3;
+        return r < N && c < M ? mat[r][c] : 0.0;
+      }
+
+      // Evaluate a component of a vector valued data
+      template <class T, int N>
+      double evaluateImpl (int comp, FieldVector<T,N> const& vec) const
+      {
+        return comp < N ? vec[comp] : 0.0;
+      }
+
+      // Return the scalar values
+      template <class T>
+      double evaluateImpl (int comp, T const& value) const
+      {
+        assert(comp == 0);
+        return value;
+      }
 
-  private:
-    LocalFunction localFct_;
-  };
+    private:
+      LocalFunction localFct_;
+    };
 
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/filereader.hh b/dune/vtk/filereader.hh
index 08fb98686ce057777a79acd09dbb6bb58ccd508d..22e0634b4243cc8e6dcc6c7f00aa0a55a2f2ac70 100644
--- a/dune/vtk/filereader.hh
+++ b/dune/vtk/filereader.hh
@@ -10,70 +10,75 @@
 
 namespace Dune
 {
-  template <class Grid, class FilerReaderImp>
-  class FileReader
+
+  namespace Vtk
   {
-  private:
-    // type of underlying implementation, for internal use only
-    using Implementation = FilerReaderImp;
 
-    /// \brief An accessor class to call protected members of reader implementations.
-    struct Accessor : public Implementation
+    template <class Grid, class FilerReaderImp>
+    class FileReader
     {
+    private:
+      // type of underlying implementation, for internal use only
+      using Implementation = FilerReaderImp;
+
+      /// \brief An accessor class to call protected members of reader implementations.
+      struct Accessor : public Implementation
+      {
+        template <class... Args>
+        static std::unique_ptr<Grid> createGridFromFileImpl (Args&&... args)
+        {
+          return Implementation::createGridFromFileImpl(std::forward<Args>(args)...);
+        }
+
+        template <class... Args>
+        static void fillFactoryImpl (Args&&... args)
+        {
+          return Implementation::fillFactoryImpl(std::forward<Args>(args)...);
+        }
+      };
+
+    public:
+      /// Reads the grid from a file with filename and returns a unique_ptr to the created grid.
+      /// Redirects to concrete implementation of derivated class.
       template <class... Args>
-      static std::unique_ptr<Grid> createGridFromFileImpl (Args&&... args)
+      static std::unique_ptr<Grid> createGridFromFile (const std::string &filename, Args&&... args)
       {
-        return Implementation::createGridFromFileImpl(std::forward<Args>(args)...);
+        return Accessor::createGridFromFileImpl(filename, std::forward<Args>(args)...);
       }
 
+      /// Reads the grid from a file with filename into a grid-factory.
+      /// Redirects to concrete implementation of derivated class.
       template <class... Args>
-      static void fillFactoryImpl (Args&&... args)
+      static void fillFactory (GridFactory<Grid> &factory,
+                              const std::string &filename,
+                              Args&&... args)
       {
-        return Implementation::fillFactoryImpl(std::forward<Args>(args)...);
+        Accessor::fillFactoryImpl(factory, filename, std::forward<Args>(args)...);
       }
-    };
-
-  public:
-    /// Reads the grid from a file with filename and returns a unique_ptr to the created grid.
-    /// Redirects to concrete implementation of derivated class.
-    template <class... Args>
-    static std::unique_ptr<Grid> createGridFromFile (const std::string &filename, Args&&... args)
-    {
-      return Accessor::createGridFromFileImpl(filename, std::forward<Args>(args)...);
-    }
-
-    /// Reads the grid from a file with filename into a grid-factory.
-    /// Redirects to concrete implementation of derivated class.
-    template <class... Args>
-    static void fillFactory (GridFactory<Grid> &factory,
-                             const std::string &filename,
-                             Args&&... args)
-    {
-      Accessor::fillFactoryImpl(factory, filename, std::forward<Args>(args)...);
-    }
 
-  protected: // default implementations
+    protected: // default implementations
 
-    // Default implementation, redirects to fillFactory implementation.
-    template <class... Args>
-    static std::unique_ptr<Grid> createGridFromFileImpl (const std::string &filename,
-                                                         Args&&... args)
-    {
-      GridFactory<Grid> factory;
-      fillFactory(factory, filename, std::forward<Args>(args)...);
+      // Default implementation, redirects to fillFactory implementation.
+      template <class... Args>
+      static std::unique_ptr<Grid> createGridFromFileImpl (const std::string &filename,
+                                                          Args&&... args)
+      {
+        GridFactory<Grid> factory;
+        fillFactory(factory, filename, std::forward<Args>(args)...);
 
-      return std::unique_ptr<Grid>{ factory.createGrid() };
-    }
+        return std::unique_ptr<Grid>{ factory.createGrid() };
+      }
 
-    // Default implementation for reading into grid-factory: produces a runtime-error.
-    template <class... Args>
-    static void fillFactoryImpl (GridFactory<Grid> &/*factory*/,
-                                 const std::string &/*filename*/,
-                                 Args&&... /*args*/)
-    {
-      DUNE_THROW(NotImplemented,
-        "GridReader using a factory argument not implemented for concrete reader implementation.");
-    }
-  };
+      // Default implementation for reading into grid-factory: produces a runtime-error.
+      template <class... Args>
+      static void fillFactoryImpl (GridFactory<Grid> &/*factory*/,
+                                  const std::string &/*filename*/,
+                                  Args&&... /*args*/)
+      {
+        DUNE_THROW(NotImplemented,
+          "GridReader using a factory argument not implemented for concrete reader implementation.");
+      }
+    };
 
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/filewriter.hh b/dune/vtk/filewriter.hh
index 803fe2566dc5a2bc2d93183e7145508497854ef0..fc90b1871da105138225794610c763031a495905 100644
--- a/dune/vtk/filewriter.hh
+++ b/dune/vtk/filewriter.hh
@@ -7,14 +7,19 @@
 
 namespace Dune
 {
-  class FileWriter
+
+  namespace Vtk
   {
-  public:
-    /// Virtual destructor
-    virtual ~FileWriter () = default;
 
-    /// Write to file given by `filename` and (optionally) store additional data in `dataDir`
-    virtual std::string write (std::string const& filename, std::optional<std::string> dataDir = {}) const = 0;
-  };
+    class FileWriter
+    {
+    public:
+      /// Virtual destructor
+      virtual ~FileWriter () = default;
+
+      /// Write to file given by `filename` and (optionally) store additional data in `dataDir`
+      virtual std::string write (std::string const& filename, std::optional<std::string> dataDir = {}) const = 0;
+    };
 
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/forward.hh b/dune/vtk/forward.hh
index f3946eab572f9d9362f38a34d9e480788c09e838..71c6c6605092687d8e469700a3ccf39d31a19428 100644
--- a/dune/vtk/forward.hh
+++ b/dune/vtk/forward.hh
@@ -4,81 +4,82 @@
 
 namespace Dune
 {
-  // forward declaration of all classes in dune-vtk
+  namespace Vtk
+  {
 
-  template <class GridView, class Derived, class Partition = Partitions::InteriorBorder>
-  class DataCollectorInterface;
+    // forward declaration of all classes in dune-vtk
 
-  // @{ datacollectors
-  template <class GridView, class Derived, class Partition = Partitions::InteriorBorder>
-  class UnstructuredDataCollectorInterface;
+    template <class GridView, class Derived, class Partition = Partitions::InteriorBorder>
+    class DataCollectorInterface;
 
-  // @{ unstructured-datacollectors
-  template <class GridView, class Partition = Partitions::InteriorBorder>
-  class ContinuousDataCollector;
+    // @{ datacollectors
+    template <class GridView, class Derived, class Partition = Partitions::InteriorBorder>
+    class UnstructuredDataCollectorInterface;
 
-  template <class GridView, class Partition = Partitions::InteriorBorder>
-  class DiscontinuousDataCollector;
+    // @{ unstructured-datacollectors
+    template <class GridView, class Partition = Partitions::InteriorBorder>
+    class ContinuousDataCollector;
 
-  template <class GridView>
-  class QuadraticDataCollector;
-  // @} unstructured-datacollectors
+    template <class GridView, class Partition = Partitions::InteriorBorder>
+    class DiscontinuousDataCollector;
 
-  template <class GridView, class Derived>
-  class StructuredDataCollectorInterface;
+    template <class GridView>
+    class QuadraticDataCollector;
+    // @} unstructured-datacollectors
 
-  namespace Impl
-  {
-    // Should be specialized for concrete structured grid
-    template <class GridView, class Grid>
-    struct StructuredDataCollectorImpl;
-  }
+    template <class GridView, class Derived>
+    class StructuredDataCollectorInterface;
 
-  template <class GridView>
-  using StructuredDataCollector = typename Impl::StructuredDataCollectorImpl<GridView, typename GridView::Grid>::type;
+    namespace Impl
+    {
+      // Should be specialized for concrete structured grid
+      template <class GridView, class Grid>
+      struct StructuredDataCollectorImpl;
+    }
 
-  // @{ structured-datacollectors
-  template <class GridView>
-  class SPDataCollector;
+    template <class GridView>
+    using StructuredDataCollector = typename Impl::StructuredDataCollectorImpl<GridView, typename GridView::Grid>::type;
 
-  template <class GridView>
-  class YaspDataCollector;
-  // @} structured-datacollectors
+    // @{ structured-datacollectors
+    template <class GridView>
+    class SPDataCollector;
 
-  // @} datacollectors
+    template <class GridView>
+    class YaspDataCollector;
+    // @} structured-datacollectors
 
-  template <class Grid, class Derived>
-  class GridCreatorInterface;
+    // @} datacollectors
 
-  template <class GridCreator, class Derived>
-  struct DerivedGridCreator;
+    template <class Grid, class Derived>
+    class GridCreatorInterface;
 
-  // @{ gridcreators
-  template <class Grid>
-  struct ContinuousGridCreator;
+    template <class GridCreator, class Derived>
+    struct DerivedGridCreator;
 
-  template <class Grid>
-  struct DiscontinuousGridCreator;
+    // @{ gridcreators
+    template <class Grid>
+    struct ContinuousGridCreator;
 
-  template <class Grid>
-  struct ParallelGridCreator;
+    template <class Grid>
+    struct DiscontinuousGridCreator;
 
-  template <class Grid>
-  struct SerialGridCreator;
+    template <class Grid>
+    struct ParallelGridCreator;
 
-  template <class Grid>
-  class LagrangeGridCreator;
-  // @} gridcreators
+    template <class Grid>
+    struct SerialGridCreator;
 
+    template <class Grid>
+    class LagrangeGridCreator;
+    // @} gridcreators
 
-  template <class Grid, class FilerReaderImp>
-  class FileReader;
 
-  template <class Grid, class GridCreator = ContinuousGridCreator<Grid>>
-  class VtkReader;
+    template <class Grid, class FilerReaderImp>
+    class FileReader;
 
+    class FileWriter;
 
-  class FileWriter;
+  } //end namespace Vtk
 
   // @{ filewriters
   template <class VtkWriter>
@@ -91,19 +92,22 @@ namespace Dune
   class VtkWriterInterface;
 
   // @{ vtkwriters
-  template <class GridView, class DataCollector = StructuredDataCollector<GridView>>
+  template <class GridView, class DataCollector = Vtk::StructuredDataCollector<GridView>>
   class VtkImageDataWriter;
 
-  template <class GridView, class DataCollector = StructuredDataCollector<GridView>>
+  template <class GridView, class DataCollector = Vtk::StructuredDataCollector<GridView>>
   class VtkRectilinearGridWriter;
 
-  template <class GridView, class DataCollector = StructuredDataCollector<GridView>>
+  template <class GridView, class DataCollector = Vtk::StructuredDataCollector<GridView>>
   class VtkStructuredGridWriter;
 
-  template <class GridView, class DataCollector = ContinuousDataCollector<GridView>>
+  template <class GridView, class DataCollector = Vtk::ContinuousDataCollector<GridView>>
   class VtkUnstructuredGridWriter;
   // @} vtkwriters
 
   // @} filewriters
 
+  template <class Grid, class GridCreator = Vtk::ContinuousGridCreator<Grid>>
+  class VtkReader;
+
 } // end namespace Dune
diff --git a/dune/vtk/function.hh b/dune/vtk/function.hh
new file mode 100644
index 0000000000000000000000000000000000000000..060c6ef17052165b4092c84e8bbafacab36ba0e1
--- /dev/null
+++ b/dune/vtk/function.hh
@@ -0,0 +1,132 @@
+#pragma once
+
+#include <optional>
+#include <type_traits>
+
+#include <dune/common/std/type_traits.hh>
+
+#include "localfunction.hh"
+#include "types.hh"
+
+namespace Dune
+{
+
+  template <class T, int N>
+  class FieldVector;
+
+  template <class T, int N, int M>
+  class FieldMatrix;
+
+  namespace Vtk
+  {
+
+    /// Wrapper class for functions allowing local evaluations.
+    template <class GridView>
+    class Function
+    {
+      template <class F>
+      using LocalFunction = decltype(localFunction(std::declval<F>()));
+
+      using Domain = typename GridView::template Codim<0>::Entity::Geometry::LocalCoordinate;
+
+      template <class F>
+      using Range = std::decay_t<std::result_of_t<F(Domain)>>;
+
+    private:
+
+      template <class T, int N>
+      static auto sizeOfImpl (FieldVector<T,N> const&)
+        -> std::integral_constant<int, N> { return {}; }
+
+      template <class T, int N, int M>
+      static auto sizeOfImpl (FieldMatrix<T,N,M> const&)
+        -> std::integral_constant<int, N*M> { return {}; }
+
+      static auto sizeOfImpl (...)
+        -> std::integral_constant<int, 1> { return {}; }
+
+      template <class T>
+      static constexpr int sizeOf () { return decltype(sizeOfImpl(std::declval<T>()))::value; }
+
+    public:
+      /// Constructor VtkFunction from legacy VTKFunction
+      /**
+      * \param fct   The VTKFunction to wrap
+      * \param type  The VTK datatype how to write the function values to the output [Vtk::FLOAT64]
+      **/
+      Function (std::shared_ptr<VTKFunction<GridView> const> const& fct,
+                  std::optional<Vtk::DataTypes> type = {})
+        : localFct_(fct)
+        , name_(fct->name())
+        , ncomps_(fct->ncomps())
+        , type_(type ? *type : Vtk::FLOAT64)
+      {}
+
+      /// Construct VtkFunction from dune-functions GridFunction with Signature
+      // NOTE: Stores the localFunction(fct) by value.
+      /**
+      * \param fct     A Grid(View)-function, providing a `localFunction(fct)`
+      * \param name    The name to use component identification in the VTK file
+      * \param ncomps  Number of components of the pointwise data. Is extracted
+      *                from the range type of the GridFunction if not given.
+      * \param type    The \ref Vtk::DataTypes used in the output. E.g. FLOAT32,
+      *                or FLOAT64. Is extracted from the range type of the
+      *                GridFunction if not given.
+      **/
+      template <class F,
+        class = void_t<LocalFunction<F>> >
+      Function (F&& fct, std::string name,
+                  std::optional<int> ncomps = {},
+                  std::optional<Vtk::DataTypes> type = {})
+        : localFct_(localFunction(std::forward<F>(fct)))
+        , name_(std::move(name))
+      {
+        using R = Range<LocalFunction<F>>;
+
+        ncomps_ = ncomps ? *ncomps : sizeOf<R>();
+        type_ = type ? *type : Vtk::Map::type<R>();
+      }
+
+      /// Constructor that forward the number of components and data type to the other constructor
+      template <class F,
+        class = void_t<LocalFunction<F>> >
+      Function (F&& fct, Vtk::FieldInfo fieldInfo,
+                  std::optional<Vtk::DataTypes> type = {})
+        : Function(std::forward<F>(fct), fieldInfo.name(), fieldInfo.ncomps(), type)
+      {}
+
+      Function () = default;
+
+      /// Create a LocalFunction
+      friend Vtk::LocalFunction<GridView> localFunction (Function const& self)
+      {
+        return self.localFct_;
+      }
+
+      /// Return a name associated with the function
+      std::string const& name () const
+      {
+        return name_;
+      }
+
+      /// Return the number of components of the Range
+      int ncomps () const
+      {
+        return ncomps_ > 3 ? 9 : ncomps_ > 1 ? 3 : 1; // tensor, vector, scalar
+      }
+
+      /// Return the VTK Datatype associated with the functions range type
+      Vtk::DataTypes type () const
+      {
+        return type_;
+      }
+
+    private:
+      Vtk::LocalFunction<GridView> localFct_;
+      std::string name_;
+      int ncomps_ = 1;
+      Vtk::DataTypes type_ = Vtk::FLOAT32;
+    };
+
+  } // end namespace Vtk
+} // end namespace Dune
diff --git a/dune/vtk/gridcreatorinterface.hh b/dune/vtk/gridcreatorinterface.hh
index 4a6b941a73c5a9dcbacc7abeb9591a6218f153af..774e818c8bf68079a2afa5acd613ecb5026272d7 100644
--- a/dune/vtk/gridcreatorinterface.hh
+++ b/dune/vtk/gridcreatorinterface.hh
@@ -13,107 +13,111 @@
 namespace Dune
 {
 
-  /// Base class for grid creators in a CRTP style.
-  /**
-   * Construct a grid from data read from VTK files.
-   *
-   * \tparam GridType         Model of Dune::Grid
-   * \tparam GlobalCoordType  Type of the global coordinates.
-   * \tparam DerivedType      Implementation of a concrete GridCreator.
-   **/
-  template <class GridType, class DerivedType>
-  class GridCreatorInterface
+  namespace Vtk
   {
-  public:
-    using Grid = GridType;
-    using GlobalCoordinate = typename Grid::template Codim<0>::Entity::Geometry::GlobalCoordinate;
-    using Derived = DerivedType;
-
-  public:
-    /// Constructor. Stores a reference to the passed GridFactory
-    GridCreatorInterface (GridFactory<Grid>& factory)
-      : factory_(&factory)
-    {}
-
-    /// Insert all points as vertices into the factory
-    void insertVertices (std::vector<GlobalCoordinate> const& points,
-                         std::vector<std::uint64_t> const& point_ids)
-    {
-      asDerived().insertVerticesImpl(points, point_ids);
-    }
-
-    /// Create elements based on type and connectivity description
-    void insertElements (std::vector<std::uint8_t> const& types,
-                         std::vector<std::int64_t> const& offsets,
-                         std::vector<std::int64_t> const& connectivity)
-    {
-      asDerived().insertElementsImpl(types, offsets, connectivity);
-    }
-
-    /// Insert part of a grid stored in file into factory
-    void insertPieces (std::vector<std::string> const& pieces)
-    {
-      asDerived().insertPiecesImpl(pieces);
-    }
-
-    /// Construct the actual grid using the GridFactory
-    std::unique_ptr<Grid> createGrid () const
-    {
-      return factory_->createGrid();
-    }
 
-    /// Return the associated GridFactory
-    GridFactory<Grid>& factory ()
+    /// Base class for grid creators in a CRTP style.
+    /**
+    * Construct a grid from data read from VTK files.
+    *
+    * \tparam GridType         Model of Dune::Grid
+    * \tparam GlobalCoordType  Type of the global coordinates.
+    * \tparam DerivedType      Implementation of a concrete GridCreator.
+    **/
+    template <class GridType, class DerivedType>
+    class GridCreatorInterface
     {
-      return *factory_;
-    }
-
-    /// Return the associated (const) GridFactory
-    GridFactory<Grid> const& factory () const
-    {
-      return *factory_;
-    }
-
-    /// Return the mpi collective communicator
-    auto comm () const
-    {
-      return MPIHelper::getCollectiveCommunication();
-    }
-
-  protected: // cast to derived type
-
-    Derived& asDerived ()
-    {
-      return static_cast<Derived&>(*this);
-    }
-
-    const Derived& asDerived () const
-    {
-      return static_cast<const Derived&>(*this);
-    }
-
-  public: // default implementations
-
-    void insertVerticesImpl (std::vector<GlobalCoordinate> const&,
-                             std::vector<std::uint64_t> const&)
-    {
-      /* do nothing */
-    }
-
-    void insertElementsImpl (std::vector<std::uint8_t> const&,
-                             std::vector<std::int64_t> const&,
-                             std::vector<std::int64_t> const&)
-    {
-      /* do nothing */
-    }
-
-    void insertPiecesImpl (std::vector<std::string> const&)
-    {
-      /* do nothing */;
-    }
-
-  protected:
-    GridFactory<Grid>* factory_;
-  };
-
+    public:
+      using Grid = GridType;
+      using GlobalCoordinate = typename Grid::template Codim<0>::Entity::Geometry::GlobalCoordinate;
+      using Derived = DerivedType;
+
+    public:
+      /// Constructor. Stores a reference to the passed GridFactory
+      GridCreatorInterface (GridFactory<Grid>& factory)
+        : factory_(&factory)
+      {}
+
+      /// Insert all points as vertices into the factory
+      void insertVertices (std::vector<GlobalCoordinate> const& points,
+                          std::vector<std::uint64_t> const& point_ids)
+      {
+        asDerived().insertVerticesImpl(points, point_ids);
+      }
+
+      /// Create elements based on type and connectivity description
+      void insertElements (std::vector<std::uint8_t> const& types,
+                          std::vector<std::int64_t> const& offsets,
+                          std::vector<std::int64_t> const& connectivity)
+      {
+        asDerived().insertElementsImpl(types, offsets, connectivity);
+      }
+
+      /// Insert part of a grid stored in file into factory
+      void insertPieces (std::vector<std::string> const& pieces)
+      {
+        asDerived().insertPiecesImpl(pieces);
+      }
+
+      /// Construct the actual grid using the GridFactory
+      std::unique_ptr<Grid> createGrid () const
+      {
+        return factory_->createGrid();
+      }
+
+      /// Return the associated GridFactory
+      GridFactory<Grid>& factory ()
+      {
+        return *factory_;
+      }
+
+      /// Return the associated (const) GridFactory
+      GridFactory<Grid> const& factory () const
+      {
+        return *factory_;
+      }
+
+      /// Return the mpi collective communicator
+      auto comm () const
+      {
+        return MPIHelper::getCollectiveCommunication();
+      }
+
+    protected: // cast to derived type
+
+      Derived& asDerived ()
+      {
+        return static_cast<Derived&>(*this);
+      }
+
+      const Derived& asDerived () const
+      {
+        return static_cast<const Derived&>(*this);
+      }
+
+    public: // default implementations
+
+      void insertVerticesImpl (std::vector<GlobalCoordinate> const&,
+                              std::vector<std::uint64_t> const&)
+      {
+        /* do nothing */
+      }
+
+      void insertElementsImpl (std::vector<std::uint8_t> const&,
+                              std::vector<std::int64_t> const&,
+                              std::vector<std::int64_t> const&)
+      {
+        /* do nothing */
+      }
+
+      void insertPiecesImpl (std::vector<std::string> const&)
+      {
+        /* do nothing */;
+      }
+
+    protected:
+      GridFactory<Grid>* factory_;
+    };
+
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/gridcreators/common.hh b/dune/vtk/gridcreators/common.hh
index e47eb4dcd50ee57d37b7bd746972ab52c694b00f..6031cec47d2c8ee006fb38791408b39007970b63 100644
--- a/dune/vtk/gridcreators/common.hh
+++ b/dune/vtk/gridcreators/common.hh
@@ -4,19 +4,23 @@
 
 namespace Dune
 {
-  template <class Factory, class... Args>
-  using HasInsertVertex = decltype( std::declval<Factory>().insertVertex(std::declval<Args>()...) );
-
-  namespace Impl
+  namespace Vtk
   {
-    template <class GF, class = void>
-    struct VertexIdType { using type = unsigned int; };
 
-    template <class GF>
-    struct VertexIdType<GF, typename GF::VertexId> { using type = typename GF::VertexId; };
-  }
+    template <class Factory, class... Args>
+    using HasInsertVertex = decltype( std::declval<Factory>().insertVertex(std::declval<Args>()...) );
+
+    namespace Impl
+    {
+      template <class GF, class = void>
+      struct VertexIdType { using type = unsigned int; };
 
-  template <class GF>
-  using VertexId_t = typename Impl::VertexIdType<GF>::type;
+      template <class GF>
+      struct VertexIdType<GF, typename GF::VertexId> { using type = typename GF::VertexId; };
+    }
+
+    template <class GF>
+    using VertexId_t = typename Impl::VertexIdType<GF>::type;
 
+  } //end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/gridcreators/continuousgridcreator.hh b/dune/vtk/gridcreators/continuousgridcreator.hh
index 3db1d0cabcf7ce53d23f2506309253b0517e7821..48715a94d3294bd476a72f626596188ac9c73441 100644
--- a/dune/vtk/gridcreators/continuousgridcreator.hh
+++ b/dune/vtk/gridcreators/continuousgridcreator.hh
@@ -9,62 +9,67 @@
 #include <dune/common/hybridutilities.hh>
 #include <dune/grid/common/gridfactory.hh>
 
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
 #include <dune/vtk/gridcreatorinterface.hh>
 
 namespace Dune
 {
-  // Create a grid where the input points and connectivity is already
-  // connected correctly.
-  template <class Grid>
-  struct ContinuousGridCreator
-      : public GridCreatorInterface<Grid, ContinuousGridCreator<Grid>>
+
+  namespace Vtk
   {
-    using Self = ContinuousGridCreator;
-    using Super = GridCreatorInterface<Grid, ContinuousGridCreator>;
-    using GlobalCoordinate = typename Super::GlobalCoordinate;
 
-    ContinuousGridCreator (GridFactory<Grid>& factory)
-      : Super(factory)
-    {}
+    // Create a grid where the input points and connectivity is already
+    // connected correctly.
+    template <class Grid>
+    struct ContinuousGridCreator
+        : public GridCreatorInterface<Grid, ContinuousGridCreator<Grid>>
+    {
+      using Self = ContinuousGridCreator;
+      using Super = GridCreatorInterface<Grid, ContinuousGridCreator>;
+      using GlobalCoordinate = typename Super::GlobalCoordinate;
 
-    using Super::factory;
+      ContinuousGridCreator (GridFactory<Grid>& factory)
+        : Super(factory)
+      {}
 
-    void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
-                             std::vector<std::uint64_t> const& /*point_ids*/)
-    {
-      for (auto const& p : points)
-        factory().insertVertex(p);
-    }
+      using Super::factory;
 
-    void insertElementsImpl (std::vector<std::uint8_t> const& types,
-                             std::vector<std::int64_t> const& offsets,
-                             std::vector<std::int64_t> const& connectivity)
-    {
-      std::size_t idx = 0;
-      for (std::size_t i = 0; i < types.size(); ++i) {
-        auto type = Vtk::to_geometry(types[i]);
-        Vtk::CellType cellType{type};
-        auto refElem = referenceElement<double,Grid::dimension>(type);
+      void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
+                              std::vector<std::uint64_t> const& /*point_ids*/)
+      {
+        for (auto const& p : points)
+          factory().insertVertex(p);
+      }
 
-        int nNodes = offsets[i] - (i == 0 ? 0 : offsets[i-1]);
-        assert(nNodes == refElem.size(Grid::dimension));
-        std::vector<unsigned int> vtk_cell; vtk_cell.reserve(nNodes);
-        for (int j = 0; j < nNodes; ++j)
-          vtk_cell.push_back( connectivity[idx++] );
+      void insertElementsImpl (std::vector<std::uint8_t> const& types,
+                              std::vector<std::int64_t> const& offsets,
+                              std::vector<std::int64_t> const& connectivity)
+      {
+        std::size_t idx = 0;
+        for (std::size_t i = 0; i < types.size(); ++i) {
+          auto type = Vtk::to_geometry(types[i]);
+          Vtk::CellType cellType{type};
+          auto refElem = referenceElement<double,Grid::dimension>(type);
 
-        if (cellType.noPermutation())
-          factory().insertElement(type,vtk_cell);
-        else {
-          // apply index permutation
-          std::vector<unsigned int> cell(nNodes);
+          int nNodes = offsets[i] - (i == 0 ? 0 : offsets[i-1]);
+          assert(nNodes == refElem.size(Grid::dimension));
+          std::vector<unsigned int> vtk_cell; vtk_cell.reserve(nNodes);
           for (int j = 0; j < nNodes; ++j)
-            cell[j] = vtk_cell[cellType.permutation(j)];
+            vtk_cell.push_back( connectivity[idx++] );
+
+          if (cellType.noPermutation())
+            factory().insertElement(type,vtk_cell);
+          else {
+            // apply index permutation
+            std::vector<unsigned int> cell(nNodes);
+            for (int j = 0; j < nNodes; ++j)
+              cell[j] = vtk_cell[cellType.permutation(j)];
 
-          factory().insertElement(type,cell);
+            factory().insertElement(type,cell);
+          }
         }
       }
-    }
-  };
+    };
 
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/gridcreators/derivedgridcreator.hh b/dune/vtk/gridcreators/derivedgridcreator.hh
index 91fa8463f028872f9160a8cce9638945615399c3..9522d4daaf2e6dc6b73095b2b5e810a08be74827 100644
--- a/dune/vtk/gridcreators/derivedgridcreator.hh
+++ b/dune/vtk/gridcreators/derivedgridcreator.hh
@@ -12,40 +12,45 @@
 
 namespace Dune
 {
-  template <class GridCreator, class Derived>
-  struct DerivedGridCreator
-      : public GridCreatorInterface<typename GridCreator::Grid, Derived>
-  {
-    using Self = DerivedGridCreator;
-    using Super = GridCreatorInterface<typename GridCreator::Grid, Derived>;
-    using Grid = typename GridCreator::Grid;
-    using GlobalCoordinate = typename Super::GlobalCoordinate;
-
-    DerivedGridCreator (GridFactory<Grid>& factory)
-      : Super(factory)
-      , gridCreator_(factory)
-    {}
-
-    void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
-                             std::vector<std::uint64_t> const& point_ids)
-    {
-      gridCreator_.insertVertices(points, point_ids);
-    }
 
-    void insertElementsImpl (std::vector<std::uint8_t> const& types,
-                             std::vector<std::int64_t> const& offsets,
-                             std::vector<std::int64_t> const& connectivity)
-    {
-      gridCreator_.insertElements(types, offsets, connectivity);
-    }
+  namespace Vtk
+  {
 
-    void insertPiecesImpl (std::vector<std::string> const& pieces)
+    template <class GridCreator, class Derived>
+    struct DerivedGridCreator
+        : public GridCreatorInterface<typename GridCreator::Grid, Derived>
     {
-      gridCreator_.insertPieces(pieces);
-    }
-
-  private:
-    GridCreator gridCreator_;
-  };
-
+      using Self = DerivedGridCreator;
+      using Super = GridCreatorInterface<typename GridCreator::Grid, Derived>;
+      using Grid = typename GridCreator::Grid;
+      using GlobalCoordinate = typename Super::GlobalCoordinate;
+
+      DerivedGridCreator (GridFactory<Grid>& factory)
+        : Super(factory)
+        , gridCreator_(factory)
+      {}
+
+      void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
+                              std::vector<std::uint64_t> const& point_ids)
+      {
+        gridCreator_.insertVertices(points, point_ids);
+      }
+
+      void insertElementsImpl (std::vector<std::uint8_t> const& types,
+                              std::vector<std::int64_t> const& offsets,
+                              std::vector<std::int64_t> const& connectivity)
+      {
+        gridCreator_.insertElements(types, offsets, connectivity);
+      }
+
+      void insertPiecesImpl (std::vector<std::string> const& pieces)
+      {
+        gridCreator_.insertPieces(pieces);
+      }
+
+    private:
+      GridCreator gridCreator_;
+    };
+
+  } // end namespace Vtk;
 } // end namespace Dune
diff --git a/dune/vtk/gridcreators/discontinuousgridcreator.hh b/dune/vtk/gridcreators/discontinuousgridcreator.hh
index c285ba51a44036d26cb8adb99372a4479e7c88f5..71dd4000b8dd676247a43494ada2bb67d639b0d6 100644
--- a/dune/vtk/gridcreators/discontinuousgridcreator.hh
+++ b/dune/vtk/gridcreators/discontinuousgridcreator.hh
@@ -9,89 +9,93 @@
 #include <dune/common/hybridutilities.hh>
 #include <dune/grid/common/gridfactory.hh>
 
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
 #include <dune/vtk/gridcreatorinterface.hh>
 namespace Dune
 {
-  // Create a grid where the input points are not connected and the connectivity
-  // describes separated elements.
-  template <class Grid>
-  struct DiscontinuousGridCreator
-      : public GridCreatorInterface<Grid, DiscontinuousGridCreator<Grid>>
+  namespace Vtk
   {
-    using Self = DiscontinuousGridCreator;
-    using Super = GridCreatorInterface<Grid, DiscontinuousGridCreator>;
-    using GlobalCoordinate = typename Super::GlobalCoordinate;
 
-    struct CoordLess
+    // Create a grid where the input points are not connected and the connectivity
+    // describes separated elements.
+    template <class Grid>
+    struct DiscontinuousGridCreator
+        : public GridCreatorInterface<Grid, DiscontinuousGridCreator<Grid>>
     {
-      template <class T, int N>
-      bool operator() (FieldVector<T,N> const& lhs, FieldVector<T,N> const& rhs) const
+      using Self = DiscontinuousGridCreator;
+      using Super = GridCreatorInterface<Grid, DiscontinuousGridCreator>;
+      using GlobalCoordinate = typename Super::GlobalCoordinate;
+
+      struct CoordLess
       {
-        for (int i = 0; i < N; ++i) {
-          if (std::abs(lhs[i] - rhs[i]) < std::numeric_limits<T>::epsilon())
-            continue;
-          return lhs[i] < rhs[i];
+        template <class T, int N>
+        bool operator() (FieldVector<T,N> const& lhs, FieldVector<T,N> const& rhs) const
+        {
+          for (int i = 0; i < N; ++i) {
+            if (std::abs(lhs[i] - rhs[i]) < std::numeric_limits<T>::epsilon())
+              continue;
+            return lhs[i] < rhs[i];
+          }
+          return false;
         }
-        return false;
-      }
-    };
+      };
 
-    DiscontinuousGridCreator (GridFactory<Grid>& factory)
-      : Super(factory)
-    {}
+      DiscontinuousGridCreator (GridFactory<Grid>& factory)
+        : Super(factory)
+      {}
 
-    using Super::factory;
-    void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
-                             std::vector<std::uint64_t> const& /*point_ids*/)
-    {
-      points_ = &points;
-      uniquePoints_.clear();
-      std::size_t idx = 0;
+      using Super::factory;
+      void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
+                              std::vector<std::uint64_t> const& /*point_ids*/)
+      {
+        points_ = &points;
+        uniquePoints_.clear();
+        std::size_t idx = 0;
 
-      for (auto const& p : points) {
-        auto b = uniquePoints_.emplace(std::make_pair(p,idx));
-        if (b.second) {
-          factory().insertVertex(p);
-          ++idx;
+        for (auto const& p : points) {
+          auto b = uniquePoints_.emplace(std::make_pair(p,idx));
+          if (b.second) {
+            factory().insertVertex(p);
+            ++idx;
+          }
         }
       }
-    }
 
-    void insertElementsImpl (std::vector<std::uint8_t> const& types,
-                             std::vector<std::int64_t> const& offsets,
-                             std::vector<std::int64_t> const& connectivity)
-    {
-      assert(points_ != nullptr);
-      std::size_t idx = 0;
-      for (std::size_t i = 0; i < types.size(); ++i) {
-        auto type = Vtk::to_geometry(types[i]);
-        Vtk::CellType cellType{type};
+      void insertElementsImpl (std::vector<std::uint8_t> const& types,
+                              std::vector<std::int64_t> const& offsets,
+                              std::vector<std::int64_t> const& connectivity)
+      {
+        assert(points_ != nullptr);
+        std::size_t idx = 0;
+        for (std::size_t i = 0; i < types.size(); ++i) {
+          auto type = Vtk::to_geometry(types[i]);
+          Vtk::CellType cellType{type};
 
-        int nNodes = offsets[i] - (i == 0 ? 0 : offsets[i-1]);
-        assert(nNodes > 0);
-        std::vector<unsigned int> vtk_cell; vtk_cell.reserve(nNodes);
-        for (int j = 0; j < nNodes; ++j) {
-          std::size_t v_j = connectivity[idx++];
-          std::size_t new_idx = uniquePoints_[(*points_)[v_j]];
-          vtk_cell.push_back(new_idx);
-        }
+          int nNodes = offsets[i] - (i == 0 ? 0 : offsets[i-1]);
+          assert(nNodes > 0);
+          std::vector<unsigned int> vtk_cell; vtk_cell.reserve(nNodes);
+          for (int j = 0; j < nNodes; ++j) {
+            std::size_t v_j = connectivity[idx++];
+            std::size_t new_idx = uniquePoints_[(*points_)[v_j]];
+            vtk_cell.push_back(new_idx);
+          }
 
-        if (cellType.noPermutation()) {
-          factory().insertElement(type,vtk_cell);
-        } else {
-          // apply index permutation
-          std::vector<unsigned int> cell(nNodes);
-          for (int j = 0; j < nNodes; ++j)
-            cell[j] = vtk_cell[cellType.permutation(j)];
-          factory().insertElement(type,cell);
+          if (cellType.noPermutation()) {
+            factory().insertElement(type,vtk_cell);
+          } else {
+            // apply index permutation
+            std::vector<unsigned int> cell(nNodes);
+            for (int j = 0; j < nNodes; ++j)
+              cell[j] = vtk_cell[cellType.permutation(j)];
+            factory().insertElement(type,cell);
+          }
         }
       }
-    }
 
-  private:
-    std::vector<GlobalCoordinate> const* points_ = nullptr;
-    std::map<GlobalCoordinate, std::size_t, CoordLess> uniquePoints_;
-  };
+    private:
+      std::vector<GlobalCoordinate> const* points_ = nullptr;
+      std::map<GlobalCoordinate, std::size_t, CoordLess> uniquePoints_;
+    };
 
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/gridcreators/lagrangegridcreator.hh b/dune/vtk/gridcreators/lagrangegridcreator.hh
index 8b6b39a2a7709538b0e0185dccfb46961fd21fde..f6d1f5c754fa79d36b0e12f5ea1c694bb06c11c3 100644
--- a/dune/vtk/gridcreators/lagrangegridcreator.hh
+++ b/dune/vtk/gridcreators/lagrangegridcreator.hh
@@ -9,377 +9,382 @@
 #include <dune/common/exceptions.hh>
 #include <dune/common/hybridutilities.hh>
 #include <dune/geometry/utility/typefromvertexcount.hh>
+#include <dune/geometry/multilineargeometry.hh>
 #include <dune/localfunctions/lagrange.hh>
 #include <dune/grid/common/gridfactory.hh>
 
 #include <dune/vtk/forward.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
 #include <dune/vtk/gridcreatorinterface.hh>
 #include <dune/vtk/utility/lagrangepoints.hh>
 
 namespace Dune
 {
-  // \brief Create a grid from data that represents higher (lagrange) cells.
-  /**
-   * The grid is created from the first nodes of a cell parametrization, representing
-   * the  corner vertices. Thus a piecewise "flat" grid is constructed. The
-   * parametrization is 1. passed as a local element parametrization to the
-   * `insertElement()` function of a gridFactory to allow the grid itself to handle the
-   * parametrization and 2. is stored internally that can be accessed by using this
-   * GridCreator object as a grid function, or by extracting locally the parametrization
-   * on each existing grid element after creation of the grid.
-   *
-   * So, the LagrangeGridCreator models both, a `GridCreator` and a `GridFunction`.
-   **/
-  template <class GridType>
-  struct LagrangeGridCreator
-      : public GridCreatorInterface<GridType, LagrangeGridCreator<GridType>>
+  namespace Vtk
   {
-    using Self = LagrangeGridCreator;
-    using Super = GridCreatorInterface<GridType, Self>;
-    using GlobalCoordinate = typename Super::GlobalCoordinate;
 
-    using Nodes = std::vector<GlobalCoordinate>;
-
-    struct ElementParametrization
+    // \brief Create a grid from data that represents higher (lagrange) cells.
+    /**
+    * The grid is created from the first nodes of a cell parametrization, representing
+    * the  corner vertices. Thus a piecewise "flat" grid is constructed. The
+    * parametrization is 1. passed as a local element parametrization to the
+    * `insertElement()` function of a gridFactory to allow the grid itself to handle the
+    * parametrization and 2. is stored internally that can be accessed by using this
+    * GridCreator object as a grid function, or by extracting locally the parametrization
+    * on each existing grid element after creation of the grid.
+    *
+    * So, the LagrangeGridCreator models both, a `GridCreator` and a `GridFunction`.
+    **/
+    template <class GridType>
+    struct LagrangeGridCreator
+        : public GridCreatorInterface<GridType, LagrangeGridCreator<GridType>>
     {
-      GeometryType type;                  //< Geometry type of the element
-      std::vector<std::int64_t> nodes;    //< Indices of the w.r.t. `nodes_` vector
-      std::vector<unsigned int> corners;  //< Insertion-indices of the element corner nodes
-    };
-
-    using Parametrization = std::vector<ElementParametrization>;
-    using Element = typename GridType::template Codim<0>::Entity;
-    using LocalCoordinate = typename Element::Geometry::LocalCoordinate;
-
-    class LocalParametrization;
-    class LocalFunction;
-
-  public:
-    using Super::factory;
-
-    LagrangeGridCreator (GridFactory<GridType>& factory)
-      : Super(factory)
-    {}
+      using Self = LagrangeGridCreator;
+      using Super = GridCreatorInterface<GridType, Self>;
+      using GlobalCoordinate = typename Super::GlobalCoordinate;
+
+      using Nodes = std::vector<GlobalCoordinate>;
+
+      struct ElementParametrization
+      {
+        GeometryType type;                  //< Geometry type of the element
+        std::vector<std::int64_t> nodes;    //< Indices of the w.r.t. `nodes_` vector
+        std::vector<unsigned int> corners;  //< Insertion-indices of the element corner nodes
+      };
+
+      using Parametrization = std::vector<ElementParametrization>;
+      using Element = typename GridType::template Codim<0>::Entity;
+      using LocalCoordinate = typename Element::Geometry::LocalCoordinate;
+
+      class LocalParametrization;
+      class LocalFunction;
+
+    public:
+      using Super::factory;
+
+      LagrangeGridCreator (GridFactory<GridType>& factory)
+        : Super(factory)
+      {}
+
+      /// Implementation of the interface function `insertVertices()`
+      void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
+                              std::vector<std::uint64_t> const& /*point_ids*/)
+      {
+        // store point coordinates in member variable
+        nodes_ = points;
+      }
 
-    /// Implementation of the interface function `insertVertices()`
-    void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
-                             std::vector<std::uint64_t> const& /*point_ids*/)
-    {
-      // store point coordinates in member variable
-      nodes_ = points;
-    }
-
-    template <class F>
-    using HasParametrizedElements = decltype(std::declval<F>().insertElement(std::declval<GeometryType>(),
-      std::declval<std::vector<unsigned int> const&>(), std::declval<std::function<GlobalCoordinate(LocalCoordinate)>>()));
-
-    /// Implementation of the interface function `insertElements()`
-    void insertElementsImpl (std::vector<std::uint8_t> const& types,
-                             std::vector<std::int64_t> const& offsets,
-                             std::vector<std::int64_t> const& connectivity)
-    {
-      assert(nodes_.size() > 0);
-
-      // mapping of node index to element-vertex index
-      std::vector<std::int64_t> elementVertices(nodes_.size(), -1);
-      parametrization_.reserve(types.size());
-
-      std::int64_t vertexIndex = 0;
-      for (std::size_t i = 0; i < types.size(); ++i) {
-        auto type = Vtk::to_geometry(types[i]);
-        if (type.dim() != GridType::dimension)
-          continue;
-
-        Vtk::CellType cellType{type};
-        auto refElem = referenceElement<double,GridType::dimension>(type);
-
-        std::int64_t shift = (i == 0 ? 0 : offsets[i-1]);
-        int nNodes = offsets[i] - shift;
-        int nVertices = refElem.size(GridType::dimension);
-
-        // insert vertices into grid and construct element vertices
-        std::vector<unsigned int> element(nVertices);
-        for (int j = 0; j < nVertices; ++j) {
-          auto index = connectivity.at(shift + j);
-          auto& vertex = elementVertices.at(index);
-          if (vertex < 0) {
-            factory().insertVertex(nodes_.at(index));
-            vertex = vertexIndex++;
+      template <class F>
+      using HasParametrizedElements = decltype(std::declval<F>().insertElement(std::declval<GeometryType>(),
+        std::declval<std::vector<unsigned int> const&>(), std::declval<std::function<GlobalCoordinate(LocalCoordinate)>>()));
+
+      /// Implementation of the interface function `insertElements()`
+      void insertElementsImpl (std::vector<std::uint8_t> const& types,
+                              std::vector<std::int64_t> const& offsets,
+                              std::vector<std::int64_t> const& connectivity)
+      {
+        assert(nodes_.size() > 0);
+
+        // mapping of node index to element-vertex index
+        std::vector<std::int64_t> elementVertices(nodes_.size(), -1);
+        parametrization_.reserve(types.size());
+
+        std::int64_t vertexIndex = 0;
+        for (std::size_t i = 0; i < types.size(); ++i) {
+          auto type = Vtk::to_geometry(types[i]);
+          if (type.dim() != GridType::dimension)
+            continue;
+
+          Vtk::CellType cellType{type};
+          auto refElem = referenceElement<double,GridType::dimension>(type);
+
+          std::int64_t shift = (i == 0 ? 0 : offsets[i-1]);
+          int nNodes = offsets[i] - shift;
+          int nVertices = refElem.size(GridType::dimension);
+
+          // insert vertices into grid and construct element vertices
+          std::vector<unsigned int> element(nVertices);
+          for (int j = 0; j < nVertices; ++j) {
+            auto index = connectivity.at(shift + j);
+            auto& vertex = elementVertices.at(index);
+            if (vertex < 0) {
+              factory().insertVertex(nodes_.at(index));
+              vertex = vertexIndex++;
+            }
+            element[j] = vertex;
           }
-          element[j] = vertex;
-        }
 
-        // permute element indices
-        if (!cellType.noPermutation()) {
-          // apply index permutation
-          std::vector<unsigned int> cell(element.size());
-          for (int j = 0; j < element.size(); ++j)
-            cell[j] = element[cellType.permutation(j)];
-          std::swap(element, cell);
-        }
+          // permute element indices
+          if (!cellType.noPermutation()) {
+            // apply index permutation
+            std::vector<unsigned int> cell(element.size());
+            for (int j = 0; j < element.size(); ++j)
+              cell[j] = element[cellType.permutation(j)];
+            std::swap(element, cell);
+          }
 
-        // fill vector of element parametrizations
-        parametrization_.push_back(ElementParametrization{type});
-        auto& param = parametrization_.back();
-
-        param.nodes.resize(nNodes);
-        for (int j = 0; j < nNodes; ++j)
-          param.nodes[j] = connectivity.at(shift + j);
-        param.corners = element;
-
-        // try to create element with parametrization
-        if constexpr (Std::is_detected_v<HasParametrizedElements, GridFactory<GridType>>) {
-          try {
-            factory().insertElement(type, element,
-              localParametrization(parametrization_.size()-1));
-          } catch (Dune::GridError const& /* notImplemented */) {
+          // fill vector of element parametrizations
+          parametrization_.push_back(ElementParametrization{type});
+          auto& param = parametrization_.back();
+
+          param.nodes.resize(nNodes);
+          for (int j = 0; j < nNodes; ++j)
+            param.nodes[j] = connectivity.at(shift + j);
+          param.corners = element;
+
+          // try to create element with parametrization
+          if constexpr (Std::is_detected_v<HasParametrizedElements, GridFactory<GridType>>) {
+            try {
+              factory().insertElement(type, element,
+                localParametrization(parametrization_.size()-1));
+            } catch (Dune::GridError const& /* notImplemented */) {
+              factory().insertElement(type, element);
+            }
+          } else {
             factory().insertElement(type, element);
           }
-        } else {
-          factory().insertElement(type, element);
         }
       }
-    }
-
-    /// \brief Construct an element parametrization
-    /**
-     * The returned LocalParametrization is a mapping `GlobalCoordinate(LocalCoordinate)`
-     * where `LocalCoordinate is w.r.t. the local coordinate system in an element with
-     * given `insertionIndex` (defined by the inserted corner vertices) and
-     * `GlobalCoordinate` a world coordinate in the parametrized grid.
-     **/
-    LocalParametrization localParametrization (unsigned int insertionIndex) const
-    {
-      assert(!nodes_.empty() && !parametrization_.empty());
-      auto const& localParam = parametrization_.at(insertionIndex);
-      return LocalParametrization{nodes_, localParam, order(localParam)};
-    }
 
-    /// \brief Construct an element parametrization
-    /**
-     * The returned LocalParametrization is a mapping `GlobalCoordinate(LocalCoordinate)`
-     * where `LocalCoordinate is w.r.t. the local coordinate system in the passed element
-     * and `GlobalCoordinate` a world coordinate in the parametrized grid.
-     *
-     * Note, when an element is passed, it might have a different local coordinate system
-     * than the coordinate system used to defined the element parametrization. Thus
-     * coordinate transform is internally chained to the evaluation of the local
-     * parametrization. This local geometry transform is obtained by figuring out the
-     * permutation of corners in the element corresponding to the inserted corner
-     * vertices.
-     **/
-    LocalParametrization localParametrization (Element const& element) const
-    {
-      assert(!nodes_.empty() && !parametrization_.empty());
-
-      unsigned int insertionIndex = factory().insertionIndex(element);
-      auto const& localParam = parametrization_.at(insertionIndex);
-      assert(element.type() == localParam.type);
-
-      // collect indices of vertices
-      std::vector<unsigned int> indices(element.subEntities(GridType::dimension));
-      for (unsigned int i = 0; i < element.subEntities(GridType::dimension); ++i)
-        indices[i] = factory().insertionIndex(element.template subEntity<GridType::dimension>(i));
-
-      // calculate permutation vector
-      std::vector<unsigned int> permutation(indices.size());
-      for (std::size_t i = 0; i < indices.size(); ++i) {
-        auto it = std::find(localParam.corners.begin(), localParam.corners.end(), indices[i]);
-        assert(it != localParam.corners.end());
-        permutation[i] = std::distance(localParam.corners.begin(), it);
+      /// \brief Construct an element parametrization
+      /**
+      * The returned LocalParametrization is a mapping `GlobalCoordinate(LocalCoordinate)`
+      * where `LocalCoordinate is w.r.t. the local coordinate system in an element with
+      * given `insertionIndex` (defined by the inserted corner vertices) and
+      * `GlobalCoordinate` a world coordinate in the parametrized grid.
+      **/
+      LocalParametrization localParametrization (unsigned int insertionIndex) const
+      {
+        assert(!nodes_.empty() && !parametrization_.empty());
+        auto const& localParam = parametrization_.at(insertionIndex);
+        return LocalParametrization{nodes_, localParam, order(localParam)};
       }
 
-      return LocalParametrization{nodes_, localParam, order(localParam), permutation};
-    }
-
-    /// Determine lagrange order from number of points
-    template <class LocalParam>
-    unsigned int order (LocalParam const& localParam) const
-    {
-      GeometryType type = localParam.type;
-      std::size_t nNodes = localParam.nodes.size();
-      for (unsigned int o = 1; o <= nNodes; ++o)
-        if (numLagrangePoints(type.id(), type.dim(), o) == nNodes)
-          return o;
-
-      return 1;
-    }
-
-    /// Determine lagrange order from number of points from the first element parametrization
-    unsigned int order () const
-    {
-      assert(!parametrization_.empty());
-      return order(parametrization_.front());
-    }
+      /// \brief Construct an element parametrization
+      /**
+      * The returned LocalParametrization is a mapping `GlobalCoordinate(LocalCoordinate)`
+      * where `LocalCoordinate is w.r.t. the local coordinate system in the passed element
+      * and `GlobalCoordinate` a world coordinate in the parametrized grid.
+      *
+      * Note, when an element is passed, it might have a different local coordinate system
+      * than the coordinate system used to defined the element parametrization. Thus
+      * coordinate transform is internally chained to the evaluation of the local
+      * parametrization. This local geometry transform is obtained by figuring out the
+      * permutation of corners in the element corresponding to the inserted corner
+      * vertices.
+      **/
+      LocalParametrization localParametrization (Element const& element) const
+      {
+        assert(!nodes_.empty() && !parametrization_.empty());
+
+        unsigned int insertionIndex = factory().insertionIndex(element);
+        auto const& localParam = parametrization_.at(insertionIndex);
+        assert(element.type() == localParam.type);
+
+        // collect indices of vertices
+        std::vector<unsigned int> indices(element.subEntities(GridType::dimension));
+        for (unsigned int i = 0; i < element.subEntities(GridType::dimension); ++i)
+          indices[i] = factory().insertionIndex(element.template subEntity<GridType::dimension>(i));
+
+        // calculate permutation vector
+        std::vector<unsigned int> permutation(indices.size());
+        for (std::size_t i = 0; i < indices.size(); ++i) {
+          auto it = std::find(localParam.corners.begin(), localParam.corners.end(), indices[i]);
+          assert(it != localParam.corners.end());
+          permutation[i] = std::distance(localParam.corners.begin(), it);
+        }
 
-  public:
-    /// \brief Local function representing the parametrization of the grid.
-    /**
-     * The returned object models Functions::Concept::LocalFunction
-     * and can thus be bound to an element of the created grid and evaluated in
-     * the local coordinates of the bound element.
-     *
-     * It is implemented in terms of the \ref LocalParametrization function
-     * returned by the method \ref localParametrization(element). See comments
-     * there for further details.
-     *
-     * Note, this methods requires the GridCreator to be based by
-     * lvalue-reference. This is necessary, since we want to guarantee that all
-     * internal storage is preserved while evaluating the local function.
-     **/
-    friend LocalFunction localFunction (LagrangeGridCreator& gridCreator)
-    {
-      return LocalFunction{gridCreator};
-    }
+        return LocalParametrization{nodes_, localParam, order(localParam), permutation};
+      }
 
-    friend LocalFunction localFunction (LagrangeGridCreator const& gridCreator)
-    {
-      return LocalFunction{gridCreator};
-    }
+      /// Determine lagrange order from number of points
+      template <class LocalParam>
+      unsigned int order (LocalParam const& localParam) const
+      {
+        GeometryType type = localParam.type;
+        std::size_t nNodes = localParam.nodes.size();
+        for (unsigned int o = 1; o <= nNodes; ++o)
+          if (numLagrangePoints(type.id(), type.dim(), o) == nNodes)
+            return o;
+
+        return 1;
+      }
 
-    friend LocalFunction localFunction (LagrangeGridCreator&& gridCreator)
-    {
-      DUNE_THROW(Dune::Exception, "Cannot pass temporary LagrangeGridCreator to localFunction(). Pass an lvalue-reference instead.");
-      return LocalFunction{gridCreator};
-    }
+      /// Determine lagrange order from number of points from the first element parametrization
+      unsigned int order () const
+      {
+        assert(!parametrization_.empty());
+        return order(parametrization_.front());
+      }
 
-    struct EntitySet
-    {
-      using Grid = GridType;
-      using GlobalCoordinate = typename Self::GlobalCoordinate;
-    };
+    public:
+      /// \brief Local function representing the parametrization of the grid.
+      /**
+      * The returned object models Functions::Concept::LocalFunction
+      * and can thus be bound to an element of the created grid and evaluated in
+      * the local coordinates of the bound element.
+      *
+      * It is implemented in terms of the \ref LocalParametrization function
+      * returned by the method \ref localParametrization(element). See comments
+      * there for further details.
+      *
+      * Note, this methods requires the GridCreator to be based by
+      * lvalue-reference. This is necessary, since we want to guarantee that all
+      * internal storage is preserved while evaluating the local function.
+      **/
+      friend LocalFunction localFunction (LagrangeGridCreator& gridCreator)
+      {
+        return LocalFunction{gridCreator};
+      }
 
-    /// Dummy function returning a placeholder entityset
-    EntitySet entitySet () const
-    {
-      assert(false && "Should not be used!");
-      return EntitySet{};
-    }
+      friend LocalFunction localFunction (LagrangeGridCreator const& gridCreator)
+      {
+        return LocalFunction{gridCreator};
+      }
 
-    /// Dummy function returning a placeholder entityset
-    GlobalCoordinate operator() (GlobalCoordinate const&) const
-    {
-      assert(false && "Should not be used!");
-      return GlobalCoordinate{};
-    }
+      friend LocalFunction localFunction (LagrangeGridCreator&& gridCreator)
+      {
+        DUNE_THROW(Dune::Exception, "Cannot pass temporary LagrangeGridCreator to localFunction(). Pass an lvalue-reference instead.");
+        return LocalFunction{gridCreator};
+      }
 
-  private:
-    /// All point coordinates inclusing the higher-order lagrange points
-    Nodes nodes_;
+      struct EntitySet
+      {
+        using Grid = GridType;
+        using GlobalCoordinate = typename Self::GlobalCoordinate;
+      };
+
+      /// Dummy function returning a placeholder entityset
+      EntitySet entitySet () const
+      {
+        assert(false && "Should not be used!");
+        return EntitySet{};
+      }
 
-    /// Parametrization for all elements
-    Parametrization parametrization_;
-  };
+      /// Dummy function returning a placeholder entityset
+      GlobalCoordinate operator() (GlobalCoordinate const&) const
+      {
+        assert(false && "Should not be used!");
+        return GlobalCoordinate{};
+      }
 
+    private:
+      /// All point coordinates inclusing the higher-order lagrange points
+      Nodes nodes_;
 
-  template <class Grid>
-  class LagrangeGridCreator<Grid>::LocalParametrization
-  {
-    using ctype = typename Grid::ctype;
-
-    using GlobalCoordinate = typename Grid::template Codim<0>::Entity::Geometry::GlobalCoordinate;
-    using LocalCoordinate = typename Grid::template Codim<0>::Entity::Geometry::LocalCoordinate;
-    using LocalGeometry = MultiLinearGeometry<ctype,Grid::dimension,Grid::dimension>;
-
-    using LocalFE = LagrangeLocalFiniteElement<VtkLagrangePointSet, Grid::dimension, ctype, ctype>;
-    using LocalBasis = typename LocalFE::Traits::LocalBasisType;
-    using LocalBasisTraits = typename LocalBasis::Traits;
-
-  public:
-    /// Construct a local element parametrization
-    template <class Nodes, class LocalParam>
-    LocalParametrization (Nodes const& nodes, LocalParam const& param, unsigned int order)
-      : localFE_(param.type, order)
-      , localNodes_(param.nodes.size())
-    {
-      for (std::size_t i = 0; i < localNodes_.size(); ++i)
-        localNodes_[i] = nodes[param.nodes[i]];
-    }
-
-    /// Construct a local element parametrization for elements with permuted corners
-    template <class Nodes, class LocalParam, class Permutation>
-    LocalParametrization (Nodes const& nodes, LocalParam const& param, unsigned int order, Permutation const& permutation)
-      : LocalParametrization(nodes, param, order)
-    {
-      auto refElem = referenceElement<ctype,Grid::dimension>(param.type);
-      std::vector<LocalCoordinate> corners(permutation.size());
-      for (std::size_t i = 0; i < permutation.size(); ++i)
-        corners[i] = refElem.position(permutation[i], Grid::dimension);
+      /// Parametrization for all elements
+      Parametrization parametrization_;
+    };
 
-      localGeometry_.emplace(param.type, corners);
-    }
 
-    /// Evaluate the local parametrization in local coordinates
-    template <class LocalCoordinate>
-    GlobalCoordinate operator() (LocalCoordinate const& local) const
+    template <class Grid>
+    class LagrangeGridCreator<Grid>::LocalParametrization
     {
-      // map coordinates if element corners are permuted
-      LocalCoordinate x = localGeometry_ ? localGeometry_->global(local) : local;
-
-      LocalBasis const& localBasis = localFE_.localBasis();
-      localBasis.evaluateFunction(x, shapeValues_);
-      assert(shapeValues_.size() == localNodes_.size());
+      using ctype = typename Grid::ctype;
+
+      using GlobalCoordinate = typename Grid::template Codim<0>::Entity::Geometry::GlobalCoordinate;
+      using LocalCoordinate = typename Grid::template Codim<0>::Entity::Geometry::LocalCoordinate;
+      using LocalGeometry = MultiLinearGeometry<ctype,Grid::dimension,Grid::dimension>;
+
+      using LocalFE = LagrangeLocalFiniteElement<Vtk::LagrangePointSet, Grid::dimension, ctype, ctype>;
+      using LocalBasis = typename LocalFE::Traits::LocalBasisType;
+      using LocalBasisTraits = typename LocalBasis::Traits;
+
+    public:
+      /// Construct a local element parametrization
+      template <class Nodes, class LocalParam>
+      LocalParametrization (Nodes const& nodes, LocalParam const& param, unsigned int order)
+        : localFE_(param.type, order)
+        , localNodes_(param.nodes.size())
+      {
+        for (std::size_t i = 0; i < localNodes_.size(); ++i)
+          localNodes_[i] = nodes[param.nodes[i]];
+      }
 
-      GlobalCoordinate out(0);
-      for (std::size_t i = 0; i < shapeValues_.size(); ++i)
-        out.axpy(shapeValues_[i], localNodes_[i]);
+      /// Construct a local element parametrization for elements with permuted corners
+      template <class Nodes, class LocalParam, class Permutation>
+      LocalParametrization (Nodes const& nodes, LocalParam const& param, unsigned int order, Permutation const& permutation)
+        : LocalParametrization(nodes, param, order)
+      {
+        auto refElem = referenceElement<ctype,Grid::dimension>(param.type);
+        std::vector<LocalCoordinate> corners(permutation.size());
+        for (std::size_t i = 0; i < permutation.size(); ++i)
+          corners[i] = refElem.position(permutation[i], Grid::dimension);
+
+        localGeometry_.emplace(param.type, corners);
+      }
 
-      return out;
-    }
+      /// Evaluate the local parametrization in local coordinates
+      template <class LocalCoordinate>
+      GlobalCoordinate operator() (LocalCoordinate const& local) const
+      {
+        // map coordinates if element corners are permuted
+        LocalCoordinate x = localGeometry_ ? localGeometry_->global(local) : local;
 
-  private:
-    LocalFE localFE_;
-    std::vector<GlobalCoordinate> localNodes_;
-    std::optional<LocalGeometry> localGeometry_;
+        LocalBasis const& localBasis = localFE_.localBasis();
+        localBasis.evaluateFunction(x, shapeValues_);
+        assert(shapeValues_.size() == localNodes_.size());
 
-    mutable std::vector<typename LocalBasisTraits::RangeType> shapeValues_;
-  };
+        GlobalCoordinate out(0);
+        for (std::size_t i = 0; i < shapeValues_.size(); ++i)
+          out.axpy(shapeValues_[i], localNodes_[i]);
 
+        return out;
+      }
 
-  template <class Grid>
-  class LagrangeGridCreator<Grid>::LocalFunction
-  {
-    using ctype = typename Grid::ctype;
-    using LocalContext = typename Grid::template Codim<0>::Entity;
-    using GlobalCoordinate = typename LocalContext::Geometry::GlobalCoordinate;
-    using LocalCoordinate = typename LocalContext::Geometry::LocalCoordinate;
-    using LocalParametrization = typename LagrangeGridCreator::LocalParametrization;
+    private:
+      LocalFE localFE_;
+      std::vector<GlobalCoordinate> localNodes_;
+      std::optional<LocalGeometry> localGeometry_;
 
-  public:
-    explicit LocalFunction (LagrangeGridCreator const& gridCreator)
-      : gridCreator_(&gridCreator)
-    {}
+      mutable std::vector<typename LocalBasisTraits::RangeType> shapeValues_;
+    };
 
-    explicit LocalFunction (LagrangeGridCreator&& gridCreator) = delete;
 
-    /// Collect a local parametrization on the element
-    void bind (LocalContext const& element)
+    template <class Grid>
+    class LagrangeGridCreator<Grid>::LocalFunction
     {
-      localContext_ = element;
-      localParametrization_.emplace(gridCreator_->localParametrization(element));
-    }
+      using ctype = typename Grid::ctype;
+      using LocalContext = typename Grid::template Codim<0>::Entity;
+      using GlobalCoordinate = typename LocalContext::Geometry::GlobalCoordinate;
+      using LocalCoordinate = typename LocalContext::Geometry::LocalCoordinate;
+      using LocalParametrization = typename LagrangeGridCreator::LocalParametrization;
+
+    public:
+      explicit LocalFunction (LagrangeGridCreator const& gridCreator)
+        : gridCreator_(&gridCreator)
+      {}
+
+      explicit LocalFunction (LagrangeGridCreator&& gridCreator) = delete;
+
+      /// Collect a local parametrization on the element
+      void bind (LocalContext const& element)
+      {
+        localContext_ = element;
+        localParametrization_.emplace(gridCreator_->localParametrization(element));
+      }
 
-    void unbind () { /* do nothing */ }
+      void unbind () { /* do nothing */ }
 
-    /// Evaluate the local parametrization in local coordinates
-    GlobalCoordinate operator() (LocalCoordinate const& local) const
-    {
-      assert(!!localParametrization_);
-      return (*localParametrization_)(local);
-    }
+      /// Evaluate the local parametrization in local coordinates
+      GlobalCoordinate operator() (LocalCoordinate const& local) const
+      {
+        assert(!!localParametrization_);
+        return (*localParametrization_)(local);
+      }
 
-    /// Return the bound element
-    LocalContext const& localContext () const
-    {
-      return localContext_;
-    }
+      /// Return the bound element
+      LocalContext const& localContext () const
+      {
+        return localContext_;
+      }
 
-  private:
-    LagrangeGridCreator const* gridCreator_;
+    private:
+      LagrangeGridCreator const* gridCreator_;
 
-    LocalContext localContext_;
-    std::optional<LocalParametrization> localParametrization_;
-  };
+      LocalContext localContext_;
+      std::optional<LocalParametrization> localParametrization_;
+    };
 
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/gridcreators/parallelgridcreator.hh b/dune/vtk/gridcreators/parallelgridcreator.hh
index 9f4632c0d98139c852074113ca28b0965faea0a1..926bb7a0979b9ec041a87f745794dc09b2321a0f 100644
--- a/dune/vtk/gridcreators/parallelgridcreator.hh
+++ b/dune/vtk/gridcreators/parallelgridcreator.hh
@@ -13,38 +13,42 @@
 
 namespace Dune
 {
-  // create a distributed grid in parallel. Currently only supported by ALUGrid
-  template <class Grid>
-  struct ParallelGridCreator
-      : public DerivedGridCreator<ContinuousGridCreator<Grid>, ParallelGridCreator<Grid>>
+  namespace Vtk
   {
-    using Self = ParallelGridCreator;
-    using Super = DerivedGridCreator<ContinuousGridCreator<Grid>, Self>;
-    using GlobalCoordinate = typename Super::GlobalCoordinate;
-    using VertexId = VertexId_t<GridFactory<Grid>>;
 
-    // The GridFactory must support insertion of global vertex IDs
-    static_assert(Std::is_detected<HasInsertVertex, GridFactory<Grid>, GlobalCoordinate, VertexId>{}, "");
-
-    ParallelGridCreator (GridFactory<Grid>& factory)
-      : Super(factory)
-    {}
-
-    void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
-                             std::vector<std::uint64_t> const& point_ids)
+    // create a distributed grid in parallel. Currently only supported by ALUGrid
+    template <class Grid>
+    struct ParallelGridCreator
+        : public DerivedGridCreator<ContinuousGridCreator<Grid>, ParallelGridCreator<Grid>>
     {
-      assert(point_ids.size() == points.size());
-      for (std::size_t i = 0; i < points.size(); ++i)
-        this->factory().insertVertex(points[i], VertexId(point_ids[i]));
-    }
+      using Self = ParallelGridCreator;
+      using Super = DerivedGridCreator<ContinuousGridCreator<Grid>, Self>;
+      using GlobalCoordinate = typename Super::GlobalCoordinate;
+      using VertexId = VertexId_t<GridFactory<Grid>>;
+
+      // The GridFactory must support insertion of global vertex IDs
+      static_assert(Std::is_detected<HasInsertVertex, GridFactory<Grid>, GlobalCoordinate, VertexId>{}, "");
+
+      ParallelGridCreator (GridFactory<Grid>& factory)
+        : Super(factory)
+      {}
+
+      void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
+                              std::vector<std::uint64_t> const& point_ids)
+      {
+        assert(point_ids.size() == points.size());
+        for (std::size_t i = 0; i < points.size(); ++i)
+          this->factory().insertVertex(points[i], VertexId(point_ids[i]));
+      }
 
-    void insertPiecesImpl (std::vector<std::string> const& pieces)
-    {
-      if (int(pieces.size()) == this->comm().size()) {
-        VtkReader<Grid, Self> pieceReader(this->factory());
-        pieceReader.read(pieces[this->comm().rank()], true);
+      void insertPiecesImpl (std::vector<std::string> const& pieces)
+      {
+        if (int(pieces.size()) == this->comm().size()) {
+          VtkReader<Grid, Self> pieceReader(this->factory());
+          pieceReader.read(pieces[this->comm().rank()], true);
+        }
       }
-    }
-  };
+    };
 
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/gridcreators/serialgridcreator.hh b/dune/vtk/gridcreators/serialgridcreator.hh
index 2d99777dd8b5214edb207bb0432e554d1c3c5048..f2d2cfe934a75c122e6e5e2efbb672e8cf3e4c3e 100644
--- a/dune/vtk/gridcreators/serialgridcreator.hh
+++ b/dune/vtk/gridcreators/serialgridcreator.hh
@@ -11,64 +11,69 @@
 
 namespace Dune
 {
-  // create a distributed grid on rank 0. Needs to be load balanced afterwards.
-  template <class Grid>
-  struct SerialGridCreator
-      : public GridCreatorInterface<Grid, SerialGridCreator<Grid>>
-  {
-    using Self = SerialGridCreator;
-    using Super = GridCreatorInterface<Grid, Self>;
-    using GlobalCoordinate = typename Super::GlobalCoordinate;
 
-    SerialGridCreator (GridFactory<Grid>& factory)
-      : Super(factory)
-    {}
+  namespace Vtk
+  {
 
-    void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
-                             std::vector<std::uint64_t> const& /*point_ids*/)
+    // create a distributed grid on rank 0. Needs to be load balanced afterwards.
+    template <class Grid>
+    struct SerialGridCreator
+        : public GridCreatorInterface<Grid, SerialGridCreator<Grid>>
     {
-      shift_.push_back(points_.size());
-      points_.reserve(points_.size() + points.size());
-      points_.insert(points_.end(), points.begin(), points.end());
-    }
+      using Self = SerialGridCreator;
+      using Super = GridCreatorInterface<Grid, Self>;
+      using GlobalCoordinate = typename Super::GlobalCoordinate;
 
-    void insertElementsImpl (std::vector<std::uint8_t> const& types,
-                             std::vector<std::int64_t> const& offsets,
-                             std::vector<std::int64_t> const& connectivity)
-    {
-      types_.reserve(types_.size() + types.size());
-      types_.insert(types_.end(), types.begin(), types.end());
+      SerialGridCreator (GridFactory<Grid>& factory)
+        : Super(factory)
+      {}
 
-      offsets_.reserve(offsets_.size() + offsets.size());
-      std::transform(offsets.begin(), offsets.end(), std::back_inserter(offsets_),
-        [shift=offsets_.empty() ? 0 : offsets_.back()](std::int64_t o) { return o + shift; });
+      void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
+                              std::vector<std::uint64_t> const& /*point_ids*/)
+      {
+        shift_.push_back(points_.size());
+        points_.reserve(points_.size() + points.size());
+        points_.insert(points_.end(), points.begin(), points.end());
+      }
 
-      connectivity_.reserve(connectivity_.size() + connectivity.size());
-      std::transform(connectivity.begin(), connectivity.end(), std::back_inserter(connectivity_),
-        [shift=shift_.back()](std::int64_t idx) { return idx + shift; });
-    }
+      void insertElementsImpl (std::vector<std::uint8_t> const& types,
+                              std::vector<std::int64_t> const& offsets,
+                              std::vector<std::int64_t> const& connectivity)
+      {
+        types_.reserve(types_.size() + types.size());
+        types_.insert(types_.end(), types.begin(), types.end());
 
-    void insertPiecesImpl (std::vector<std::string> const& pieces)
-    {
-      if (this->comm().rank() == 0) {
-        VtkReader<Grid, Self> pieceReader(*this);
-        for (std::string const& piece : pieces) {
-          pieceReader.read(piece, false);
-          pieceReader.fillGridCreator(false);
-        }
+        offsets_.reserve(offsets_.size() + offsets.size());
+        std::transform(offsets.begin(), offsets.end(), std::back_inserter(offsets_),
+          [shift=offsets_.empty() ? 0 : offsets_.back()](std::int64_t o) { return o + shift; });
+
+        connectivity_.reserve(connectivity_.size() + connectivity.size());
+        std::transform(connectivity.begin(), connectivity.end(), std::back_inserter(connectivity_),
+          [shift=shift_.back()](std::int64_t idx) { return idx + shift; });
+      }
 
-        DiscontinuousGridCreator<Grid> creator(this->factory());
-        creator.insertVertices(points_, {});
-        creator.insertElements(types_, offsets_, connectivity_);
+      void insertPiecesImpl (std::vector<std::string> const& pieces)
+      {
+        if (this->comm().rank() == 0) {
+          VtkReader<Grid, Self> pieceReader(*this);
+          for (std::string const& piece : pieces) {
+            pieceReader.read(piece, false);
+            pieceReader.fillGridCreator(false);
+          }
+
+          DiscontinuousGridCreator<Grid> creator(this->factory());
+          creator.insertVertices(points_, {});
+          creator.insertElements(types_, offsets_, connectivity_);
+        }
       }
-    }
 
-  private:
-    std::vector<GlobalCoordinate> points_;
-    std::vector<std::uint8_t> types_;
-    std::vector<std::int64_t> offsets_;
-    std::vector<std::int64_t> connectivity_;
-    std::vector<std::int64_t> shift_;
-  };
+    private:
+      std::vector<GlobalCoordinate> points_;
+      std::vector<std::uint8_t> types_;
+      std::vector<std::int64_t> offsets_;
+      std::vector<std::int64_t> connectivity_;
+      std::vector<std::int64_t> shift_;
+    };
 
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/legacyvtkfunction.hh b/dune/vtk/legacyvtkfunction.hh
index c0339e8496b11f1cce61be9dc3666ca45a6e650f..c05a1229581de28f21322653713ee1eaae042a0a 100644
--- a/dune/vtk/legacyvtkfunction.hh
+++ b/dune/vtk/legacyvtkfunction.hh
@@ -4,46 +4,51 @@
 
 #include <dune/grid/io/file/vtk/function.hh>
 
-#include "vtklocalfunctioninterface.hh"
+#include "localfunctioninterface.hh"
 
 namespace Dune
 {
-  /// Type erasure for Legacy VTKFunction
-  template <class GridView>
-  class VTKLocalFunctionWrapper final
-      : public VtkLocalFunctionInterface<GridView>
-  {
-    using Interface = VtkLocalFunctionInterface<GridView>;
-    using Entity = typename Interface::Entity;
-    using LocalCoordinate = typename Interface::LocalCoordinate;
-
-  public:
-    /// Constructor. Stores a shared pointer to the passed Dune::VTKFunction
-    VTKLocalFunctionWrapper (std::shared_ptr<VTKFunction<GridView> const> const& fct)
-      : fct_(fct)
-    {}
-
-    /// Stores a pointer to the passed entity
-    virtual void bind (Entity const& entity) override
-    {
-      entity_ = &entity;
-    }
 
-    /// Unsets the stored entity pointer
-    virtual void unbind () override
-    {
-      entity_ = nullptr;
-    }
+  namespace Vtk
+  {
 
-    /// Evaluate the Dune::VTKFunction in LocalCoordinates on the stored Entity
-    virtual double evaluate (int comp, LocalCoordinate const& xi) const override
+    /// Type erasure for Legacy VTKFunction
+    template <class GridView>
+    class VTKLocalFunctionWrapper final
+        : public LocalFunctionInterface<GridView>
     {
-      return fct_->evaluate(comp, *entity_, xi);
-    }
-
-  private:
-    std::shared_ptr<VTKFunction<GridView> const> fct_;
-    Entity const* entity_;
-  };
-
+      using Interface = LocalFunctionInterface<GridView>;
+      using Entity = typename Interface::Entity;
+      using LocalCoordinate = typename Interface::LocalCoordinate;
+
+    public:
+      /// Constructor. Stores a shared pointer to the passed Dune::VTKFunction
+      VTKLocalFunctionWrapper (std::shared_ptr<VTKFunction<GridView> const> const& fct)
+        : fct_(fct)
+      {}
+
+      /// Stores a pointer to the passed entity
+      virtual void bind (Entity const& entity) override
+      {
+        entity_ = &entity;
+      }
+
+      /// Unsets the stored entity pointer
+      virtual void unbind () override
+      {
+        entity_ = nullptr;
+      }
+
+      /// Evaluate the Dune::VTKFunction in LocalCoordinates on the stored Entity
+      virtual double evaluate (int comp, LocalCoordinate const& xi) const override
+      {
+        return fct_->evaluate(comp, *entity_, xi);
+      }
+
+    private:
+      std::shared_ptr<VTKFunction<GridView> const> fct_;
+      Entity const* entity_;
+    };
+
+  } //end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/localfunction.hh b/dune/vtk/localfunction.hh
new file mode 100644
index 0000000000000000000000000000000000000000..022fa9b040bbb6f08517596b5c4276982e4f0088
--- /dev/null
+++ b/dune/vtk/localfunction.hh
@@ -0,0 +1,76 @@
+#pragma once
+
+#include <memory>
+#include <type_traits>
+
+#include <dune/common/std/type_traits.hh>
+
+#include "localfunctioninterface.hh"
+#include "legacyvtkfunction.hh"
+#include "defaultvtkfunction.hh"
+
+namespace Dune
+{
+
+  namespace Vtk
+  {
+
+    /// \brief A Vtk::LocalFunction is a function-like object that can be bound to a grid element
+    /// an that provides an evaluate method with a component argument.
+    /**
+    * Stores internally a Vtk::LocalFunctionInterface object for the concrete evaluation.
+    **/
+    template <class GridView>
+    class LocalFunction
+    {
+      using Self = LocalFunction;
+      using Entity = typename GridView::template Codim<0>::Entity;
+      using LocalCoordinate = typename Entity::Geometry::LocalCoordinate;
+
+      template <class LF, class E>
+      using HasBind = decltype(std::declval<LF>().bind(std::declval<E>()));
+
+    public:
+      /// Construct the Vtk::LocalFunction from any function object that has a bind(element) method.
+      template <class LF,
+        disableCopyMove<Self, LF> = 0,
+        class = void_t<HasBind<LF,Entity>> >
+      LocalFunction (LF&& lf)
+        : localFct_(std::make_shared<LocalFunctionWrapper<GridView,LF>>(std::forward<LF>(lf)))
+      {}
+
+      /// Construct a Vtk::LocalFunction from a legacy VTKFunction
+      LocalFunction (std::shared_ptr<VTKFunction<GridView> const> const& lf)
+        : localFct_(std::make_shared<VTKLocalFunctionWrapper<GridView>>(lf))
+      {}
+
+      /// Allow the default construction of a Vtk::LocalFunction
+      LocalFunction () = default;
+
+      /// Bind the function to the grid entity
+      void bind (Entity const& entity)
+      {
+        assert(bool(localFct_));
+        localFct_->bind(entity);
+      }
+
+      /// Unbind from the currently bound entity
+      void unbind ()
+      {
+        assert(bool(localFct_));
+        localFct_->unbind();
+      }
+
+      /// Evaluate the `comp` component of the Range value at local coordinate `xi`
+      double evaluate (int comp, LocalCoordinate const& xi) const
+      {
+        assert(bool(localFct_));
+        return localFct_->evaluate(comp, xi);
+      }
+
+    private:
+      std::shared_ptr<LocalFunctionInterface<GridView>> localFct_ = nullptr;
+    };
+
+  } // end namespace Vtk
+} // end namespace Dune
diff --git a/dune/vtk/localfunctioninterface.hh b/dune/vtk/localfunctioninterface.hh
new file mode 100644
index 0000000000000000000000000000000000000000..900b069fb052bbb25371cbac6256aff6ebf8e2a7
--- /dev/null
+++ b/dune/vtk/localfunctioninterface.hh
@@ -0,0 +1,32 @@
+#pragma once
+
+namespace Dune
+{
+
+  namespace Vtk
+  {
+
+    /// \brief An abstract base class for LocalFunctions that can be bound to an element and
+    /// evaluated in local coordinates w.r.t. to a component of its value.
+    template <class GridView>
+    class LocalFunctionInterface
+    {
+    public:
+      using Entity = typename GridView::template Codim<0>::Entity;
+      using LocalCoordinate = typename Entity::Geometry::LocalCoordinate;
+
+      /// Bind the function to the grid entity
+      virtual void bind (Entity const& entity) = 0;
+
+      /// Unbind from the currently bound entity
+      virtual void unbind () = 0;
+
+      /// Evaluate single component comp in the entity at local coordinates xi
+      virtual double evaluate (int comp, LocalCoordinate const& xi) const = 0;
+
+      /// Virtual destructor
+      virtual ~LocalFunctionInterface () = default;
+    };
+
+  } // end namespace Vtk
+} // end namespace Dune
diff --git a/dune/vtk/pvdwriter.hh b/dune/vtk/pvdwriter.hh
index 7b256cc802e8f292219e46177a555a1f11ffa1cf..a3d78ecb7250208247be225154fb9982112a86d0 100644
--- a/dune/vtk/pvdwriter.hh
+++ b/dune/vtk/pvdwriter.hh
@@ -6,7 +6,7 @@
 #include <vector>
 #include <tuple>
 
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
 #include <dune/vtk/filewriter.hh>
 #include <dune/vtk/forward.hh>
 
@@ -15,7 +15,7 @@ namespace Dune
   /// File-Writer for ParaView .pvd files
   template <class VtkWriter>
   class PvdWriter
-      : public FileWriter
+      : public Vtk::FileWriter
   {
     using Self = PvdWriter;
 
@@ -53,7 +53,7 @@ namespace Dune
      **/
     virtual std::string write (std::string const& fn, std::optional<std::string> dir = {}) const override;
 
-    /// Attach point data to the writer, \see VtkFunction for possible arguments
+    /// Attach point data to the writer, \see Vtk::Function for possible arguments
     template <class Function, class... Args>
     PvdWriter& addPointData (Function const& fct, Args&&... args)
     {
@@ -61,7 +61,7 @@ namespace Dune
       return *this;
     }
 
-    /// Attach cell data to the writer, \see VtkFunction for possible arguments
+    /// Attach cell data to the writer, \see Vtk::Function for possible arguments
     template <class Function, class... Args>
     PvdWriter& addCellData (Function const& fct, Args&&... args)
     {
diff --git a/dune/vtk/pvdwriter.impl.hh b/dune/vtk/pvdwriter.impl.hh
index 7f84bd1ffdbd945fb99941010355c4eb3f410084..e7b7a563540a85f18ff7956ba16175ec791d70a5 100644
--- a/dune/vtk/pvdwriter.impl.hh
+++ b/dune/vtk/pvdwriter.impl.hh
@@ -12,13 +12,13 @@ template <class W>
 void PvdWriter<W>
   ::writeTimestep (double time, std::string const& fn, std::optional<std::string> dir, bool writeCollection) const
 {
-  auto p = filesystem::path(fn);
+  auto p = Vtk::path(fn);
   auto name = p.stem();
   p.remove_filename();
 
-  filesystem::path fn_dir = p;
-  filesystem::path data_dir = dir ? filesystem::path(*dir) : fn_dir;
-  filesystem::path rel_dir = filesystem::relative(data_dir, fn_dir);
+  Vtk::path fn_dir = p;
+  Vtk::path data_dir = dir ? Vtk::path(*dir) : fn_dir;
+  Vtk::path rel_dir = Vtk::relative(data_dir, fn_dir);
 
   std::string pvd_fn = fn_dir.string() + '/' + name.string();
   std::string seq_fn = data_dir.string() + '/' + name.string() + "_t" + std::to_string(timesteps_.size());
@@ -52,7 +52,7 @@ template <class W>
 std::string PvdWriter<W>
   ::write (std::string const& fn, std::optional<std::string> /*dir*/) const
 {
-  auto p = filesystem::path(fn);
+  auto p = Vtk::path(fn);
   auto name = p.stem();
   p.remove_filename();
   p /= name.string();
diff --git a/dune/vtk/vtkreader.hh b/dune/vtk/reader.hh
similarity index 97%
rename from dune/vtk/vtkreader.hh
rename to dune/vtk/reader.hh
index ccb06e5564ab3fb73cf6f1c6a5318f326faafbea..417d125a7bf647a5f23075c6b424d7e1673ba4a6 100644
--- a/dune/vtk/vtkreader.hh
+++ b/dune/vtk/reader.hh
@@ -10,7 +10,7 @@
 
 #include <dune/vtk/filereader.hh>
 #include <dune/vtk/forward.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
 
 // default GridCreator
 #include <dune/vtk/gridcreators/continuousgridcreator.hh>
@@ -26,7 +26,7 @@ namespace Dune
    **/
   template <class Grid, class GridCreator>
   class VtkReader
-      : public FileReader<Grid, VtkReader<Grid, GridCreator>>
+      : public Vtk::FileReader<Grid, VtkReader<Grid, GridCreator>>
   {
     // Sections visited during the xml parsing
     enum Sections {
@@ -203,7 +203,7 @@ namespace Dune
   // deduction guides
   template <class Grid>
   VtkReader (GridFactory<Grid>&)
-    -> VtkReader<Grid, ContinuousGridCreator<Grid>>;
+    -> VtkReader<Grid, Vtk::ContinuousGridCreator<Grid>>;
 
   template <class GridCreator,
     class = std::void_t<typename GridCreator::Grid>>
@@ -217,4 +217,4 @@ namespace Dune
 
 } // end namespace Dune
 
-#include "vtkreader.impl.hh"
+#include "reader.impl.hh"
diff --git a/dune/vtk/vtkreader.impl.hh b/dune/vtk/reader.impl.hh
similarity index 98%
rename from dune/vtk/vtkreader.impl.hh
rename to dune/vtk/reader.impl.hh
index dff794ae03dbd595d5efa75522673ace153a84f9..542f8371976b89338bdbc57d52aaf1d1b53766a3 100644
--- a/dune/vtk/vtkreader.impl.hh
+++ b/dune/vtk/reader.impl.hh
@@ -19,13 +19,13 @@ template <class Grid, class Creator>
 void VtkReader<Grid,Creator>::read (std::string const& filename, bool fillCreator)
 {
   // check whether file exists!
-  if (!filesystem::exists(filename))
+  if (!Vtk::exists(filename))
     DUNE_THROW(IOError, "File " << filename << " does not exist!");
 
   std::ifstream input(filename, std::ios_base::in | std::ios_base::binary);
   assert(input.is_open());
 
-  std::string ext = filesystem::path(filename).extension().string();
+  std::string ext = Vtk::path(filename).extension().string();
   if (ext == ".vtu") {
     readSerialFileFromStream(input, fillCreator);
     pieces_.push_back(filename);
@@ -49,7 +49,7 @@ void VtkReader<Grid,Creator>::readSerialFileFromStream (std::ifstream& input, bo
 
   Sections section = NO_SECTION;
   for (std::string line; std::getline(input, line); ) {
-    ltrim(line);
+    Vtk::ltrim(line);
 
     if (isSection(line, "VTKFile", section)) {
       bool closed = false;
@@ -109,7 +109,7 @@ void VtkReader<Grid,Creator>::readSerialFileFromStream (std::ifstream& input, bo
       data_type = Vtk::Map::to_datatype[attr["type"]];
 
       if (!attr["Name"].empty())
-        data_name = to_lower(attr["Name"]);
+        data_name = Vtk::to_lower(attr["Name"]);
       else if (section == POINTS)
         data_name = "points";
 
@@ -118,7 +118,7 @@ void VtkReader<Grid,Creator>::readSerialFileFromStream (std::ifstream& input, bo
         data_components = std::stoul(attr["NumberOfComponents"]);
 
       // determine FormatType
-      data_format = to_lower(attr["format"]);
+      data_format = Vtk::to_lower(attr["format"]);
       if (data_format == "appended") {
         format_ = !compressor.empty() ? Vtk::COMPRESSED : Vtk::BINARY;
       } else {
@@ -139,7 +139,7 @@ void VtkReader<Grid,Creator>::readSerialFileFromStream (std::ifstream& input, bo
       if (data_format == "appended") {
         if (!closed) {
           while (std::getline(input, line)) {
-            ltrim(line);
+            Vtk::ltrim(line);
             if (line.substr(1,10) == "/DataArray")
               break;
           }
@@ -237,7 +237,7 @@ void VtkReader<Grid,Creator>::readParallelFileFromStream (std::ifstream& input,
 
   Sections section = NO_SECTION;
   for (std::string line; std::getline(input, line); ) {
-    ltrim(line);
+    Vtk::ltrim(line);
 
     if (isSection(line, "VTKFile", section)) {
       bool closed = false;
@@ -297,7 +297,7 @@ Sections readDataArray (IStream& input, std::vector<T>& values, std::size_t max_
 
   std::size_t idx = 0;
   for (std::string line; std::getline(input, line);) {
-    trim(line);
+    Vtk::trim(line);
     if (line.substr(1,10) == "/DataArray")
       return parent_section;
     if (line[0] == '<')
@@ -318,7 +318,7 @@ template <class IStream, class Sections>
 Sections skipRestOfDataArray (IStream& input, Sections section, Sections parent_section)
 {
   for (std::string line; std::getline(input, line);) {
-    ltrim(line);
+    Vtk::ltrim(line);
     if (line.substr(1,10) == "/DataArray")
       return parent_section;
   }
diff --git a/dune/vtk/test/test-typededuction.cc b/dune/vtk/test/test-typededuction.cc
index 0cc21c23430ac1bcdce0bf343f38209da0324ffe..164cfdbd178c9310fecfd78bb398e5fbc9817dc2 100644
--- a/dune/vtk/test/test-typededuction.cc
+++ b/dune/vtk/test/test-typededuction.cc
@@ -5,8 +5,8 @@
 #include "config.h" // autoconf defines, needed by the dune headers
 #endif
 
-#include <dune/vtk/vtkreader.hh>
-#include <dune/vtk/vtkwriter.hh>
+#include <dune/vtk/reader.hh>
+#include <dune/vtk/writer.hh>
 
 #if HAVE_UG
   #include <dune/grid/uggrid.hh>
@@ -28,7 +28,7 @@ int main (int argc, char** argv)
   VtkUnstructuredGridWriter writer1(grid->leafGridView());
 
   // 2. construct writer from datacollector
-  ContinuousDataCollector dataCollector1(grid->leafGridView());
+  Vtk::ContinuousDataCollector dataCollector1(grid->leafGridView());
   VtkUnstructuredGridWriter writer2(dataCollector1);
   VtkUnstructuredGridWriter writer3(stackobject_to_shared_ptr(dataCollector1));
 
@@ -40,6 +40,6 @@ int main (int argc, char** argv)
   VtkReader reader1(factory);
 
   // 5. construct reader from grid-creator
-  ContinuousGridCreator creator(factory);
+  Vtk::ContinuousGridCreator creator(factory);
   VtkReader reader2(creator);
 }
\ No newline at end of file
diff --git a/dune/vtk/test/test-vtkwriter.cc b/dune/vtk/test/test-vtkwriter.cc
index 721b61366a41c206c976dfaac266461348a89f9e..8cc2200012f6f3c3309f120be8642cf84c4ca2f3 100644
--- a/dune/vtk/test/test-vtkwriter.cc
+++ b/dune/vtk/test/test-vtkwriter.cc
@@ -20,7 +20,7 @@
   #include <dune/grid/uggrid.hh>
   #include <dune/grid/utility/structuredgridfactory.hh>
 #endif
-#include <dune/vtk/vtkwriter.hh>
+#include <dune/vtk/writer.hh>
 
 #include "checkvtkfile.hh"
 
diff --git a/dune/vtk/vtktimeserieswriter.hh b/dune/vtk/timeserieswriter.hh
similarity index 92%
rename from dune/vtk/vtktimeserieswriter.hh
rename to dune/vtk/timeserieswriter.hh
index bbac2b301ff0b4f167c0b7814b863ddee83f3a9e..ba500528b2cd15922825d297ee58062242a0b83d 100644
--- a/dune/vtk/vtktimeserieswriter.hh
+++ b/dune/vtk/timeserieswriter.hh
@@ -8,7 +8,7 @@
 
 #include <dune/vtk/filewriter.hh>
 #include <dune/vtk/forward.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
 #include <dune/vtk/utility/filesystem.hh>
 #include <dune/vtk/utility/uid.hh>
 
@@ -22,7 +22,7 @@ namespace Dune
    **/
   template <class VtkWriter>
   class VtkTimeseriesWriter
-      : public FileWriter
+      : public Vtk::FileWriter
   {
   protected:
     using Self = VtkTimeseriesWriter;
@@ -47,7 +47,7 @@ namespace Dune
       std::srand(std::time(nullptr));
 
       // put temporary file to /tmp directory
-      tmpDir_ = filesystem::path("/tmp/vtktimeserieswriter_" + uid(10) + "/");
+      tmpDir_ = Vtk::path("/tmp/vtktimeserieswriter_" + Vtk::uid(10) + "/");
     }
 
     /// Remove all written intermediate files and remove temporary directory
@@ -75,7 +75,7 @@ namespace Dune
      *
      * \param fn   Filename of the Timeseries file. May contain a directory and any file extension.
      * \param dir  The optional parameter specifies the directory of the partition files.
-     * 
+     *
      * \returns File name of the actual written file
      **/
     virtual std::string write (std::string const& fn, std::optional<std::string> dir = {}) const override;
@@ -88,7 +88,7 @@ namespace Dune
       return *this;
     }
 
-    /// Attach cell data to the writer, \see VtkFunction for possible arguments
+    /// Attach cell data to the writer, \see Vtk::Function for possible arguments
     template <class Function, class... Args>
     VtkTimeseriesWriter& addCellData (Function const& fct, Args&&... args)
     {
@@ -98,7 +98,7 @@ namespace Dune
 
   protected:
     VtkWriter vtkWriter_;
-    filesystem::path tmpDir_;
+    Vtk::path tmpDir_;
 
     mutable bool initialized_ = false;
 
@@ -111,4 +111,4 @@ namespace Dune
 
 } // end namespace Dune
 
-#include "vtktimeserieswriter.impl.hh"
+#include "timeserieswriter.impl.hh"
diff --git a/dune/vtk/vtktimeserieswriter.impl.hh b/dune/vtk/timeserieswriter.impl.hh
similarity index 90%
rename from dune/vtk/vtktimeserieswriter.impl.hh
rename to dune/vtk/timeserieswriter.impl.hh
index f3822daa7fc832cbe157a3e91fe71b06a55146fe..3a3ef66bc67c8ecaf5663c6a28aaccb482ed6a47 100644
--- a/dune/vtk/vtktimeserieswriter.impl.hh
+++ b/dune/vtk/timeserieswriter.impl.hh
@@ -41,8 +41,8 @@ template <class W>
 void VtkTimeseriesWriter<W>
   ::writeTimestep (double time, std::string const& fn, std::optional<std::string> tmpDir, bool writeCollection) const
 {
-  auto name = filesystem::path(fn).stem();
-  auto tmpBase = tmpDir ? filesystem::path(*tmpDir) : tmpDir_;
+  auto name = Vtk::path(fn).stem();
+  auto tmpBase = tmpDir ? Vtk::path(*tmpDir) : tmpDir_;
   auto tmp = tmpBase;
   tmp /= name.string();
 
@@ -54,7 +54,7 @@ void VtkTimeseriesWriter<W>
     filenameBase = tmp.string() + "_p" + std::to_string(vtkWriter_.comm().rank());
 
   if (!initialized_) {
-    filesystem::create_directories(tmpBase);
+    Vtk::create_directories(tmpBase);
 
     // write points and cells only once
     filenameMesh_ = filenameBase + ".mesh.vtkdata";
@@ -79,13 +79,13 @@ std::string VtkTimeseriesWriter<W>
 {
   assert( initialized_ );
 
-  auto p = filesystem::path(fn);
+  auto p = Vtk::path(fn);
   auto name = p.stem();
   p.remove_filename();
 
-  filesystem::path fn_dir = p;
-  filesystem::path data_dir = dir ? filesystem::path(*dir) : fn_dir;
-  filesystem::path rel_dir = filesystem::relative(data_dir, fn_dir);
+  Vtk::path fn_dir = p;
+  Vtk::path data_dir = dir ? Vtk::path(*dir) : fn_dir;
+  Vtk::path rel_dir = Vtk::relative(data_dir, fn_dir);
 
   std::string serial_fn = fn_dir.string() + '/' + name.string() + "_ts";
   std::string parallel_fn = data_dir.string() + '/' + name.string() + "_ts";
diff --git a/dune/vtk/types.cc b/dune/vtk/types.cc
new file mode 100644
index 0000000000000000000000000000000000000000..866f3410fd1c367c72f5572c468ca30ceab55512
--- /dev/null
+++ b/dune/vtk/types.cc
@@ -0,0 +1,204 @@
+#include "types.hh"
+
+#include <iostream>
+
+#include <dune/common/exceptions.hh>
+
+namespace Dune
+{
+
+  namespace Vtk
+  {
+
+    std::string to_string (FormatTypes type)
+    {
+      switch (type) {
+        case ASCII:      return "ascii";
+        case BINARY:     return "binary";
+        case COMPRESSED: return "compressed";
+        case APPENDED:   return "appended";
+        default:
+          DUNE_THROW(RangeError, "FormatType not found.");
+          std::abort();
+      }
+    }
+
+    std::string to_string (DataTypes type)
+    {
+      switch (type) {
+        case INT8:    return "Int8";
+        case UINT8:   return "UInt8";
+        case INT16:   return "Int16";
+        case UINT16:  return "UInt16";
+        case INT32:   return "Int32";
+        case UINT32:  return "UInt32";
+        case INT64:   return "Int64";
+        case UINT64:  return "UInt64";
+        case FLOAT32: return "Float32";
+        case FLOAT64: return "Float64";
+        default:
+          DUNE_THROW(RangeError, "DataType not found.");
+          std::abort();
+      }
+    }
+
+    GeometryType to_geometry (std::uint8_t cell)
+    {
+      switch (cell) {
+        case VERTEX:     return GeometryTypes::vertex;
+        case LINE:       return GeometryTypes::line;
+        case TRIANGLE:   return GeometryTypes::triangle;
+        case QUAD:       return GeometryTypes::quadrilateral;
+        case TETRA:      return GeometryTypes::tetrahedron;
+        case HEXAHEDRON: return GeometryTypes::hexahedron;
+        case WEDGE:      return GeometryTypes::prism;
+        case PYRAMID:    return GeometryTypes::pyramid;
+
+        // Quadratic VTK cell types
+        case QUADRATIC_EDGE:        return GeometryTypes::line;
+        case QUADRATIC_TRIANGLE:    return GeometryTypes::triangle;
+        case QUADRATIC_QUAD:        return GeometryTypes::quadrilateral;
+        case QUADRATIC_TETRA:       return GeometryTypes::tetrahedron;
+        case QUADRATIC_HEXAHEDRON:  return GeometryTypes::hexahedron;
+
+        // Arbitrary order Lagrange elements
+        case LAGRANGE_CURVE:        return GeometryTypes::line;
+        case LAGRANGE_TRIANGLE:     return GeometryTypes::triangle;
+        case LAGRANGE_QUADRILATERAL:return GeometryTypes::quadrilateral;
+        case LAGRANGE_TETRAHEDRON:  return GeometryTypes::tetrahedron;
+        case LAGRANGE_HEXAHEDRON:   return GeometryTypes::hexahedron;
+        case LAGRANGE_WEDGE:        return GeometryTypes::prism;
+        default:
+          DUNE_THROW(RangeError, "CellType does not map to GeometryType.");
+          std::abort();
+      }
+    }
+
+
+    std::map<std::string, DataTypes> Map::to_datatype = {
+      {"Int8",    INT8},
+      {"UInt8",   UINT8},
+      {"Int16",   INT16},
+      {"UInt16",  UINT16},
+      {"Int32",   INT32},
+      {"UInt32",  UINT32},
+      {"Int64",   INT64},
+      {"UInt64",  UINT64},
+      {"Float32", FLOAT32},
+      {"Float64", FLOAT64}
+    };
+
+
+
+    CellType::CellType (GeometryType const& t, CellParametrization parametrization)
+      : noPermutation_(true)
+    {
+      if (parametrization == LINEAR) {
+        if (t.isVertex()) {
+          type_ = VERTEX;
+          permutation_ = {0};
+        }
+        else if (t.isLine()) {
+          type_ = LINE;
+          permutation_ = {0,1};
+        }
+        else if (t.isTriangle()) {
+          type_ = TRIANGLE;
+          permutation_ = {0,1,2};
+        }
+        else if (t.isQuadrilateral()) {
+          type_ = QUAD;
+          permutation_ = {0,1,3,2};
+          noPermutation_ = false;
+        }
+        else if (t.isTetrahedron()) {
+          type_ = TETRA;
+          permutation_ = {0,1,2,3};
+        }
+        else if (t.isHexahedron()) {
+          type_ = HEXAHEDRON;
+          permutation_ = {0,1,3,2,4,5,7,6};
+          noPermutation_ = false;
+        }
+        else if (t.isPrism()) {
+          type_ = WEDGE;
+          permutation_ = {0,2,1,3,5,4};
+          noPermutation_ = false;
+        }
+        else if (t.isPyramid()) {
+          type_ = PYRAMID;
+          permutation_ = {0,1,3,2,4};
+          noPermutation_ = false;
+        }
+        else if (t.isNone() && t.dim() == 1) {
+          type_ = LINE;
+          permutation_ = {0,1};
+        }
+        else if (t.isNone() && t.dim() == 2) {
+          type_ = POLYGON;
+          permutation_ = {0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19};
+        }
+        else {
+          std::cerr << "Geometry Type not supported by VTK!\n";
+          std::abort();
+        }
+      } else if (parametrization == QUADRATIC) {
+        if (t.isLine()) {
+          type_ = QUADRATIC_EDGE;
+          permutation_ = {0,1, 0};
+        }
+        else if (t.isTriangle()) {
+          type_ = QUADRATIC_TRIANGLE;
+          permutation_ = {0,1,2, 0,2,1};
+          noPermutation_ = false;
+        }
+        else if (t.isQuadrilateral()) {
+          type_ = QUADRATIC_QUAD;
+          permutation_ = {0,1,3,2, 2,1,3,0};
+          noPermutation_ = false;
+        }
+        else if (t.isTetrahedron()) {
+          type_ = QUADRATIC_TETRA;
+          permutation_ = {0,1,2,3, 0,2,1,3,4,5};
+          noPermutation_ = false;
+        }
+        else if (t.isHexahedron()) {
+          type_ = QUADRATIC_HEXAHEDRON;
+          permutation_ = {0,1,3,2,4,5,7,6, 6,5,7,4,10,9,11,8,0,1,3,2};
+          noPermutation_ = false;
+        }
+        else {
+          std::cerr << "Geometry Type not supported by VTK!\n";
+          std::abort();
+        }
+      } else if (parametrization == LAGRANGE) {
+        if (t.isLine()) {
+          type_ = LAGRANGE_CURVE;
+        }
+        else if (t.isTriangle()) {
+          type_ = LAGRANGE_TRIANGLE;
+        }
+        else if (t.isQuadrilateral()) {
+          type_ = LAGRANGE_QUADRILATERAL;
+        }
+        else if (t.isTetrahedron()) {
+          type_ = LAGRANGE_TETRAHEDRON;
+        }
+        else if (t.isHexahedron()) {
+          type_ = LAGRANGE_HEXAHEDRON;
+        }
+        else if (t.isPrism()) {
+          type_ = LAGRANGE_WEDGE;
+        }
+        else if (t.isPyramid()) {
+          type_ = LAGRANGE_PYRAMID;
+        }
+        else {
+          std::cerr << "Geometry Type not supported by VTK!\n";
+          std::abort();
+        }
+      }
+    }
+
+  } // end namespace Vtk
+} // end namespace Dune
diff --git a/dune/vtk/vtktypes.hh b/dune/vtk/types.hh
similarity index 100%
rename from dune/vtk/vtktypes.hh
rename to dune/vtk/types.hh
diff --git a/dune/vtk/utility/enum.hh b/dune/vtk/utility/enum.hh
index 6c0b3fc16f71845d9215c2b25069ea336921b00d..025719e97b22448052c3de059483d4a474000d1c 100644
--- a/dune/vtk/utility/enum.hh
+++ b/dune/vtk/utility/enum.hh
@@ -4,11 +4,15 @@
 
 namespace Dune
 {
-  template <class E, class Integer,
-    std::enable_if_t<std::is_enum<E>::value, int> = 0>
-  constexpr bool is_a(E a, Integer b)
+  namespace Vtk
   {
-    return (int(a) & int(b)) != 0;
-  }
 
+    template <class E, class Integer,
+      std::enable_if_t<std::is_enum<E>::value, int> = 0>
+    constexpr bool is_a(E a, Integer b)
+    {
+      return (int(a) & int(b)) != 0;
+    }
+
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/utility/filesystem.cc b/dune/vtk/utility/filesystem.cc
index e8c4c06edf5e02dbd14d11ebb93f4fa694fd4fdc..96459982fa42768d1c83de671fbb4883c3f85278 100644
--- a/dune/vtk/utility/filesystem.cc
+++ b/dune/vtk/utility/filesystem.cc
@@ -21,7 +21,7 @@
 template <class... Args>
 void inline _ignore_(Args&&...) {}
 
-namespace Dune { namespace filesystem {
+namespace Dune { namespace Vtk {
 
 std::string path::string() const
 {
@@ -42,7 +42,7 @@ void path::split(std::string p)
   bool relative = true;
 
   trim(p);
-  Dune::split(p.begin(), p.end(), separators.begin(), separators.end(),
+  Dune::Vtk::split(p.begin(), p.end(), separators.begin(), separators.end(),
     [this,&relative](auto first, auto end)
     {
       auto token = std::string(first, end);
@@ -187,4 +187,4 @@ path relative(path const& a, path const& b)
   return rel;
 }
 
-} } // end namespace Dune::filesystem
+} } // end namespace Dune::Vtk
diff --git a/dune/vtk/utility/filesystem.hh b/dune/vtk/utility/filesystem.hh
index 12aacfd31fe6a614b079e3ba149f51fd0cf1c9c4..ab7e148b7e7b6446ab6d0a163ddbedd95ab0692d 100644
--- a/dune/vtk/utility/filesystem.hh
+++ b/dune/vtk/utility/filesystem.hh
@@ -7,8 +7,9 @@
 
 namespace Dune
 {
-  namespace filesystem
+  namespace Vtk
   {
+
     // A minimalistic filesystem class
     class path
         : public std::vector<std::string>
@@ -130,5 +131,5 @@ namespace Dune
     /// Find the path of `a` relative to directory of `b`
     path relative(path const& a, path const& b);
 
-  } // end namespace filesystem
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/utility/lagrangepoints.hh b/dune/vtk/utility/lagrangepoints.hh
index b598c89527500007b91af4a8a4cfc46316678826..065c121d517f7e41718a75861b2995bf719deb94 100644
--- a/dune/vtk/utility/lagrangepoints.hh
+++ b/dune/vtk/utility/lagrangepoints.hh
@@ -7,22 +7,26 @@
 #include <dune/geometry/type.hh>
 #include <dune/localfunctions/lagrange/equidistantpoints.hh>
 
-namespace Dune {
+namespace Dune
+{
+
+namespace Vtk
+{
 
 namespace Impl {
   // forward declaration
   template <class K, unsigned int dim>
-  class VtkLagrangePointSetBuilder;
+  class LagrangePointSetBuilder;
 }
 
 
-/// \brief A set of lagrange points compatible with the numbering of VTK
+/// \brief A set of lagrange points compatible with the numbering of VTK and Gmsh
 /**
  * \tparam K    Field-type for the coordinates
  * \tparam dim  Dimension of the coordinates
  **/
 template <class K, unsigned int dim>
-class VtkLagrangePointSet
+class LagrangePointSet
     : public EmptyPointSet<K, dim>
 {
   using Super = EmptyPointSet<K, dim>;
@@ -30,7 +34,7 @@ class VtkLagrangePointSet
 public:
   static const unsigned int dimension = dim;
 
-  VtkLagrangePointSet (std::size_t order)
+  LagrangePointSet (std::size_t order)
     : Super(order)
   {
     assert(order > 0);
@@ -63,7 +67,7 @@ public:
 
 private:
   using Super::points_;
-  Impl::VtkLagrangePointSetBuilder<K,dim> builder_;
+  Impl::LagrangePointSetBuilder<K,dim> builder_;
 };
 
 
@@ -72,7 +76,7 @@ namespace Impl {
 // Build for lagrange point sets in different dimensions
 // Specialized for dim=1,2,3
 template <class K, unsigned int dim>
-class VtkLagrangePointSetBuilder
+class LagrangePointSetBuilder
 {
 public:
   template <class Points>
@@ -85,13 +89,13 @@ public:
 
 /**
  *  The implementation of the point set builder is directly derived from VTK.
- *  Modification are a change in data-types and naming scheme. Additionally 
+ *  Modification are a change in data-types and naming scheme. Additionally
  *  a LocalKey is created to reflect the concept of a Dune PointSet.
- * 
+ *
  *  Included is the license of the BSD 3-clause License included in the VTK
  *  source code from 2020/04/13 in commit b90dad558ce28f6d321420e4a6b17e23f5227a1c
  *  of git repository https://gitlab.kitware.com/vtk/vtk.
- * 
+ *
     Program:   Visualization Toolkit
     Module:    Copyright.txt
 
@@ -127,7 +131,7 @@ public:
 
 // Lagrange points on point geometries
 template <class K>
-class VtkLagrangePointSetBuilder<K,0>
+class LagrangePointSetBuilder<K,0>
 {
   static constexpr int dim = 0;
   using LP = LagrangePoint<K,dim>;
@@ -146,7 +150,7 @@ public:
 
 // Lagrange points on line geometries
 template <class K>
-class VtkLagrangePointSetBuilder<K,1>
+class LagrangePointSetBuilder<K,1>
 {
   static constexpr int dim = 1;
   using LP = LagrangePoint<K,dim>;
@@ -178,14 +182,14 @@ public:
 
 // Lagrange points on 2d geometries
 template <class K>
-class VtkLagrangePointSetBuilder<K,2>
+class LagrangePointSetBuilder<K,2>
 {
   static constexpr int dim = 2;
   using LP = LagrangePoint<K,dim>;
   using Vec = typename LP::Vector;
   using Key = LocalKey;
 
-  friend class VtkLagrangePointSetBuilder<K,3>;
+  friend class LagrangePointSetBuilder<K,3>;
 
 public:
   template <class Points>
@@ -236,7 +240,7 @@ private:
     };
 
     std::array<int,3> bindex;
-    
+
     double order_d = double(order);
     for (std::size_t p = 0; p < nPoints; ++p) {
       barycentricIndex(p, bindex, order);
@@ -315,7 +319,7 @@ private:
     const bool ibdy = (i == 0 || i == order[0]);
     const bool jbdy = (j == 0 || j == order[1]);
     const int nbdy = (ibdy ? 1 : 0) + (jbdy ? 1 : 0); // How many boundaries do we lie on at once?
-    
+
     int dof = 0;
     unsigned int entityIndex = 0;
     unsigned int index = 0;
@@ -347,7 +351,7 @@ private:
 
     // nbdy == 0: Face DOF
     dof = offset + (i - 1) + (order[0]-1) * ((j - 1));
-    Key innerKey = VtkLagrangePointSetBuilder<K,2>::calcQuadKey(i-1,j-1,{order[0]-2, order[1]-2}).second;
+    Key innerKey = LagrangePointSetBuilder<K,2>::calcQuadKey(i-1,j-1,{order[0]-2, order[1]-2}).second;
     return std::make_pair(dof, Key{0, dim-2, innerKey.index()});
   }
 };
@@ -355,7 +359,7 @@ private:
 
 // Lagrange points on 3d geometries
 template <class K>
-class VtkLagrangePointSetBuilder<K,3>
+class LagrangePointSetBuilder<K,3>
 {
   static constexpr int dim = 3;
   using LP = LagrangePoint<K,dim>;
@@ -419,7 +423,7 @@ private:
     };
 
     std::array<int,4> bindex;
-    
+
     double order_d = double(order);
     for (std::size_t p = 0; p < nPoints; ++p) {
       barycentricIndex(p, bindex, order);
@@ -482,11 +486,11 @@ private:
       if (order == 3)
         projectedBIndex[0] = projectedBIndex[1] = projectedBIndex[2] = 0;
       else
-        VtkLagrangePointSetBuilder<K,2>::barycentricIndex(vertexId, projectedBIndex, order-3);
+        LagrangePointSetBuilder<K,2>::barycentricIndex(vertexId, projectedBIndex, order-3);
 
       for (int i = 0; i < 3; i++)
         bindex[faceBCoords[faceId][i]] = (min + 1 + projectedBIndex[i]);
-        
+
       bindex[faceMinCoord[faceId]] = min;
     }
   }
@@ -501,7 +505,7 @@ private:
     points.resize(nPoints);
 
     std::array<int,3> orders{order, order, order};
-    std::array<Vec,8> nodes{Vec{0., 0., 0.}, Vec{1., 0., 0.}, Vec{1., 1., 0.}, Vec{0., 1., 0.}, 
+    std::array<Vec,8> nodes{Vec{0., 0., 0.}, Vec{1., 0., 0.}, Vec{1., 1., 0.}, Vec{0., 1., 0.},
                             Vec{0., 0., 1.}, Vec{1., 0., 1.}, Vec{1., 1., 1.}, Vec{0., 1., 1.}};
 
     for (int p = 0; p <= orders[2]; ++p) {
@@ -555,7 +559,7 @@ private:
         index = j;
         entityIndex += (j ? 3 : 2);
       }
-      else 
+      else
       { // !kbdy, On k axis
         offset += 4 * (order[0]-1) + 4 * (order[1]-1);
         dof = (k - 1) + (order[2]-1) * (i ? (j ? 3 : 1) : (j ? 2 : 0)) + offset;
@@ -573,7 +577,7 @@ private:
       {
         dof = (j - 1) + ((order[1]-1) * (k - 1)) + (i ? (order[1]-1) * (order[2]-1) : 0) + offset;
         entityIndex = (i ? 1 : 0);
-        faceKey = VtkLagrangePointSetBuilder<K,2>::calcQuadKey(j-1,k-1,{order[1]-2, order[2]-2}).second;
+        faceKey = LagrangePointSetBuilder<K,2>::calcQuadKey(j-1,k-1,{order[1]-2, order[2]-2}).second;
       }
       else {
         offset += 2 * (order[1] - 1) * (order[2] - 1);
@@ -581,14 +585,14 @@ private:
         {
           dof = (i - 1) + ((order[0]-1) * (k - 1)) + (j ? (order[2]-1) * (order[0]-1) : 0) + offset;
           entityIndex = (j ? 3 : 2);
-          faceKey = VtkLagrangePointSetBuilder<K,2>::calcQuadKey(i-1,k-1,{order[0]-2, order[2]-2}).second;
+          faceKey = LagrangePointSetBuilder<K,2>::calcQuadKey(i-1,k-1,{order[0]-2, order[2]-2}).second;
         }
-        else 
+        else
         { // kbdy, On k-normal face
           offset += 2 * (order[2]-1) * (order[0]-1);
           dof = (i - 1) + ((order[0]-1) * (j - 1)) + (k ? (order[0]-1) * (order[1]-1) : 0) + offset;
           entityIndex = (k ? 5 : 4);
-          faceKey = VtkLagrangePointSetBuilder<K,2>::calcQuadKey(i-1,j-1,{order[0]-2, order[1]-2}).second;
+          faceKey = LagrangePointSetBuilder<K,2>::calcQuadKey(i-1,j-1,{order[0]-2, order[1]-2}).second;
         }
       }
       return std::make_pair(dof, Key{entityIndex, dim-2, faceKey.index()});
@@ -597,9 +601,9 @@ private:
     // nbdy == 0: Body DOF
     offset += 2 * ((order[1]-1) * (order[2]-1) + (order[2]-1) * (order[0]-1) + (order[0]-1) * (order[1]-1));
     dof = offset + (i - 1) + (order[0]-1) * ((j - 1) + (order[1]-1) * ((k - 1)));
-    Key innerKey = VtkLagrangePointSetBuilder<K,3>::calcHexKey(i-1,j-1,k-1,{order[0]-2, order[1]-2, order[2]-2}).second;
+    Key innerKey = LagrangePointSetBuilder<K,3>::calcHexKey(i-1,j-1,k-1,{order[0]-2, order[1]-2, order[2]-2}).second;
     return std::make_pair(dof, Key{0, dim-3, innerKey.index()});
   }
 };
 
-}} // end namespace Dune::Impl
+}}} // end namespace Dune::Vtk::Impl
diff --git a/dune/vtk/utility/string.hh b/dune/vtk/utility/string.hh
index 7724911a9df42db7474b297b0f52424913242695..4cd7125a36045dfcb89d343dd1bcbd54c45c9a51 100644
--- a/dune/vtk/utility/string.hh
+++ b/dune/vtk/utility/string.hh
@@ -8,113 +8,118 @@
 
 namespace Dune
 {
-  /// convert all characters in a string to upper case
-  inline std::string to_upper(std::string input)
-  {
-    for (auto& c : input)
-      c = toupper(c);
-    return input;
-  }
 
-  /// convert all characters in a string to upper case
-  inline std::string to_lower(std::string input)
+  namespace Vtk
   {
-    for (auto& c : input)
-      c = tolower(c);
-    return input;
-  }
 
-  /// trim a string from the left
-  inline std::string& ltrim(std::string& str)
-  {
-    auto it =  std::find_if(str.begin(), str.end(), [](char ch)
+    /// convert all characters in a string to upper case
+    inline std::string to_upper(std::string input)
     {
-      return !std::isspace<char>(ch, std::locale::classic());
-    });
-    str.erase(str.begin() , it);
-    return str;
-  }
-
-  /// trim a string from the right
-  inline std::string& rtrim(std::string& str)
-  {
-    auto it =  std::find_if(str.rbegin(), str.rend(), [](char ch)
+      for (auto& c : input)
+        c = toupper(c);
+      return input;
+    }
+
+    /// convert all characters in a string to upper case
+    inline std::string to_lower(std::string input)
     {
-      return !std::isspace<char>(ch, std::locale::classic());
-    });
-    str.erase(it.base(), str.end());
-    return str;
-  }
-
-  /// trim a string from both sides
-  inline std::string& trim(std::string& str)
-  {
-    return ltrim(rtrim(str));
-  }
+      for (auto& c : input)
+        c = tolower(c);
+      return input;
+    }
 
-  /// trim a (copy of the) string from both sides
-  inline std::string trim_copy(std::string const& str)
-  {
-    auto s = str;
-    return trim(s);
-  }
+    /// trim a string from the left
+    inline std::string& ltrim(std::string& str)
+    {
+      auto it =  std::find_if(str.begin(), str.end(), [](char ch)
+      {
+        return !std::isspace<char>(ch, std::locale::classic());
+      });
+      str.erase(str.begin() , it);
+      return str;
+    }
 
+    /// trim a string from the right
+    inline std::string& rtrim(std::string& str)
+    {
+      auto it =  std::find_if(str.rbegin(), str.rend(), [](char ch)
+      {
+        return !std::isspace<char>(ch, std::locale::classic());
+      });
+      str.erase(it.base(), str.end());
+      return str;
+    }
 
-  template <class InputIter, class T, class Func>
-  void split(InputIter first, InputIter end, T const& t, Func f)
-  {
-    if (first == end)
-      return;
-
-    while (true) {
-      InputIter found = std::find(first, end, t);
-      f(first, found);
-      if (found == end)
-        break;
-      first = ++found;
+    /// trim a string from both sides
+    inline std::string& trim(std::string& str)
+    {
+      return ltrim(rtrim(str));
     }
-  }
 
-  template <class InputIter, class SeparatorIter, class Func>
-  void split(InputIter first, InputIter end, SeparatorIter s_first, SeparatorIter s_end, Func f)
-  {
-    if (first == end)
-      return;
-
-    while (true) {
-      InputIter found = std::find_first_of(first, end, s_first, s_end);
-      f(first, found);
-      if (found == end)
-        break;
-      first = ++found;
+    /// trim a (copy of the) string from both sides
+    inline std::string trim_copy(std::string const& str)
+    {
+      auto s = str;
+      return trim(s);
     }
-  }
 
-  /// Replace all occurences of substring `from` with `to` in source `str`.
-  inline void replaceAll(std::string& str, std::string const& from, std::string const& to)
-  {
-    if (from.empty())
-      return;
-    std::size_t start_pos = 0;
-    while ((start_pos = str.find(from, start_pos)) != std::string::npos)
+
+    template <class InputIter, class T, class Func>
+    void split(InputIter first, InputIter end, T const& t, Func f)
     {
-      str.replace(start_pos, from.length(), to);
-      start_pos += to.length();
+      if (first == end)
+        return;
+
+      while (true) {
+        InputIter found = std::find(first, end, t);
+        f(first, found);
+        if (found == end)
+          break;
+        first = ++found;
+      }
     }
-  }
 
+    template <class InputIter, class SeparatorIter, class Func>
+    void split(InputIter first, InputIter end, SeparatorIter s_first, SeparatorIter s_end, Func f)
+    {
+      if (first == end)
+        return;
+
+      while (true) {
+        InputIter found = std::find_first_of(first, end, s_first, s_end);
+        f(first, found);
+        if (found == end)
+          break;
+        first = ++found;
+      }
+    }
 
-  template <class InputIter>
-  std::string join (InputIter first, InputIter end, std::string sep = " ")
-  {
-    if (first == end)
-      return "";
+    /// Replace all occurences of substring `from` with `to` in source `str`.
+    inline void replaceAll(std::string& str, std::string const& from, std::string const& to)
+    {
+      if (from.empty())
+        return;
+      std::size_t start_pos = 0;
+      while ((start_pos = str.find(from, start_pos)) != std::string::npos)
+      {
+        str.replace(start_pos, from.length(), to);
+        start_pos += to.length();
+      }
+    }
 
-    std::ostringstream os;
-    os << *first++;
-    while (first != end)
-      os << sep << *first++;
-    return os.str();
-  }
 
+    template <class InputIter>
+    std::string join (InputIter first, InputIter end, std::string sep = " ")
+    {
+      if (first == end)
+        return "";
+
+      std::ostringstream os;
+      os << *first++;
+      while (first != end)
+        os << sep << *first++;
+      return os.str();
+    }
+
+  } // end namespace Vtk
 } // end namspace Dune
diff --git a/dune/vtk/utility/test/test-lagrange.cc b/dune/vtk/utility/test/test-lagrange.cc
index 62ed50a11e6a44cb4bc73ace86ccffc31a2b2fdc..f25c2cb55826ebfc2c5435e6d14c6904bbf69dd3 100644
--- a/dune/vtk/utility/test/test-lagrange.cc
+++ b/dune/vtk/utility/test/test-lagrange.cc
@@ -75,8 +75,8 @@ bool test(const Basis &basis, const Points &points, bool verbose)
 template <class Topology>
 bool test(unsigned int order, bool verbose = false)
 {
-  typedef Dune::LagrangeBasisFactory<Dune::VtkLagrangePointSet,Topology::dimension,StorageField,ComputeField> BasisFactory;
-  typedef Dune::LagrangeCoefficientsFactory< Dune::VtkLagrangePointSet,  Topology::dimension,double > LagrangeCoefficientsFactory;
+  typedef Dune::LagrangeBasisFactory<Dune::Vtk::LagrangePointSet,Topology::dimension,StorageField,ComputeField> BasisFactory;
+  typedef Dune::LagrangeCoefficientsFactory< Dune::Vtk::LagrangePointSet,  Topology::dimension,double > LagrangeCoefficientsFactory;
 
   bool ret = true;
 
@@ -154,14 +154,14 @@ int main ( int argc, char **argv )
 #endif
 
 #ifdef CHECKDIM2
-  tests &= test<Prism<Prism<Point> > > (order); // quad
-  tests &= test<Pyramid<Pyramid<Point> > >(order); // triangle
+  //tests &= test<Prism<Prism<Point> > > (order); // quad
+  //tests &= test<Pyramid<Pyramid<Point> > >(order); // triangle
 #endif
 
 #ifdef CHECKDIM3
-  tests &= test<Prism<Prism<Prism<Point> > > >(order); // hexahedron
+  // tests &= test<Prism<Prism<Prism<Point> > > >(order); // hexahedron
   // tests &= test<Prism<Pyramid<Pyramid<Point> > > >(order);
-  tests &= test<Pyramid<Pyramid<Pyramid<Point> > > >(order); // tetrahedron
+  // tests &= test<Pyramid<Pyramid<Pyramid<Point> > > >(order); // tetrahedron
 #endif
 
   return (tests ? 0 : 1);
diff --git a/dune/vtk/utility/uid.hh b/dune/vtk/utility/uid.hh
index 22c2bb7e7bf4eb72c4aa2062cea325b39c428825..6a686507403c621b07a337fb3ee2ffd5bd19f087 100644
--- a/dune/vtk/utility/uid.hh
+++ b/dune/vtk/utility/uid.hh
@@ -7,16 +7,20 @@
 
 namespace Dune
 {
-  inline std::string uid (std::size_t len = 8)
+  namespace Vtk
   {
-    static const auto digits = "0123456789abcdefghijklmnopqrstuvwxyzABCDEFGHIJKLMNOPQRSTUVWXYZ";
-    static const int N = std::strlen(digits);
 
-    std::string id(len,' ');
-    for (std::size_t i = 0; i < len; ++i)
-      id[i] = digits[std::rand()%N];
+    inline std::string uid (std::size_t len = 8)
+    {
+      static const auto digits = "0123456789abcdefghijklmnopqrstuvwxyzABCDEFGHIJKLMNOPQRSTUVWXYZ";
+      static const int N = std::strlen(digits);
 
-    return id;
-  }
+      std::string id(len,' ');
+      for (std::size_t i = 0; i < len; ++i)
+        id[i] = digits[std::rand()%N];
 
+      return id;
+    }
+
+  } // end namespace Vtk
 } // end namespace Dune
diff --git a/dune/vtk/vtkfunction.hh b/dune/vtk/vtkfunction.hh
deleted file mode 100644
index d51f64c9af1a69c6affc6ff4481031a289fdeba9..0000000000000000000000000000000000000000
--- a/dune/vtk/vtkfunction.hh
+++ /dev/null
@@ -1,128 +0,0 @@
-#pragma once
-
-#include <optional>
-#include <type_traits>
-
-#include <dune/common/std/type_traits.hh>
-
-#include "vtklocalfunction.hh"
-#include "vtktypes.hh"
-
-namespace Dune
-{
-  template <class T, int N>
-  class FieldVector;
-
-  template <class T, int N, int M>
-  class FieldMatrix;
-
-
-  /// Wrapper class for functions allowing local evaluations.
-  template <class GridView>
-  class VtkFunction
-  {
-    template <class F>
-    using LocalFunction = decltype(localFunction(std::declval<F>()));
-
-    using Domain = typename GridView::template Codim<0>::Entity::Geometry::LocalCoordinate;
-
-    template <class F>
-    using Range = std::decay_t<std::result_of_t<F(Domain)>>;
-
-  private:
-
-    template <class T, int N>
-    static auto sizeOfImpl (FieldVector<T,N> const&)
-      -> std::integral_constant<int, N> { return {}; }
-
-    template <class T, int N, int M>
-    static auto sizeOfImpl (FieldMatrix<T,N,M> const&)
-      -> std::integral_constant<int, N*M> { return {}; }
-
-    static auto sizeOfImpl (...)
-      -> std::integral_constant<int, 1> { return {}; }
-
-    template <class T>
-    static constexpr int sizeOf () { return decltype(sizeOfImpl(std::declval<T>()))::value; }
-
-  public:
-    /// Constructor VtkFunction from legacy VTKFunction
-    /**
-     * \param fct   The VTKFunction to wrap
-     * \param type  The VTK datatype how to write the function values to the output [Vtk::FLOAT64]
-     **/
-    VtkFunction (std::shared_ptr<VTKFunction<GridView> const> const& fct,
-                 std::optional<Vtk::DataTypes> type = {})
-      : localFct_(fct)
-      , name_(fct->name())
-      , ncomps_(fct->ncomps())
-      , type_(type ? *type : Vtk::FLOAT64)
-    {}
-
-    /// Construct VtkFunction from dune-functions GridFunction with Signature
-    // NOTE: Stores the localFunction(fct) by value.
-    /**
-     * \param fct     A Grid(View)-function, providing a `localFunction(fct)`
-     * \param name    The name to use component identification in the VTK file
-     * \param ncomps  Number of components of the pointwise data. Is extracted
-     *                from the range type of the GridFunction if not given.
-     * \param type    The \ref Vtk::DataTypes used in the output. E.g. FLOAT32,
-     *                or FLOAT64. Is extracted from the range type of the
-     *                GridFunction if not given.
-     **/
-    template <class F,
-      class = void_t<LocalFunction<F>> >
-    VtkFunction (F&& fct, std::string name,
-                 std::optional<int> ncomps = {},
-                 std::optional<Vtk::DataTypes> type = {})
-      : localFct_(localFunction(std::forward<F>(fct)))
-      , name_(std::move(name))
-    {
-      using R = Range<LocalFunction<F>>;
-
-      ncomps_ = ncomps ? *ncomps : sizeOf<R>();
-      type_ = type ? *type : Vtk::Map::type<R>();
-    }
-
-    /// Constructor that forward the number of components and data type to the other constructor
-    template <class F,
-      class = void_t<LocalFunction<F>> >
-    VtkFunction (F&& fct, Vtk::FieldInfo fieldInfo,
-                 std::optional<Vtk::DataTypes> type = {})
-      : VtkFunction(std::forward<F>(fct), fieldInfo.name(), fieldInfo.ncomps(), type)
-    {}
-
-    VtkFunction () = default;
-
-    /// Create a LocalFunction
-    friend VtkLocalFunction<GridView> localFunction (VtkFunction const& self)
-    {
-      return self.localFct_;
-    }
-
-    /// Return a name associated with the function
-    std::string const& name () const
-    {
-      return name_;
-    }
-
-    /// Return the number of components of the Range
-    int ncomps () const
-    {
-      return ncomps_ > 3 ? 9 : ncomps_ > 1 ? 3 : 1; // tensor, vector, scalar
-    }
-
-    /// Return the VTK Datatype associated with the functions range type
-    Vtk::DataTypes type () const
-    {
-      return type_;
-    }
-
-  private:
-    VtkLocalFunction<GridView> localFct_;
-    std::string name_;
-    int ncomps_ = 1;
-    Vtk::DataTypes type_ = Vtk::FLOAT32;
-  };
-
-} // end namespace Dune
diff --git a/dune/vtk/vtklocalfunction.hh b/dune/vtk/vtklocalfunction.hh
deleted file mode 100644
index e636b0b4396adb03db94ab4366227ed35cce7c3d..0000000000000000000000000000000000000000
--- a/dune/vtk/vtklocalfunction.hh
+++ /dev/null
@@ -1,71 +0,0 @@
-#pragma once
-
-#include <memory>
-#include <type_traits>
-
-#include <dune/common/std/type_traits.hh>
-
-#include "vtklocalfunctioninterface.hh"
-#include "legacyvtkfunction.hh"
-#include "defaultvtkfunction.hh"
-
-namespace Dune
-{
-  /// \brief A VtkLocalFunction is a function-like object that can be bound to a grid element
-  /// an that provides an evaluate method with a component argument.
-  /**
-   * Stores internally a VtkLocalFunctionInterface object for the concrete evaluation.
-   **/
-  template <class GridView>
-  class VtkLocalFunction
-  {
-    using Self = VtkLocalFunction;
-    using Entity = typename GridView::template Codim<0>::Entity;
-    using LocalCoordinate = typename Entity::Geometry::LocalCoordinate;
-
-    template <class LF, class E>
-    using HasBind = decltype(std::declval<LF>().bind(std::declval<E>()));
-
-  public:
-    /// Construct the VtkLocalFunction from any function object that has a bind(element) method.
-    template <class LF,
-      disableCopyMove<Self, LF> = 0,
-      class = void_t<HasBind<LF,Entity>> >
-    VtkLocalFunction (LF&& lf)
-      : localFct_(std::make_shared<LocalFunctionWrapper<GridView,LF>>(std::forward<LF>(lf)))
-    {}
-
-    /// Construct a VtkLocalFunction from a legacy VTKFunction
-    VtkLocalFunction (std::shared_ptr<VTKFunction<GridView> const> const& lf)
-      : localFct_(std::make_shared<VTKLocalFunctionWrapper<GridView>>(lf))
-    {}
-
-    /// Allow the default construction of a VtkLocalFunction
-    VtkLocalFunction () = default;
-
-    /// Bind the function to the grid entity
-    void bind (Entity const& entity)
-    {
-      assert(bool(localFct_));
-      localFct_->bind(entity);
-    }
-
-    /// Unbind from the currently bound entity
-    void unbind ()
-    {
-      assert(bool(localFct_));
-      localFct_->unbind();
-    }
-
-    /// Evaluate the `comp` component of the Range value at local coordinate `xi`
-    double evaluate (int comp, LocalCoordinate const& xi) const
-    {
-      assert(bool(localFct_));
-      return localFct_->evaluate(comp, xi);
-    }
-
-  private:
-    std::shared_ptr<VtkLocalFunctionInterface<GridView>> localFct_ = nullptr;
-  };
-
-} // end namespace Dune
diff --git a/dune/vtk/vtklocalfunctioninterface.hh b/dune/vtk/vtklocalfunctioninterface.hh
deleted file mode 100644
index 9ba2afeccea47eb143edfab28d2d84e35cd0bc1e..0000000000000000000000000000000000000000
--- a/dune/vtk/vtklocalfunctioninterface.hh
+++ /dev/null
@@ -1,27 +0,0 @@
-#pragma once
-
-namespace Dune
-{
-  /// \brief An abstract base class for LocalFunctions that can be bound to an element and
-  /// evaluated in local coordinates w.r.t. to a component of its value.
-  template <class GridView>
-  class VtkLocalFunctionInterface
-  {
-  public:
-    using Entity = typename GridView::template Codim<0>::Entity;
-    using LocalCoordinate = typename Entity::Geometry::LocalCoordinate;
-
-    /// Bind the function to the grid entity
-    virtual void bind (Entity const& entity) = 0;
-
-    /// Unbind from the currently bound entity
-    virtual void unbind () = 0;
-
-    /// Evaluate single component comp in the entity at local coordinates xi
-    virtual double evaluate (int comp, LocalCoordinate const& xi) const = 0;
-
-    /// Virtual destructor
-    virtual ~VtkLocalFunctionInterface () = default;
-  };
-
-} // end namespace Dune
diff --git a/dune/vtk/vtktypes.cc b/dune/vtk/vtktypes.cc
deleted file mode 100644
index 90f993ec79f5442974039d02f6bac52435e850db..0000000000000000000000000000000000000000
--- a/dune/vtk/vtktypes.cc
+++ /dev/null
@@ -1,200 +0,0 @@
-#include "vtktypes.hh"
-
-#include <iostream>
-
-#include <dune/common/exceptions.hh>
-
-namespace Dune {
-namespace Vtk {
-
-std::string to_string (FormatTypes type)
-{
-  switch (type) {
-    case ASCII:      return "ascii";
-    case BINARY:     return "binary";
-    case COMPRESSED: return "compressed";
-    case APPENDED:   return "appended";
-    default:
-      DUNE_THROW(RangeError, "FormatType not found.");
-      std::abort();
-  }
-}
-
-std::string to_string (DataTypes type)
-{
-  switch (type) {
-    case INT8:    return "Int8";
-    case UINT8:   return "UInt8";
-    case INT16:   return "Int16";
-    case UINT16:  return "UInt16";
-    case INT32:   return "Int32";
-    case UINT32:  return "UInt32";
-    case INT64:   return "Int64";
-    case UINT64:  return "UInt64";
-    case FLOAT32: return "Float32";
-    case FLOAT64: return "Float64";
-    default:
-      DUNE_THROW(RangeError, "DataType not found.");
-      std::abort();
-  }
-}
-
-GeometryType to_geometry (std::uint8_t cell)
-{
-  switch (cell) {
-    case VERTEX:     return GeometryTypes::vertex;
-    case LINE:       return GeometryTypes::line;
-    case TRIANGLE:   return GeometryTypes::triangle;
-    case QUAD:       return GeometryTypes::quadrilateral;
-    case TETRA:      return GeometryTypes::tetrahedron;
-    case HEXAHEDRON: return GeometryTypes::hexahedron;
-    case WEDGE:      return GeometryTypes::prism;
-    case PYRAMID:    return GeometryTypes::pyramid;
-
-    // Quadratic VTK cell types
-    case QUADRATIC_EDGE:        return GeometryTypes::line;
-    case QUADRATIC_TRIANGLE:    return GeometryTypes::triangle;
-    case QUADRATIC_QUAD:        return GeometryTypes::quadrilateral;
-    case QUADRATIC_TETRA:       return GeometryTypes::tetrahedron;
-    case QUADRATIC_HEXAHEDRON:  return GeometryTypes::hexahedron;
-    
-    // Arbitrary order Lagrange elements
-    case LAGRANGE_CURVE:        return GeometryTypes::line;
-    case LAGRANGE_TRIANGLE:     return GeometryTypes::triangle;
-    case LAGRANGE_QUADRILATERAL:return GeometryTypes::quadrilateral;
-    case LAGRANGE_TETRAHEDRON:  return GeometryTypes::tetrahedron;
-    case LAGRANGE_HEXAHEDRON:   return GeometryTypes::hexahedron;
-    case LAGRANGE_WEDGE:        return GeometryTypes::prism;
-    default:
-      DUNE_THROW(RangeError, "CellType does not map to GeometryType.");
-      std::abort();
-  }
-}
-
-
-std::map<std::string, DataTypes> Map::to_datatype = {
-  {"Int8",    INT8},
-  {"UInt8",   UINT8},
-  {"Int16",   INT16},
-  {"UInt16",  UINT16},
-  {"Int32",   INT32},
-  {"UInt32",  UINT32},
-  {"Int64",   INT64},
-  {"UInt64",  UINT64},
-  {"Float32", FLOAT32},
-  {"Float64", FLOAT64}
-};
-
-
-
-CellType::CellType (GeometryType const& t, CellParametrization parametrization)
-  : noPermutation_(true)
-{
-  if (parametrization == LINEAR) {
-    if (t.isVertex()) {
-      type_ = VERTEX;
-      permutation_ = {0};
-    }
-    else if (t.isLine()) {
-      type_ = LINE;
-      permutation_ = {0,1};
-    }
-    else if (t.isTriangle()) {
-      type_ = TRIANGLE;
-      permutation_ = {0,1,2};
-    }
-    else if (t.isQuadrilateral()) {
-      type_ = QUAD;
-      permutation_ = {0,1,3,2};
-      noPermutation_ = false;
-    }
-    else if (t.isTetrahedron()) {
-      type_ = TETRA;
-      permutation_ = {0,1,2,3};
-    }
-    else if (t.isHexahedron()) {
-      type_ = HEXAHEDRON;
-      permutation_ = {0,1,3,2,4,5,7,6};
-      noPermutation_ = false;
-    }
-    else if (t.isPrism()) {
-      type_ = WEDGE;
-      permutation_ = {0,2,1,3,5,4};
-      noPermutation_ = false;
-    }
-    else if (t.isPyramid()) {
-      type_ = PYRAMID;
-      permutation_ = {0,1,3,2,4};
-      noPermutation_ = false;
-    }
-    else if (t.isNone() && t.dim() == 1) {
-      type_ = LINE;
-      permutation_ = {0,1};
-    }
-    else if (t.isNone() && t.dim() == 2) {
-      type_ = POLYGON;
-      permutation_ = {0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19};
-    }
-    else {
-      std::cerr << "Geometry Type not supported by VTK!\n";
-      std::abort();
-    }
-  } else if (parametrization == QUADRATIC) {
-    if (t.isLine()) {
-      type_ = QUADRATIC_EDGE;
-      permutation_ = {0,1, 0};
-    }
-    else if (t.isTriangle()) {
-      type_ = QUADRATIC_TRIANGLE;
-      permutation_ = {0,1,2, 0,2,1};
-      noPermutation_ = false;
-    }
-    else if (t.isQuadrilateral()) {
-      type_ = QUADRATIC_QUAD;
-      permutation_ = {0,1,3,2, 2,1,3,0};
-      noPermutation_ = false;
-    }
-    else if (t.isTetrahedron()) {
-      type_ = QUADRATIC_TETRA;
-      permutation_ = {0,1,2,3, 0,2,1,3,4,5};
-      noPermutation_ = false;
-    }
-    else if (t.isHexahedron()) {
-      type_ = QUADRATIC_HEXAHEDRON;
-      permutation_ = {0,1,3,2,4,5,7,6, 6,5,7,4,10,9,11,8,0,1,3,2};
-      noPermutation_ = false;
-    }
-    else {
-      std::cerr << "Geometry Type not supported by VTK!\n";
-      std::abort();
-    }
-  } else if (parametrization == LAGRANGE) {
-    if (t.isLine()) {
-      type_ = LAGRANGE_CURVE;
-    }
-    else if (t.isTriangle()) {
-      type_ = LAGRANGE_TRIANGLE;
-    }
-    else if (t.isQuadrilateral()) {
-      type_ = LAGRANGE_QUADRILATERAL;
-    }
-    else if (t.isTetrahedron()) {
-      type_ = LAGRANGE_TETRAHEDRON;
-    }
-    else if (t.isHexahedron()) {
-      type_ = LAGRANGE_HEXAHEDRON;
-    }
-    else if (t.isPrism()) {
-      type_ = LAGRANGE_WEDGE;
-    }
-    else if (t.isPyramid()) {
-      type_ = LAGRANGE_PYRAMID;
-    }
-    else {
-      std::cerr << "Geometry Type not supported by VTK!\n";
-      std::abort();
-    }
-  }
-}
-
-}} // end namespace Dune::Vtk
diff --git a/dune/vtk/vtkwriter.hh b/dune/vtk/writer.hh
similarity index 81%
rename from dune/vtk/vtkwriter.hh
rename to dune/vtk/writer.hh
index b360f86ea7e657a8da9ea1f290842c4d5a694ecb..dd350212a854ada52a06b4d13a6e460c12fd2b4a 100644
--- a/dune/vtk/vtkwriter.hh
+++ b/dune/vtk/writer.hh
@@ -1,9 +1,9 @@
 #pragma once
 
-#include <dune/vtk/writers/vtkimagedatawriter.hh>
-#include <dune/vtk/writers/vtkrectilineargridwriter.hh>
-#include <dune/vtk/writers/vtkstructuredgridwriter.hh>
-#include <dune/vtk/writers/vtkunstructuredgridwriter.hh>
+#include <dune/vtk/writers/imagedatawriter.hh>
+#include <dune/vtk/writers/rectilineargridwriter.hh>
+#include <dune/vtk/writers/structuredgridwriter.hh>
+#include <dune/vtk/writers/unstructuredgridwriter.hh>
 
 #if HAVE_DUNE_SPGRID
 #include <dune/grid/spgrid.hh>
@@ -40,7 +40,7 @@ namespace Dune
       : public Impl::VtkWriterImpl<GridView, typename GridView::Grid>::type
   {
     using Super = typename Impl::VtkWriterImpl<GridView, typename GridView::Grid>::type;
-    
+
   public:
     using Super::Super;
   };
@@ -58,7 +58,7 @@ namespace Dune
     template <class GridView, int dim, class Coordinates>
     struct VtkWriterImpl<GridView, YaspGrid<dim,Coordinates>>
     {
-      using type = VtkRectilinearGridWriter<GridView, YaspDataCollector<GridView>>;
+      using type = VtkRectilinearGridWriter<GridView, Vtk::YaspDataCollector<GridView>>;
     };
 
 #if HAVE_DUNE_SPGRID
@@ -66,7 +66,7 @@ namespace Dune
     template <class GridView, class ct, int dim, template <int> class Ref, class Comm>
     struct VtkWriterImpl<GridView, SPGrid<ct,dim,Ref,Comm>>
     {
-      using type = VtkRectilinearGridWriter<GridView, SPDataCollector<GridView>>;
+      using type = VtkRectilinearGridWriter<GridView, Vtk::SPDataCollector<GridView>>;
     };
 #endif
 
@@ -74,7 +74,7 @@ namespace Dune
     template <class GridView, int dim, class ct>
     struct VtkWriterImpl<GridView, YaspGrid<dim,TensorProductCoordinates<ct,dim>>>
     {
-      using type = VtkRectilinearGridWriter<GridView, YaspDataCollector<GridView>>;
+      using type = VtkRectilinearGridWriter<GridView, Vtk::YaspDataCollector<GridView>>;
     };
 
     // A transformed structured grid has structured connectivity but unstructured point
@@ -82,7 +82,7 @@ namespace Dune
     template <class GridView, int dim, class Coordinates, class CoordFunction, class Allocator>
     struct VtkWriterImpl<GridView, GeometryGrid<YaspGrid<dim,Coordinates>, CoordFunction, Allocator>>
     {
-      using type = VtkStructuredGridWriter<GridView, YaspDataCollector<GridView>>;
+      using type = VtkStructuredGridWriter<GridView, Vtk::YaspDataCollector<GridView>>;
     };
 
   } // end namespace Impl
diff --git a/dune/vtk/vtkwriterinterface.hh b/dune/vtk/writerinterface.hh
similarity index 95%
rename from dune/vtk/vtkwriterinterface.hh
rename to dune/vtk/writerinterface.hh
index 8dadb0a5ed609c1c8d9af6bb7b066935fee4c12c..7e8053821deb72dcb0bb05c4cc1fdc9f357fcee4 100644
--- a/dune/vtk/vtkwriterinterface.hh
+++ b/dune/vtk/writerinterface.hh
@@ -10,8 +10,8 @@
 #include <dune/common/parallel/mpihelper.hh>
 #include <dune/vtk/filewriter.hh>
 #include <dune/vtk/forward.hh>
-#include <dune/vtk/vtkfunction.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/function.hh>
+#include <dune/vtk/types.hh>
 
 namespace Dune
 {
@@ -22,15 +22,15 @@ namespace Dune
    **/
   template <class GridView, class DataCollector>
   class VtkWriterInterface
-      : public FileWriter
+      : public Vtk::FileWriter
   {
-    template <class> friend class VtkTimeseriesWriter;
+    template <class> friend class TimeseriesWriter;
     template <class> friend class PvdWriter;
 
   protected:
     static constexpr int dimension = GridView::dimension;
 
-    using VtkFunction = Dune::VtkFunction<GridView>;
+    using VtkFunction = Dune::Vtk::Function<GridView>;
     using Communicator = CollectiveCommunication<typename MPIHelper::MPICommunicator>;
     using pos_type = typename std::ostream::pos_type;
 
@@ -81,7 +81,7 @@ namespace Dune
     /**
      * \param fn   Filename of the VTK file. May contain a directory and any file extension.
      * \param dir  The optional parameter specifies the directory of the partition files for parallel writes.
-     * 
+     *
      * \returns File name that is actually written.
      **/
     virtual std::string write (std::string const& fn, std::optional<std::string> dir = {}) const override;
@@ -90,8 +90,8 @@ namespace Dune
     /**
      * Attach a global function to the writer that will be evaluated at grid points
      * (vertices and higher order points). The global function must be
-     * assignable to the function wrapper \ref VtkFunction. Additional argument
-     * for output datatype and number of components can be passed. See \ref VtkFunction
+     * assignable to the function wrapper \ref Vtk::Function. Additional argument
+     * for output datatype and number of components can be passed. See \ref Vtk::Function
      * Constructor for possible arguments.
      **/
     template <class Function, class... Args>
@@ -104,9 +104,9 @@ namespace Dune
     /// \brief Attach cell data to the writer
     /**
      * Attach a global function to the writer that will be evaluated at cell centers.
-     * The global function must be assignable to the function wrapper \ref VtkFunction.
+     * The global function must be assignable to the function wrapper \ref Vtk::Function.
      * Additional argument for output datatype and number of components can be passed.
-     * See \ref VtkFunction Constructor for possible arguments.
+     * See \ref Vtk::Function Constructor for possible arguments.
      **/
     template <class Function, class... Args>
     VtkWriterInterface& addCellData (Function&& fct, Args&&... args)
@@ -120,7 +120,7 @@ namespace Dune
     void setFormat (Vtk::FormatTypes format)
     {
       format_ = format;
-      
+
 #if !HAVE_VTK_ZLIB
       if (format_ == Vtk::COMPRESSED) {
         std::cout << "Dune is compiled without compression. Falling back to BINARY VTK output!\n";
@@ -247,4 +247,4 @@ namespace Dune
 
 } // end namespace Dune
 
-#include "vtkwriterinterface.impl.hh"
+#include "writerinterface.impl.hh"
diff --git a/dune/vtk/vtkwriterinterface.impl.hh b/dune/vtk/writerinterface.impl.hh
similarity index 98%
rename from dune/vtk/vtkwriterinterface.impl.hh
rename to dune/vtk/writerinterface.impl.hh
index 597774c9cbacb45aedd3269e8f6d28b7777a0e39..fbd7d5010ca1d5f98287751512dc3d7a98a41e80 100644
--- a/dune/vtk/vtkwriterinterface.impl.hh
+++ b/dune/vtk/writerinterface.impl.hh
@@ -28,13 +28,13 @@ std::string VtkWriterInterface<GV,DC>
 {
   dataCollector_->update();
 
-  auto p = filesystem::path(fn);
+  auto p = Vtk::path(fn);
   auto name = p.stem();
   p.remove_filename();
 
-  filesystem::path fn_dir = p;
-  filesystem::path data_dir = dir ? filesystem::path(*dir) : fn_dir;
-  filesystem::path rel_dir = filesystem::relative(data_dir, fn_dir);
+  Vtk::path fn_dir = p;
+  Vtk::path data_dir = dir ? Vtk::path(*dir) : fn_dir;
+  Vtk::path rel_dir = Vtk::relative(data_dir, fn_dir);
 
   std::string serial_fn = data_dir.string() + '/' + name.string();
   std::string parallel_fn = fn_dir.string() + '/' + name.string();
@@ -44,7 +44,7 @@ std::string VtkWriterInterface<GV,DC>
     serial_fn += "_p" + std::to_string(comm().rank());
 
   std::string outputFilename;
-  
+
   { // write serial file
     outputFilename = serial_fn + "." + fileExtension();
     std::ofstream serial_out(outputFilename, std::ios_base::ate | std::ios::binary);
diff --git a/dune/vtk/writers/vtkimagedatawriter.hh b/dune/vtk/writers/imagedatawriter.hh
similarity index 91%
rename from dune/vtk/writers/vtkimagedatawriter.hh
rename to dune/vtk/writers/imagedatawriter.hh
index ed7cf1fc2027909dffb6e0ac6a09241c94151734..e512cb4c75092c772f0000781d4190bfcf813167 100644
--- a/dune/vtk/writers/vtkimagedatawriter.hh
+++ b/dune/vtk/writers/imagedatawriter.hh
@@ -6,11 +6,11 @@
 
 #include <dune/vtk/filewriter.hh>
 #include <dune/vtk/forward.hh>
-#include <dune/vtk/vtkfunction.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/function.hh>
+#include <dune/vtk/types.hh>
 #include <dune/vtk/datacollectors/structureddatacollector.hh>
 
-#include <dune/vtk/vtkwriterinterface.hh>
+#include <dune/vtk/writerinterface.hh>
 
 namespace Dune
 {
@@ -62,7 +62,7 @@ namespace Dune
   template <class GridView,
     class = std::void_t<typename GridView::IndexSet>>
   VtkImageDataWriter(GridView const&, Vtk::FormatTypes = Vtk::BINARY, Vtk::DataTypes = Vtk::FLOAT32)
-    -> VtkImageDataWriter<GridView, StructuredDataCollector<GridView>>;
+    -> VtkImageDataWriter<GridView, Vtk::StructuredDataCollector<GridView>>;
 
   template <class DataCollector,
     class = std::void_t<typename DataCollector::GridView>>
@@ -76,4 +76,4 @@ namespace Dune
 
 } // end namespace Dune
 
-#include "vtkimagedatawriter.impl.hh"
+#include "imagedatawriter.impl.hh"
diff --git a/dune/vtk/writers/vtkimagedatawriter.impl.hh b/dune/vtk/writers/imagedatawriter.impl.hh
similarity index 82%
rename from dune/vtk/writers/vtkimagedatawriter.impl.hh
rename to dune/vtk/writers/imagedatawriter.impl.hh
index 898b076ca5fdfc0d43988cf5057b9b854a111cc8..03c1ea196b580874b2b7357320cdcd1109fa67bf 100644
--- a/dune/vtk/writers/vtkimagedatawriter.impl.hh
+++ b/dune/vtk/writers/imagedatawriter.impl.hh
@@ -27,13 +27,13 @@ void VtkImageDataWriter<GV,DC>
   auto const& origin = dataCollector_->origin();
   auto const& spacing = dataCollector_->spacing();
   out << "<ImageData"
-      << " WholeExtent=\"" << join(wholeExtent.begin(), wholeExtent.end()) << "\""
-      << " Origin=\"" << join(origin.begin(), origin.end()) << "\""
-      << " Spacing=\"" << join(spacing.begin(), spacing.end()) << "\""
+      << " WholeExtent=\"" << Vtk::join(wholeExtent.begin(), wholeExtent.end()) << "\""
+      << " Origin=\"" << Vtk::join(origin.begin(), origin.end()) << "\""
+      << " Spacing=\"" << Vtk::join(spacing.begin(), spacing.end()) << "\""
       << ">\n";
 
   dataCollector_->writeLocalPiece([&out](auto const& extent) {
-    out << "<Piece Extent=\"" << join(extent.begin(), extent.end()) << "\">\n";
+    out << "<Piece Extent=\"" << Vtk::join(extent.begin(), extent.end()) << "\">\n";
   });
 
   // Write data associated with grid points
@@ -67,9 +67,9 @@ void VtkImageDataWriter<GV,DC>
   auto const& spacing = dataCollector_->spacing();
   out << "<PImageData"
       << " GhostLevel=\"" << dataCollector_->ghostLevel() << "\""
-      << " WholeExtent=\"" << join(wholeExtent.begin(), wholeExtent.end()) << "\""
-      << " Origin=\"" << join(origin.begin(), origin.end()) << "\""
-      << " Spacing=\"" << join(spacing.begin(), spacing.end()) << "\""
+      << " WholeExtent=\"" << Vtk::join(wholeExtent.begin(), wholeExtent.end()) << "\""
+      << " Origin=\"" << Vtk::join(origin.begin(), origin.end()) << "\""
+      << " Spacing=\"" << Vtk::join(spacing.begin(), spacing.end()) << "\""
       << ">\n";
 
   // Write data associated with grid points
@@ -100,7 +100,7 @@ void VtkImageDataWriter<GV,DC>
     std::string piece_source = pfilename + "_p" + std::to_string(p) + "." + ext;
     out << "<Piece Source=\"" << piece_source << "\"";
     if (write_extent)
-      out << " Extent=\"" << join(extent.begin(), extent.end()) << "\"";
+      out << " Extent=\"" << Vtk::join(extent.begin(), extent.end()) << "\"";
      out << " />\n";
   });
 
diff --git a/dune/vtk/writers/vtkrectilineargridwriter.hh b/dune/vtk/writers/rectilineargridwriter.hh
similarity index 92%
rename from dune/vtk/writers/vtkrectilineargridwriter.hh
rename to dune/vtk/writers/rectilineargridwriter.hh
index f119a0501f54adff0987b4e7112fb4d386f8300a..65f4cd1e6cb5a70dae61a2908d3db6d3b5a1b5fe 100644
--- a/dune/vtk/writers/vtkrectilineargridwriter.hh
+++ b/dune/vtk/writers/rectilineargridwriter.hh
@@ -6,11 +6,11 @@
 
 #include <dune/vtk/filewriter.hh>
 #include <dune/vtk/forward.hh>
-#include <dune/vtk/vtkfunction.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/function.hh>
+#include <dune/vtk/types.hh>
 #include <dune/vtk/datacollectors/structureddatacollector.hh>
 
-#include <dune/vtk/vtkwriterinterface.hh>
+#include <dune/vtk/writerinterface.hh>
 
 namespace Dune
 {
@@ -68,7 +68,7 @@ namespace Dune
   template <class GridView,
     class = std::void_t<typename GridView::IndexSet>>
   VtkRectilinearGridWriter(GridView const&, Vtk::FormatTypes = Vtk::BINARY, Vtk::DataTypes = Vtk::FLOAT32)
-    -> VtkRectilinearGridWriter<GridView, StructuredDataCollector<GridView>>;
+    -> VtkRectilinearGridWriter<GridView, Vtk::StructuredDataCollector<GridView>>;
 
   template <class DataCollector,
     class = std::void_t<typename DataCollector::GridView>>
@@ -82,4 +82,4 @@ namespace Dune
 
 } // end namespace Dune
 
-#include "vtkrectilineargridwriter.impl.hh"
+#include "rectilineargridwriter.impl.hh"
diff --git a/dune/vtk/writers/vtkrectilineargridwriter.impl.hh b/dune/vtk/writers/rectilineargridwriter.impl.hh
similarity index 94%
rename from dune/vtk/writers/vtkrectilineargridwriter.impl.hh
rename to dune/vtk/writers/rectilineargridwriter.impl.hh
index f2c822382e96f4c5d4f0c7e4d991fbf44146e483..d26428dc8724a8731d826f298626d480ccc3cc1f 100644
--- a/dune/vtk/writers/vtkrectilineargridwriter.impl.hh
+++ b/dune/vtk/writers/rectilineargridwriter.impl.hh
@@ -25,11 +25,11 @@ void VtkRectilinearGridWriter<GV,DC>
 
   auto const& wholeExtent = dataCollector_->wholeExtent();
   out << "<RectilinearGrid"
-      << " WholeExtent=\"" << join(wholeExtent.begin(), wholeExtent.end()) << "\""
+      << " WholeExtent=\"" << Vtk::join(wholeExtent.begin(), wholeExtent.end()) << "\""
       << ">\n";
 
   dataCollector_->writeLocalPiece([&out](auto const& extent) {
-    out << "<Piece Extent=\"" << join(extent.begin(), extent.end()) << "\">\n";
+    out << "<Piece Extent=\"" << Vtk::join(extent.begin(), extent.end()) << "\">\n";
   });
 
   // Write point coordinates for x, y, and z ordinate
@@ -66,7 +66,7 @@ void VtkRectilinearGridWriter<GV,DC>
   auto const& wholeExtent = dataCollector_->wholeExtent();
   out << "<PRectilinearGrid"
       << " GhostLevel=\"" << dataCollector_->ghostLevel() << "\""
-      << " WholeExtent=\"" << join(wholeExtent.begin(), wholeExtent.end()) << "\""
+      << " WholeExtent=\"" << Vtk::join(wholeExtent.begin(), wholeExtent.end()) << "\""
       << ">\n";
 
   // Write point coordinates for x, y, and z ordinate
@@ -104,7 +104,7 @@ void VtkRectilinearGridWriter<GV,DC>
     std::string piece_source = pfilename + "_p" + std::to_string(p) + "." + ext;
     out << "<Piece Source=\"" << piece_source << "\"";
     if (write_extent)
-      out << " Extent=\"" << join(extent.begin(), extent.end()) << "\"";
+      out << " Extent=\"" << Vtk::join(extent.begin(), extent.end()) << "\"";
      out << " />\n";
   });
 
diff --git a/dune/vtk/writers/vtkstructuredgridwriter.hh b/dune/vtk/writers/structuredgridwriter.hh
similarity index 91%
rename from dune/vtk/writers/vtkstructuredgridwriter.hh
rename to dune/vtk/writers/structuredgridwriter.hh
index 21d23134e8fb166a15ac0e6afdcd93261118cb99..3d07aec07a6584b85aab9eb6b7b1b4e2f3227580 100644
--- a/dune/vtk/writers/vtkstructuredgridwriter.hh
+++ b/dune/vtk/writers/structuredgridwriter.hh
@@ -6,11 +6,11 @@
 
 #include <dune/vtk/filewriter.hh>
 #include <dune/vtk/forward.hh>
-#include <dune/vtk/vtkfunction.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/function.hh>
+#include <dune/vtk/types.hh>
 #include <dune/vtk/datacollectors/structureddatacollector.hh>
 
-#include <dune/vtk/vtkwriterinterface.hh>
+#include <dune/vtk/writerinterface.hh>
 
 namespace Dune
 {
@@ -62,7 +62,7 @@ namespace Dune
   template <class GridView,
     class = std::void_t<typename GridView::IndexSet>>
   VtkStructuredGridWriter(GridView const&, Vtk::FormatTypes = Vtk::BINARY, Vtk::DataTypes = Vtk::FLOAT32)
-    -> VtkStructuredGridWriter<GridView, StructuredDataCollector<GridView>>;
+    -> VtkStructuredGridWriter<GridView, Vtk::StructuredDataCollector<GridView>>;
 
   template <class DataCollector,
     class = std::void_t<typename DataCollector::GridView>>
@@ -76,4 +76,4 @@ namespace Dune
 
 } // end namespace Dune
 
-#include "vtkstructuredgridwriter.impl.hh"
+#include "structuredgridwriter.impl.hh"
diff --git a/dune/vtk/writers/vtkstructuredgridwriter.impl.hh b/dune/vtk/writers/structuredgridwriter.impl.hh
similarity index 90%
rename from dune/vtk/writers/vtkstructuredgridwriter.impl.hh
rename to dune/vtk/writers/structuredgridwriter.impl.hh
index ec8707990614f24175c6248950d65d0d7148686e..0d09b88d5cead0fe046c7dd0134b61443af5b506 100644
--- a/dune/vtk/writers/vtkstructuredgridwriter.impl.hh
+++ b/dune/vtk/writers/structuredgridwriter.impl.hh
@@ -24,10 +24,10 @@ void VtkStructuredGridWriter<GV,DC>
   this->writeHeader(out, "StructuredGrid");
 
   auto const& wholeExtent = dataCollector_->wholeExtent();
-  out << "<StructuredGrid WholeExtent=\"" << join(wholeExtent.begin(), wholeExtent.end()) << "\">\n";
+  out << "<StructuredGrid WholeExtent=\"" << Vtk::join(wholeExtent.begin(), wholeExtent.end()) << "\">\n";
 
   dataCollector_->writeLocalPiece([&out](auto const& extent) {
-    out << "<Piece Extent=\"" << join(extent.begin(), extent.end()) << "\">\n";
+    out << "<Piece Extent=\"" << Vtk::join(extent.begin(), extent.end()) << "\">\n";
   });
 
   // Write point coordinates
@@ -64,7 +64,7 @@ void VtkStructuredGridWriter<GV,DC>
   auto const& wholeExtent = dataCollector_->wholeExtent();
   out << "<PStructuredGrid"
       << " GhostLevel=\"" << dataCollector_->ghostLevel() << "\""
-      << " WholeExtent=\"" << join(wholeExtent.begin(), wholeExtent.end()) << "\""
+      << " WholeExtent=\"" << Vtk::join(wholeExtent.begin(), wholeExtent.end()) << "\""
       << ">\n";
 
   // Write points
@@ -103,7 +103,7 @@ void VtkStructuredGridWriter<GV,DC>
     std::string piece_source = pfilename + "_p" + std::to_string(p) + "." + ext;
     out << "<Piece Source=\"" << piece_source << "\"";
     if (write_extent)
-      out << " Extent=\"" << join(extent.begin(), extent.end()) << "\"";
+      out << " Extent=\"" << Vtk::join(extent.begin(), extent.end()) << "\"";
     out << " />\n";
   });
 
diff --git a/dune/vtk/writers/vtkunstructuredgridwriter.hh b/dune/vtk/writers/unstructuredgridwriter.hh
similarity index 94%
rename from dune/vtk/writers/vtkunstructuredgridwriter.hh
rename to dune/vtk/writers/unstructuredgridwriter.hh
index a1d557ba9a35e04c43369c123c9e257d0a598482..39d670804691115f4da664546f1c00fa1aac8f02 100644
--- a/dune/vtk/writers/vtkunstructuredgridwriter.hh
+++ b/dune/vtk/writers/unstructuredgridwriter.hh
@@ -6,11 +6,11 @@
 
 #include <dune/vtk/filewriter.hh>
 #include <dune/vtk/forward.hh>
-#include <dune/vtk/vtkfunction.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/function.hh>
+#include <dune/vtk/types.hh>
 #include <dune/vtk/datacollectors/continuousdatacollector.hh>
 
-#include <dune/vtk/vtkwriterinterface.hh>
+#include <dune/vtk/writerinterface.hh>
 
 namespace Dune
 {
@@ -94,7 +94,7 @@ namespace Dune
   template <class GridView,
     class = std::void_t<typename GridView::IndexSet>>
   VtkUnstructuredGridWriter(GridView const&, Vtk::FormatTypes = Vtk::BINARY, Vtk::DataTypes = Vtk::FLOAT32)
-    -> VtkUnstructuredGridWriter<GridView, ContinuousDataCollector<GridView>>;
+    -> VtkUnstructuredGridWriter<GridView, Vtk::ContinuousDataCollector<GridView>>;
 
   template <class DataCollector,
     class = std::void_t<typename DataCollector::GridView>>
@@ -108,4 +108,4 @@ namespace Dune
 
 } // end namespace Dune
 
-#include "vtkunstructuredgridwriter.impl.hh"
+#include "unstructuredgridwriter.impl.hh"
diff --git a/dune/vtk/writers/vtkunstructuredgridwriter.impl.hh b/dune/vtk/writers/unstructuredgridwriter.impl.hh
similarity index 100%
rename from dune/vtk/writers/vtkunstructuredgridwriter.impl.hh
rename to dune/vtk/writers/unstructuredgridwriter.impl.hh
diff --git a/src/benchmark.cc b/src/benchmark.cc
index 32a3babc84375e83c3a749e3222aac40e8500bbc..1982745799104a1f02e6c84c998f191e1995c2a5 100644
--- a/src/benchmark.cc
+++ b/src/benchmark.cc
@@ -18,7 +18,7 @@
 #include <dune/grid/io/file/vtk.hh>
 #include <dune/grid/utility/structuredgridfactory.hh>
 
-#include <dune/vtk/writers/vtkunstructuredgridwriter.hh>
+#include <dune/vtk/writers/unstructuredgridwriter.hh>
 
 using namespace Dune;
 
diff --git a/src/datacollector.cc b/src/datacollector.cc
index 9f126d2c54c19e59ecc9050ff6f8e9c606cedd30..91ac8dea703308cf7e89f5ac544af31f8f8e468f 100644
--- a/src/datacollector.cc
+++ b/src/datacollector.cc
@@ -21,12 +21,12 @@
 
 #if HAVE_UG
 #include <dune/grid/uggrid.hh>
-#endif 
+#endif
 
 #include <dune/grid/yaspgrid.hh>
 #include <dune/grid/utility/structuredgridfactory.hh>
 
-#include <dune/vtk/writers/vtkunstructuredgridwriter.hh>
+#include <dune/vtk/writers/unstructuredgridwriter.hh>
 
 #include <dune/vtk/datacollectors/continuousdatacollector.hh>
 #include <dune/vtk/datacollectors/discontinuousdatacollector.hh>
@@ -73,12 +73,12 @@ void write (std::string prefix, GridView const& gridView)
   // write analytic function
   auto p1Analytic = makeAnalyticGridViewFunction([&c](auto const& x) { return c.dot(x); }, gridView);
 
-  write_dc<ContinuousDataCollector<GridView>>(prefix + "_continuous", gridView, p1Interpol, p1Analytic);
-  write_dc<DiscontinuousDataCollector<GridView>>(prefix + "_discontinuous", gridView, p1Interpol, p1Analytic);
-  write_dc<QuadraticDataCollector<GridView>>(prefix + "_quadratic", gridView, p1Interpol, p1Analytic);
+  write_dc<Vtk::ContinuousDataCollector<GridView>>(prefix + "_continuous", gridView, p1Interpol, p1Analytic);
+  write_dc<Vtk::DiscontinuousDataCollector<GridView>>(prefix + "_discontinuous", gridView, p1Interpol, p1Analytic);
+  write_dc<Vtk::QuadraticDataCollector<GridView>>(prefix + "_quadratic", gridView, p1Interpol, p1Analytic);
 
   Hybrid::forEach(StaticIntegralRange<int,7,1>{}, [&](auto p) {
-    write_dc<LagrangeDataCollector<GridView,p>>(prefix + "_lagrange_p" + std::to_string(p), gridView, p1Interpol, p1Analytic);
+    write_dc<Vtk::LagrangeDataCollector<GridView,p>>(prefix + "_lagrange_p" + std::to_string(p), gridView, p1Interpol, p1Analytic);
   });
 }
 
diff --git a/src/geometrygrid.cc b/src/geometrygrid.cc
index 462ac2e0a8c020543c649b97b608d02897a451cc..2c0c088ac1522211255b3f41a3990d74f7cba26d 100644
--- a/src/geometrygrid.cc
+++ b/src/geometrygrid.cc
@@ -15,7 +15,7 @@
 #include <dune/grid/yaspgrid.hh>
 #include <dune/grid/geometrygrid.hh>
 
-#include <dune/vtk/vtkwriter.hh>
+#include <dune/vtk/writer.hh>
 
 using namespace Dune;
 using namespace Dune::Functions;
diff --git a/src/lagrangepoints.cc b/src/lagrangepoints.cc
index a40f83e80c6e84e0f2f324d2d6becb7032d2b39b..9839ebb40dac1f51a53422d8ba049da5a62e6d3c 100644
--- a/src/lagrangepoints.cc
+++ b/src/lagrangepoints.cc
@@ -25,7 +25,7 @@ void write (std::string prefix, index_constant<dim>)
   for (int order = 1; order < 6; ++order) {
     std::cout << "order: " << order << std::endl;
     {
-      VtkLagrangePointSet<double, dim> pointSet(order);
+      Vtk::LagrangePointSet<double, dim> pointSet(order);
       pointSet.build(GeometryTypes::cube(dim));
 
       std::size_t i = 0;
@@ -34,9 +34,9 @@ void write (std::string prefix, index_constant<dim>)
         std::cout << i++ << ") p = " << p.point() << ", key = " << p.localKey() << std::endl;
       }
 
-      using BasisF = LagrangeBasisFactory<VtkLagrangePointSet, dim, double, double>;
-      using CoefficientF = LagrangeCoefficientsFactory<VtkLagrangePointSet, dim, double>;
-      using InterpolationF = LagrangeInterpolationFactory<VtkLagrangePointSet, dim, double>;
+      using BasisF = LagrangeBasisFactory<Vtk::LagrangePointSet, dim, double, double>;
+      using CoefficientF = LagrangeCoefficientsFactory<Vtk::LagrangePointSet, dim, double>;
+      using InterpolationF = LagrangeInterpolationFactory<Vtk::LagrangePointSet, dim, double>;
       GenericLocalFiniteElement<BasisF, CoefficientF, InterpolationF> localFE(GeometryTypes::cube(dim), order);
 
       auto const& localBasis = localFE.localBasis();
@@ -45,7 +45,7 @@ void write (std::string prefix, index_constant<dim>)
     }
 
     {
-      VtkLagrangePointSet<double, dim> pointSet(order);
+      Vtk::LagrangePointSet<double, dim> pointSet(order);
       pointSet.build(GeometryTypes::simplex(dim));
 
       std::size_t i = 0;
diff --git a/src/lagrangereader.cc b/src/lagrangereader.cc
index 8743c57cd881059de3b428221d5c68ed6591a0f5..c0eaadda47d1e35c3e6a9a8345232075b0a23247 100644
--- a/src/lagrangereader.cc
+++ b/src/lagrangereader.cc
@@ -19,10 +19,10 @@
 #include <dune/geometry/multilineargeometry.hh>
 #include <dune/grid/utility/structuredgridfactory.hh>
 
-#include <dune/vtk/vtkreader.hh>
+#include <dune/vtk/reader.hh>
 #include <dune/vtk/datacollectors/lagrangedatacollector.hh>
 #include <dune/vtk/gridcreators/lagrangegridcreator.hh>
-#include <dune/vtk/writers/vtkunstructuredgridwriter.hh>
+#include <dune/vtk/writers/unstructuredgridwriter.hh>
 
 using namespace Dune;
 
@@ -42,8 +42,8 @@ int main(int argc, char** argv)
   // using GridType = Dune::ALUGrid<dim,dow,Dune::simplex,Dune::conforming>;
   using GridType = FoamGrid<dim,dow>;
   using GridView = typename GridType::LeafGridView;
-  using DataCollector = LagrangeDataCollector<GridView, order>;
-  using GridCreator = LagrangeGridCreator<GridType>;
+  using DataCollector = Vtk::LagrangeDataCollector<GridView, order>;
+  using GridCreator = Vtk::LagrangeGridCreator<GridType>;
 
   std::string filename = "triangles_" + std::to_string(dow) + "d_order" + std::to_string(order);
 
diff --git a/src/legacyvtkwriter.cc b/src/legacyvtkwriter.cc
index a5838f409e43e23814f19fd48699b7712448a45b..6e7b837ec6cc68e90243915b87fe6dc66c2a00d3 100644
--- a/src/legacyvtkwriter.cc
+++ b/src/legacyvtkwriter.cc
@@ -15,7 +15,7 @@
 #include <dune/grid/uggrid.hh>
 #include <dune/grid/yaspgrid.hh>
 
-#include <dune/vtk/writers/vtkunstructuredgridwriter.hh>
+#include <dune/vtk/writers/unstructuredgridwriter.hh>
 #include <dune/vtk/legacyvtkfunction.hh>
 
 using namespace Dune;
diff --git a/src/polygongrid.cc b/src/polygongrid.cc
index 120a80f90e930c9f890b32089e68be35297e17ff..b430211f6ba88da6e5eb4fc5311a5c8c0e0a64a6 100644
--- a/src/polygongrid.cc
+++ b/src/polygongrid.cc
@@ -16,7 +16,7 @@
 #include <dune/polygongrid/gridfactory.hh>
 
 #include <dune/vtk/legacyvtkfunction.hh>
-#include <dune/vtk/writers/vtkunstructuredgridwriter.hh>
+#include <dune/vtk/writers/unstructuredgridwriter.hh>
 
 using namespace Dune;
 
diff --git a/src/pvdwriter.cc b/src/pvdwriter.cc
index e0e448038da5d2832a0a92b5630638965741e044..1a279779cb47821eb110da71e960a8771158f533 100644
--- a/src/pvdwriter.cc
+++ b/src/pvdwriter.cc
@@ -14,7 +14,7 @@
 #include <dune/functions/gridfunctions/analyticgridviewfunction.hh>
 #include <dune/grid/yaspgrid.hh>
 #include <dune/vtk/pvdwriter.hh>
-#include <dune/vtk/writers/vtkunstructuredgridwriter.hh>
+#include <dune/vtk/writers/unstructuredgridwriter.hh>
 
 using namespace Dune;
 using namespace Dune::Functions;
diff --git a/src/structuredgridwriter.cc b/src/structuredgridwriter.cc
index 2521f58992b285b7a82087be3b583e0990e71d73..33688d2b1bd792fa8466b2768a5d81da38494588 100644
--- a/src/structuredgridwriter.cc
+++ b/src/structuredgridwriter.cc
@@ -20,10 +20,10 @@
 #include <dune/grid/uggrid.hh>
 #include <dune/grid/yaspgrid.hh>
 
-#include <dune/vtk/writers/vtkimagedatawriter.hh>
-#include <dune/vtk/writers/vtkrectilineargridwriter.hh>
-#include <dune/vtk/writers/vtkstructuredgridwriter.hh>
-#include <dune/vtk/writers/vtkunstructuredgridwriter.hh>
+#include <dune/vtk/writers/imagedatawriter.hh>
+#include <dune/vtk/writers/rectilineargridwriter.hh>
+#include <dune/vtk/writers/structuredgridwriter.hh>
+#include <dune/vtk/writers/unstructuredgridwriter.hh>
 #include <dune/vtk/datacollectors/yaspdatacollector.hh>
 #include <dune/vtk/datacollectors/spdatacollector.hh>
 
diff --git a/src/test/mixed_element_test.cc b/src/test/mixed_element_test.cc
index 6dbd9e104eb2ad5b7d7d99c0e6d07a6e00adf77e..2bbb1a363605a7ffcd91659ee5e57b7c83d51357 100644
--- a/src/test/mixed_element_test.cc
+++ b/src/test/mixed_element_test.cc
@@ -20,8 +20,8 @@
 #include <dune/grid/yaspgrid.hh>
 #include <dune/grid/utility/structuredgridfactory.hh>
 
-#include <dune/vtk/vtkreader.hh>
-#include <dune/vtk/writers/vtkunstructuredgridwriter.hh>
+#include <dune/vtk/reader.hh>
+#include <dune/vtk/writers/unstructuredgridwriter.hh>
 
 using namespace Dune;
 
diff --git a/src/test/parallel_reader_writer_test.cc b/src/test/parallel_reader_writer_test.cc
index 5a8b400af13740b14cb8bf98bb8a0dd7d1aabb8a..c15ccf602b8d392e82a21298f48b31cc7c87b18e 100644
--- a/src/test/parallel_reader_writer_test.cc
+++ b/src/test/parallel_reader_writer_test.cc
@@ -25,8 +25,8 @@
 #include <dune/alugrid/grid.hh>
 #endif
 
-#include <dune/vtk/vtkreader.hh>
-#include <dune/vtk/writers/vtkunstructuredgridwriter.hh>
+#include <dune/vtk/reader.hh>
+#include <dune/vtk/writers/unstructuredgridwriter.hh>
 #include <dune/vtk/gridcreators/parallelgridcreator.hh>
 #include <dune/vtk/gridcreators/serialgridcreator.hh>
 
@@ -79,7 +79,7 @@ using HasParallelGridFactory = Std::is_detected<HasParallelGridFactoryImpl, Grid
 
 
 template <class Test>
-void compare (Test& test, filesystem::path const& dir, filesystem::path const& name)
+void compare (Test& test, Vtk::path const& dir, Vtk::path const& name)
 {
   test.check(compare_files(dir.string() + '/' + name.string() + ".vtu",
                            dir.string() + '/' + name.string() + "_2.vtu"));
@@ -182,19 +182,19 @@ int main (int argc, char** argv)
   TestSuite test{};
 
 #if HAVE_UG
-  reader_writer_test<UGGrid<2>, SerialGridCreator<UGGrid<2>>>(mpi, test, "UGGrid<2>");
-  reader_writer_test<UGGrid<3>, SerialGridCreator<UGGrid<3>>>(mpi, test, "UGGrid<3>");
+  reader_writer_test<UGGrid<2>, Vtk::SerialGridCreator<UGGrid<2>>>(mpi, test, "UGGrid<2>");
+  reader_writer_test<UGGrid<3>, Vtk::SerialGridCreator<UGGrid<3>>>(mpi, test, "UGGrid<3>");
 #endif
 
 #if HAVE_DUNE_ALUGRID
   // Test VtkWriter for ALUGrid.
-  reader_writer_test<ALUGridType<2>, SerialGridCreator<ALUGridType<2>>>(mpi, test, "ALUGridType<2>");
-  reader_writer_test<ALUGridType<2>, ParallelGridCreator<ALUGridType<2>>>(mpi, test, "ALUGridType<2, Parallel>", false);
+  reader_writer_test<ALUGridType<2>, Vtk::SerialGridCreator<ALUGridType<2>>>(mpi, test, "ALUGridType<2>");
+  reader_writer_test<ALUGridType<2>, Vtk::ParallelGridCreator<ALUGridType<2>>>(mpi, test, "ALUGridType<2, Parallel>", false);
 
-  reader_writer_test<ALUGridType<3>, SerialGridCreator<ALUGridType<3>>>(mpi, test, "ALUGridType<3>");
+  reader_writer_test<ALUGridType<3>, Vtk::SerialGridCreator<ALUGridType<3>>>(mpi, test, "ALUGridType<3>");
   #if DUNE_VERSION_LT(DUNE_GRID,2,7)
   // Currently the 2.7 branch is not working, due to a new bisection compatibility check in 3d
-  reader_writer_test<ALUGridType<3>, ParallelGridCreator<ALUGridType<3>>>(mpi, test, "ALUGridType<3, Parallel>", false);
+  reader_writer_test<ALUGridType<3>, Vtk::ParallelGridCreator<ALUGridType<3>>>(mpi, test, "ALUGridType<3, Parallel>", false);
   #endif
 #endif
 
diff --git a/src/test/reader_writer_test.cc b/src/test/reader_writer_test.cc
index 0043bfcddd0cf37c6875fde65bec3171beee5e5b..d9e6fbbfb410527f38c8ab47b55759c731b2d0d1 100644
--- a/src/test/reader_writer_test.cc
+++ b/src/test/reader_writer_test.cc
@@ -25,8 +25,8 @@
 #include <dune/alugrid/grid.hh>
 #endif
 
-#include <dune/vtk/vtkreader.hh>
-#include <dune/vtk/writers/vtkunstructuredgridwriter.hh>
+#include <dune/vtk/reader.hh>
+#include <dune/vtk/writers/unstructuredgridwriter.hh>
 #include <dune/vtk/gridcreators/continuousgridcreator.hh>
 
 using namespace Dune;
@@ -80,7 +80,7 @@ static TestCases test_cases = {
 };
 
 template <class Test>
-void compare (Test& test, filesystem::path const& dir, filesystem::path const& name)
+void compare (Test& test, Vtk::path const& dir, Vtk::path const& name)
 {
   test.check(compare_files(dir.string() + '/' + name.string() + ".vtu",
                            dir.string() + '/' + name.string() + "_2.vtu"));
diff --git a/src/timeserieswriter.cc b/src/timeserieswriter.cc
index e8ac785add8f210d0f31b5847aee33f6e70662c5..16fa630de41ee4d748eb4b82019fce0bd092322a 100644
--- a/src/timeserieswriter.cc
+++ b/src/timeserieswriter.cc
@@ -15,8 +15,8 @@
 #include <dune/functions/gridfunctions/analyticgridviewfunction.hh>
 #include <dune/grid/yaspgrid.hh>
 
-#include <dune/vtk/vtktimeserieswriter.hh>
-#include <dune/vtk/writers/vtkunstructuredgridwriter.hh>
+#include <dune/vtk/timeserieswriter.hh>
+#include <dune/vtk/writers/unstructuredgridwriter.hh>
 
 using namespace Dune;
 using namespace Dune::Functions;
diff --git a/src/vectorwriter.cc b/src/vectorwriter.cc
index 605328e18ca5b9570dacdcdf6f2ba2eb62bfb827..d6b29b14e2270d6b724829baeede71d63381f21e 100644
--- a/src/vectorwriter.cc
+++ b/src/vectorwriter.cc
@@ -14,7 +14,7 @@
 
 #include <dune/functions/gridfunctions/analyticgridviewfunction.hh>
 #include <dune/grid/yaspgrid.hh>
-#include <dune/vtk/vtkwriter.hh>
+#include <dune/vtk/writer.hh>
 
 using namespace Dune;
 using namespace Dune::Functions;
diff --git a/src/vtkreader.cc b/src/vtkreader.cc
index 13be5109dda8a51ca976f1a418273fec185ff793..52869739dd8d6c6baf49349bd6bf78a3c4087726 100644
--- a/src/vtkreader.cc
+++ b/src/vtkreader.cc
@@ -15,10 +15,10 @@
 #include <dune/grid/uggrid.hh>
 #include <dune/grid/utility/structuredgridfactory.hh>
 
-#include <dune/vtk/vtkreader.hh>
+#include <dune/vtk/reader.hh>
 #include <dune/vtk/gridcreators/continuousgridcreator.hh>
 #include <dune/vtk/gridcreators/discontinuousgridcreator.hh>
-#include <dune/vtk/writers/vtkunstructuredgridwriter.hh>
+#include <dune/vtk/writers/unstructuredgridwriter.hh>
 
 using namespace Dune;
 
@@ -78,14 +78,14 @@ int main(int argc, char** argv)
 
   {
     std::cout << "read 'test_ascii_float32.vtu'..." << std::endl;
-    auto gridPtr = VtkReader<GridType,DiscontinuousGridCreator<GridType>>::createGridFromFile("test_ascii_float32.vtu");
+    auto gridPtr = VtkReader<GridType,Vtk::DiscontinuousGridCreator<GridType>>::createGridFromFile("test_ascii_float32.vtu");
     auto& grid = *gridPtr;
 
     VtkUnstructuredGridWriter<GridView> vtkWriter(grid.leafGridView(), Vtk::ASCII);
     vtkWriter.write("test_ascii_float32_4.vtu");
   }
 
-  if (filesystem::exists("paraview_3d.vtu")) {
+  if (Vtk::exists("paraview_3d.vtu")) {
     using GridType3d = UGGrid<3>;
     using GridView3d = typename GridType3d::LeafGridView;
 
@@ -97,7 +97,7 @@ int main(int argc, char** argv)
     vtkWriter.write("paraview_3d_ascii.vtu");
   }
 
-  if (filesystem::exists("paraview_2d.vtu")) {
+  if (Vtk::exists("paraview_2d.vtu")) {
     std::cout << "paraview_2d_ascii...\n";
     using GridType2d = UGGrid<2>;
     using GridView2d = typename GridType2d::LeafGridView;
diff --git a/src/vtkwriter.cc b/src/vtkwriter.cc
index 435c3ff1a63bb463c45c9ddf70a3227a328b8809..9c9b2306eff87a839ad0e2d82cafba2bd4012d83 100644
--- a/src/vtkwriter.cc
+++ b/src/vtkwriter.cc
@@ -20,7 +20,7 @@
 #include <dune/grid/uggrid.hh>
 #include <dune/grid/yaspgrid.hh>
 
-#include <dune/vtk/vtkwriter.hh>
+#include <dune/vtk/writer.hh>
 
 using namespace Dune;
 using namespace Dune::Functions;