diff --git a/dune/vtk/CMakeLists.txt b/dune/vtk/CMakeLists.txt
index 8effca5d89517bd44962f9c7c980fecafab93ccb..a270981487bfaa32698b952402c10de2b0ea03c6 100644
--- a/dune/vtk/CMakeLists.txt
+++ b/dune/vtk/CMakeLists.txt
@@ -3,14 +3,20 @@ dune_add_library("vtktypes" OBJECT
 
 #install headers
 install(FILES
+  defaultvtkfunction.hh
   filereader.hh
   filewriter.hh
+  gridcreator.hh
+  legacyvtkfunction.hh
   vtkfunction.hh
+  vtklocalfunction.hh
+  vtklocalfunctioninterface.hh
   vtkreader.hh
   vtkreader.impl.hh
   vtktypes.hh
   vtkwriter.hh
-  vtkwriter.impl.hh
   DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dune/vtkwriter)
 
+add_subdirectory(datacollectors)
 add_subdirectory(utility)
+add_subdirectory(writers)
diff --git a/dune/vtk/datacollectors/CMakeLists.txt b/dune/vtk/datacollectors/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..c7841ffbe015595b08f0e71d17c847be2048a243
--- /dev/null
+++ b/dune/vtk/datacollectors/CMakeLists.txt
@@ -0,0 +1,11 @@
+#install headers
+install(FILES
+  continuousdatacollector.hh
+  datacollectorinterface.hh
+  discontinuousdatacollector.hh
+  quadraticdatacollector.hh
+  spdatacollector.hh
+  structureddatacollector.hh
+  unstructureddatacollector.hh
+  yaspdatacollector.hh
+  DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dune/vtkwriter/datacollectors)
diff --git a/dune/vtk/datacollectors/continuousdatacollector.hh b/dune/vtk/datacollectors/continuousdatacollector.hh
index 82b726645a06b86a65681f525aef33b3937e89c9..1c9d4dc26931367bebc9c679efbad291ba25506a 100644
--- a/dune/vtk/datacollectors/continuousdatacollector.hh
+++ b/dune/vtk/datacollectors/continuousdatacollector.hh
@@ -1,23 +1,22 @@
 #pragma once
 
-#include <dune/vtk/datacollector.hh>
+#include "unstructureddatacollector.hh"
 
 namespace Dune { namespace experimental
 {
 
 /// Implementation of \ref DataCollector for linear cells, with continuous data.
 template <class GridView>
-class DefaultDataCollector
-    : public DataCollectorInterface<GridView, DefaultDataCollector<GridView>>
+class ContinuousDataCollector
+    : public UnstructuredDataCollectorInterface<GridView, ContinuousDataCollector<GridView>>
 {
   enum { dim = GridView::dimension };
 
-  using Self = DefaultDataCollector;
-  using Super = DataCollectorInterface<GridView, Self>;
-  using Super::gridView_;
+  using Self = ContinuousDataCollector;
+  using Super = UnstructuredDataCollectorInterface<GridView, Self>;
 
 public:
-  DefaultDataCollector (GridView const& gridView)
+  ContinuousDataCollector (GridView const& gridView)
     : Super(gridView)
   {}
 
@@ -31,7 +30,7 @@ public:
   template <class T>
   std::vector<T> pointsImpl () const
   {
-    std::vector<T> data(this->numPoints() * 3);
+    std::vector<T> data(gridView_.size(dim) * 3);
     auto const& indexSet = gridView_.indexSet();
     for (auto const& vertex : vertices(gridView_, Partitions::all)) {
       std::size_t idx = 3 * indexSet.index(vertex);
@@ -44,6 +43,12 @@ public:
     return data;
   }
 
+  /// Return number of grid cells
+  std::uint64_t numCellsImpl () const
+  {
+    return gridView_.size(0);
+  }
+
   /// Return the types, offsets and connectivity of the cells, using the same connectivity as
   /// given by the grid.
   Cells cellsImpl () const
@@ -56,9 +61,9 @@ public:
     });
 
     Cells cells;
-    cells.connectivity.reserve(this->numCells() * maxVertices);
-    cells.offsets.reserve(this->numCells());
-    cells.types.reserve(this->numCells());
+    cells.connectivity.reserve(gridView_.size(0) * maxVertices);
+    cells.offsets.reserve(gridView_.size(0));
+    cells.types.reserve(gridView_.size(0));
 
     std::int64_t old_o = 0;
     for (auto const& c : elements(gridView_, Partitions::all)) {
@@ -68,7 +73,6 @@ public:
       cells.offsets.push_back(old_o += c.subEntities(dim));
       cells.types.push_back(cellType.type());
     }
-
     return cells;
   }
 
@@ -76,7 +80,7 @@ public:
   template <class T, class GlobalFunction>
   std::vector<T> pointDataImpl (GlobalFunction const& fct) const
   {
-    std::vector<T> data(this->numPoints() * fct.ncomps());
+    std::vector<T> data(gridView_.size(dim) * fct.ncomps());
     auto const& indexSet = gridView_.indexSet();
     auto localFct = localFunction(fct);
     for (auto const& e : elements(gridView_, Partitions::all)) {
@@ -92,6 +96,9 @@ public:
     }
     return data;
   }
+
+protected:
+  using Super::gridView_;
 };
 
 }} // end namespace Dune::experimental
diff --git a/dune/vtk/datacollector.hh b/dune/vtk/datacollectors/datacollectorinterface.hh
similarity index 72%
rename from dune/vtk/datacollector.hh
rename to dune/vtk/datacollectors/datacollectorinterface.hh
index 6c384f0a4e71e026dc27c0e0effb3f5b1fdede79..5d4f16d61e8b765ec4b416166f10d309f45800e0 100644
--- a/dune/vtk/datacollector.hh
+++ b/dune/vtk/datacollectors/datacollectorinterface.hh
@@ -1,41 +1,17 @@
 #pragma once
 
-#include <algorithm>
-#include <cstdint>
-#include <vector>
-
-#include "vtktypes.hh"
+#include <dune/vtk/vtktypes.hh>
 
 namespace Dune { namespace experimental {
 
-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 DataCollectorInterface
 {
 public:
-
   DataCollectorInterface (GridView const& gridView)
     : gridView_(gridView)
   {}
 
-  /// \brief Return the number of points in the grid
-  std::uint64_t numPoints () const
-  {
-    return asDerived().numPointsImpl();
-  }
-
-  /// \brief Return the number of cells in the grid
-  std::uint64_t numCells () const
-  {
-    return asDerived().numCellsImpl();
-  }
-
   /// Update the DataCollector on the current GridView
   void update ()
   {
@@ -48,6 +24,12 @@ public:
     return asDerived().ghostLevelImpl();
   }
 
+  /// \brief Return the number of points in 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.
@@ -61,12 +43,6 @@ public:
     return asDerived().template pointsImpl<T>();
   }
 
-  /// \brief Return cell types, offsets, and connectivity. \see Cells
-  Cells cells () const
-  {
-    return asDerived().cellsImpl();
-  }
-
   /// \brief Return a flat vector of function values evaluated at the grid points.
   /**
    * In case of a vector valued function, flat the vector entries:
@@ -76,20 +52,19 @@ public:
    * [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 GlobalFunction>
-  std::vector<T> pointData (GlobalFunction const& fct) const
+  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 grid cells. \see pointData.
-  template <class T, class GlobalFunction>
-  std::vector<T> cellData (GlobalFunction const& fct) const
+  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 ()
@@ -102,13 +77,7 @@ protected: // cast to derived type
     return static_cast<const Derived&>(*this);
   }
 
-
-protected: // default implementations
-
-  std::uint64_t numCellsImpl () const
-  {
-    return gridView_.size(0);
-  }
+public: // default implementations
 
   void updateImpl ()
   {
@@ -121,10 +90,10 @@ protected: // default implementations
   }
 
   // Evaluate `fct` in center of cell
-  template <class T, class GlobalFunction>
-  std::vector<T> cellDataImpl (GlobalFunction const& fct) const
+  template <class T, class VtkFunction>
+  std::vector<T> cellDataImpl (VtkFunction const& fct) const
   {
-    std::vector<T> data(numCells() * fct.ncomps());
+    std::vector<T> data(gridView_.size(0) * fct.ncomps());
     auto const& indexSet = gridView_.indexSet();
     auto localFct = localFunction(fct);
     for (auto const& e : elements(gridView_, Partitions::all)) {
@@ -138,9 +107,7 @@ protected: // default implementations
     return data;
   }
 
-
 protected:
-
   GridView gridView_;
 };
 
diff --git a/dune/vtk/datacollectors/discontinuousdatacollector.hh b/dune/vtk/datacollectors/discontinuousdatacollector.hh
index f922eb97569aedd969ea07bafa2a6a7be1e18a3c..1816a37854fce90986d4e98f309e9e37ce97dca1 100644
--- a/dune/vtk/datacollectors/discontinuousdatacollector.hh
+++ b/dune/vtk/datacollectors/discontinuousdatacollector.hh
@@ -1,6 +1,6 @@
 #pragma once
 
-#include <dune/vtk/datacollector.hh>
+#include "unstructureddatacollector.hh"
 
 namespace Dune { namespace experimental
 {
@@ -8,20 +8,17 @@ namespace Dune { namespace experimental
 /// Implementation of \ref DataCollector for linear cells, with discontinuous data.
 template <class GridView>
 class DiscontinuousDataCollector
-    : public DataCollectorInterface<GridView, DiscontinuousDataCollector<GridView>>
+    : public UnstructuredDataCollectorInterface<GridView, DiscontinuousDataCollector<GridView>>
 {
   enum { dim = GridView::dimension };
 
   using Self = DiscontinuousDataCollector;
-  using Super = DataCollectorInterface<GridView, Self>;
-  using Super::gridView_;
+  using Super = UnstructuredDataCollectorInterface<GridView, Self>;
 
 public:
   DiscontinuousDataCollector (GridView const& gridView)
     : Super(gridView)
-  {
-    this->update();
-  }
+  {}
 
   /// Create an index map the uniquely assignes an index to each pair (element,corner)
   void updateImpl ()
@@ -47,7 +44,7 @@ public:
   template <class T>
   std::vector<T> pointsImpl () const
   {
-    std::vector<T> data(this->numPoints() * 3);
+    std::vector<T> data(numPoints_ * 3);
     auto const& indexSet = gridView_.indexSet();
     for (auto const& element : elements(gridView_, Partitions::interior)) {
       for (unsigned int i = 0; i < element.subEntities(dim); ++i) {
@@ -62,13 +59,19 @@ public:
     return data;
   }
 
+  /// Return number of grid cells
+  std::uint64_t numCellsImpl () const
+  {
+    return gridView_.size(0);
+  }
+
   /// Connect the corners of each cell. The leads to a global discontinuous grid
   Cells cellsImpl () const
   {
     Cells cells;
-    cells.connectivity.reserve(this->numPoints());
-    cells.offsets.reserve(this->numCells());
-    cells.types.reserve(this->numCells());
+    cells.connectivity.reserve(numPoints_);
+    cells.offsets.reserve(gridView_.size(0));
+    cells.types.reserve(gridView_.size(0));
 
     std::int64_t old_o = 0;
     auto const& indexSet = gridView_.indexSet();
@@ -89,7 +92,7 @@ public:
   template <class T, class GlobalFunction>
   std::vector<T> pointDataImpl (GlobalFunction const& fct) const
   {
-    std::vector<T> data(this->numPoints() * fct.ncomps());
+    std::vector<T> data(numPoints_ * fct.ncomps());
     auto const& indexSet = gridView_.indexSet();
     auto localFct = localFunction(fct);
     for (auto const& e : elements(gridView_, Partitions::interior)) {
@@ -106,7 +109,8 @@ public:
     return data;
   }
 
-private:
+protected:
+  using Super::gridView_;
   std::uint64_t numPoints_ = 0;
   std::vector<std::int64_t> indexMap_;
 };
diff --git a/dune/vtk/datacollectors/quadraticdatacollector.hh b/dune/vtk/datacollectors/quadraticdatacollector.hh
index a39042e8701622ddb19224778222e21cfc983f5d..e9c8052121d82a91966a631bd10d60f807bc820c 100644
--- a/dune/vtk/datacollectors/quadraticdatacollector.hh
+++ b/dune/vtk/datacollectors/quadraticdatacollector.hh
@@ -1,6 +1,6 @@
 #pragma once
 
-#include <dune/vtk/datacollector.hh>
+#include "unstructureddatacollector.hh"
 
 namespace Dune { namespace experimental
 {
@@ -8,13 +8,12 @@ namespace Dune { namespace experimental
 /// Implementation of \ref DataCollector for quadratic cells, with continuous data.
 template <class GridView>
 class QuadraticDataCollector
-    : public DataCollectorInterface<GridView, QuadraticDataCollector<GridView>>
+    : public UnstructuredDataCollectorInterface<GridView, QuadraticDataCollector<GridView>>
 {
   enum { dim = GridView::dimension };
 
   using Self = QuadraticDataCollector;
-  using Super = DataCollectorInterface<GridView, Self>;
-  using Super::gridView_;
+  using Super = UnstructuredDataCollectorInterface<GridView, Self>;
 
 public:
   QuadraticDataCollector (GridView const& gridView)
@@ -63,6 +62,12 @@ public:
     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,
@@ -92,7 +97,6 @@ public:
       cells.offsets.push_back(old_o += c.subEntities(dim)+c.subEntities(dim-1));
       cells.types.push_back(cellType.type());
     }
-
     return cells;
   }
 
@@ -123,6 +127,9 @@ public:
     }
     return data;
   }
+
+protected:
+  using Super::gridView_;
 };
 
 }} // end namespace Dune::extensions
diff --git a/dune/vtk/datacollectors/spdatacollector.hh b/dune/vtk/datacollectors/spdatacollector.hh
index 60545ecac5489e7835970580825776f74f64739c..cb7e2b6544c81b7591d1fdbb94dc07c192791b45 100644
--- a/dune/vtk/datacollectors/spdatacollector.hh
+++ b/dune/vtk/datacollectors/spdatacollector.hh
@@ -3,7 +3,8 @@
 #if HAVE_DUNE_SPGRID
 #include <dune/grid/spgrid.hh>
 #endif
-#include <dune/vtk/datacollectors/structureddatacollector.hh>
+
+#include "structureddatacollector.hh"
 
 namespace Dune { namespace experimental
 {
@@ -18,7 +19,6 @@ class SPDataCollector
 
   using Self = SPDataCollector;
   using Super = StructuredDataCollectorInterface<GridView, Self>;
-  using Super::gridView_;
   using ctype = typename GridView::ctype;
 
 public:
@@ -69,7 +69,8 @@ public:
     }
   }
 
-private:
+protected:
+  using Super::gridView_;
   std::array<int, 6> wholeExtent_;
   std::array<int, 6> extent_;
   FieldVector<ctype,3> spacing_;
diff --git a/dune/vtk/datacollectors/structureddatacollector.hh b/dune/vtk/datacollectors/structureddatacollector.hh
index 0f9967abb44d9dce71ba074657367bdb8890e8f0..854aa032b50baa5d04bc7442d421c81f1e102eef 100644
--- a/dune/vtk/datacollectors/structureddatacollector.hh
+++ b/dune/vtk/datacollectors/structureddatacollector.hh
@@ -1,6 +1,9 @@
 #pragma once
 
-#include <dune/vtk/datacollectors/continuousdatacollector.hh>
+#include <array>
+#include <dune/common/fvector.hh>
+
+#include "continuousdatacollector.hh"
 
 namespace Dune { namespace experimental
 {
@@ -23,14 +26,13 @@ class StructuredDataCollectorInterface
 {
 protected:
   using Super = DataCollectorInterface<GridView, Derived>;
-  using Super::gridView_;
+  using SubDataCollector = ContinuousDataCollector<GridView>;
   using ctype = typename GridView::ctype;
 
 public:
   StructuredDataCollectorInterface (GridView const& gridView)
     : Super(gridView)
-    , defaultDataCollector_(gridView)
-    , ghostLevel_(gridView.overlapSize(0))
+    , subDataCollector_(gridView)
   {}
 
   /// Sequence of Index pairs [begin, end) for the cells in each direction
@@ -45,18 +47,6 @@ public:
     return this->asDerived().extentImpl();
   }
 
-  /// 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();
-  }
-
   /// Call the `writer` with extent
   template <class Writer>
   void writeLocalPiece (Writer const& writer) const
@@ -71,12 +61,26 @@ public:
     this->asDerived().writePiecesImpl(writer);
   }
 
-  /// Return the number of overlapping elements
-  int ghostLevel () const
+  /// Interface for ImageData:
+  /// @{
+
+  /// Lower left corner of the grid
+  FieldVector<ctype, 3> origin () const
   {
-    return this->asDerived().ghostLevelImpl();
+    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
@@ -84,17 +88,15 @@ public:
     return this->asDerived().template coordinatesImpl<T>();
   }
 
-public:
-  /// Return number of grid vertices
-  std::uint64_t numPointsImpl () const
-  {
-    return gridView_.size(GridView::dimension);
-  }
+  /// @}
+
+
+public: // default implementation:
 
   /// \copyref DefaultDataCollector::update.
   void updateImpl ()
   {
-    defaultDataCollector_.update();
+    subDataCollector_.update();
 
 #if HAVE_MPI
     int rank = -1;
@@ -113,23 +115,27 @@ public:
 #endif
   }
 
+  /// Return number of grid vertices
+  std::uint64_t numPointsImpl () const
+  {
+    return subDataCollector_.numPoints();
+  }
+
   /// \copydoc DefaultDataCollector::points.
   template <class T>
   std::vector<T> pointsImpl () const
   {
-    return defaultDataCollector_.template points<T>();
+    return subDataCollector_.template points<T>();
   }
 
   /// \copydoc DefaultDataCollector::pointData
   template <class T, class GlobalFunction>
   std::vector<T> pointDataImpl (GlobalFunction const& fct) const
   {
-    return defaultDataCollector_.template pointData<T>(fct);
+    return subDataCollector_.template pointData<T>(fct);
   }
 
-
-  /// Default implementation for \ref writeLocalPiece. Calculates the extent and communicates it to
-  /// rank 0.
+  // Calculates the extent and communicates it to rank 0.
   template <class Writer>
   void writeLocalPieceImpl (Writer const& writer) const
   {
@@ -154,8 +160,7 @@ public:
     writer(extent);
   }
 
-
-  /// Receive extent from all ranks and call the `writer` with the rank's extent vector
+  // Receive extent from all ranks and call the `writer` with the rank's extent vector
   template <class Writer>
   void writePiecesImpl (Writer const& writer) const
   {
@@ -176,15 +181,21 @@ public:
 #endif
   }
 
-
-  /// Return the \ref GridView::overlapSize
-  int ghostLevelImpl () const
+  // Origin (0,0,0)
+  FieldVector<ctype, 3> originImpl () const
   {
-    return ghostLevel_;
+    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
+  // Ordinate along each axis with constant \ref spacing from the \ref origin
   template <class T>
   std::array<std::vector<T>, 3> coordinatesImpl () const
   {
@@ -206,9 +217,9 @@ public:
     return ordinates;
   }
 
-private:
-  DefaultDataCollector<GridView> defaultDataCollector_;
-  int ghostLevel_;
+protected:
+  using Super::gridView_;
+  SubDataCollector subDataCollector_;
 
 #if HAVE_MPI
   mutable std::vector<std::array<int,6>> extents_;
diff --git a/dune/vtk/datacollectors/unstructureddatacollector.hh b/dune/vtk/datacollectors/unstructureddatacollector.hh
new file mode 100644
index 0000000000000000000000000000000000000000..b831102fd6cfbda8208404d5e008c474b19dca4a
--- /dev/null
+++ b/dune/vtk/datacollectors/unstructureddatacollector.hh
@@ -0,0 +1,44 @@
+#pragma once
+
+#include <cstdint>
+#include <vector>
+
+#include "datacollectorinterface.hh"
+
+namespace Dune { namespace experimental {
+
+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 UnstructuredDataCollectorInterface
+    : public DataCollectorInterface<GridView, Derived>
+{
+  using Super = DataCollectorInterface<GridView, Derived>;
+
+public:
+  UnstructuredDataCollectorInterface (GridView const& gridView)
+    : Super(gridView)
+  {}
+
+  /// \brief Return the number of cells in the grid
+  std::uint64_t numCells () const
+  {
+    return this->asDerived().numCellsImpl();
+  }
+
+  /// \brief Return cell types, offsets, and connectivity. \see Cells
+  Cells cells () const
+  {
+    return this->asDerived().cellsImpl();
+  }
+
+protected:
+  using Super::gridView_;
+};
+
+}} // end namespace Dune::experimental
diff --git a/dune/vtk/datacollectors/yaspdatacollector.hh b/dune/vtk/datacollectors/yaspdatacollector.hh
index 6902f61ed83fb2582a73569094d8d202f93cfef3..7dbbfd4c0082fafa1fae2307ba1dba73831aa55b 100644
--- a/dune/vtk/datacollectors/yaspdatacollector.hh
+++ b/dune/vtk/datacollectors/yaspdatacollector.hh
@@ -1,7 +1,8 @@
 #pragma once
 
 #include <dune/grid/yaspgrid.hh>
-#include <dune/vtk/datacollectors/structureddatacollector.hh>
+
+#include "structureddatacollector.hh"
 
 namespace Dune { namespace experimental
 {
@@ -14,7 +15,6 @@ class YaspDataCollector
 
   using Self = YaspDataCollector;
   using Super = StructuredDataCollectorInterface<GridView, Self>;
-  using Super::gridView_;
   using ctype = typename GridView::ctype;
 
 public:
@@ -114,7 +114,8 @@ public:
     return ordinates;
   }
 
-private:
+protected:
+  using Super::gridView_;
   std::array<int, 6> wholeExtent_;
   std::array<int, 6> extent_;
   FieldVector<ctype,3> origin_;
diff --git a/dune/vtk/vtkfunction.hh b/dune/vtk/vtkfunction.hh
index 08ea87c97a5be89a7f8c96cfff6359bd2ce891b8..70428b0e6f7d92c47165e5e7d2ad25333f91c13d 100644
--- a/dune/vtk/vtkfunction.hh
+++ b/dune/vtk/vtkfunction.hh
@@ -6,6 +6,7 @@
 #include <dune/functions/common/signature.hh>
 
 #include "vtklocalfunction.hh"
+#include "vtktypes.hh"
 
 namespace Dune { namespace experimental
 {
diff --git a/dune/vtk/vtktypes.cc b/dune/vtk/vtktypes.cc
index 4cfd7ac7ad79bc388169fa7ebca67d75fad26dd6..6e2f9e6528769851bdb20d4bbafefae069f662c7 100644
--- a/dune/vtk/vtktypes.cc
+++ b/dune/vtk/vtktypes.cc
@@ -134,7 +134,7 @@ CellType::CellType (GeometryType const& t, CellParametrization parametrization)
     }
     else if (t.isQuadrilateral()) {
       type_ = QUADRATIC_QUAD;
-      permutation_ = {0,1,3,2, 3,1,0,2}; // maybe need inverse permutation as well
+      permutation_ = {0,1,3,2, 2,1,3,0};
       noPermutation_ = false;
     }
     else if (t.isTetrahedron()) {
@@ -144,7 +144,7 @@ CellType::CellType (GeometryType const& t, CellParametrization parametrization)
     }
     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}; // maybe need inverse permutation as well
+      permutation_ = {0,1,3,2,4,5,7,6, 6,5,7,4,10,9,11,8,0,1,3,2};
       noPermutation_ = false;
     }
     else {
diff --git a/dune/vtk/vtkwriter.hh b/dune/vtk/vtkwriter.hh
index b5c94d3e8e411ac465d7b51782d9d6a2840c56df..6667f8ebcace1630fd54a6b3f6c8ed1513e290fd 100644
--- a/dune/vtk/vtkwriter.hh
+++ b/dune/vtk/vtkwriter.hh
@@ -61,6 +61,6 @@ namespace Dune { namespace experimental
 
   /// Default choice for several grid types, uses the default data-collector.
   template <class GridView>
-  using VtkWriter = Impl::VtkWriterImpl<GridView, typename GridView::Grid>;
+  using VtkWriter = typename Impl::VtkWriterImpl<GridView, typename GridView::Grid>::type;
 
 }} // end namespace Dune::experimental
diff --git a/dune/vtk/writers/CMakeLists.txt b/dune/vtk/writers/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..0a9e7e241d9e1d47def09960fa7da194dc3bffe9
--- /dev/null
+++ b/dune/vtk/writers/CMakeLists.txt
@@ -0,0 +1,13 @@
+#install headers
+install(FILES
+  vtkimagedatawriter.hh
+  vtkimagedatawriter.impl.hh
+  vtkrectilineargridwriter.hh
+  vtkrectilineargridwriter.impl.hh
+  vtkstructuredgridwriter.hh
+  vtkstructuredgridwriter.impl.hh
+  vtkunstructuredgridwriter.hh
+  vtkunstructuredgridwriter.impl.hh
+  vtkwriterinterface.hh
+  vtkwriterinterface.impl.hh
+  DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dune/vtkwriter/writers)
diff --git a/dune/vtk/writers/vtkimagedatawriter.hh b/dune/vtk/writers/vtkimagedatawriter.hh
index ada9fd941bbc167880421ec03d88c0ad25c137c8..34c9fa7e8e7961e1f642b1300a2433ec42c71425 100644
--- a/dune/vtk/writers/vtkimagedatawriter.hh
+++ b/dune/vtk/writers/vtkimagedatawriter.hh
@@ -4,7 +4,6 @@
 #include <iosfwd>
 #include <map>
 
-#include <dune/vtk/datacollector.hh>
 #include <dune/vtk/filewriter.hh>
 #include <dune/vtk/vtkfunction.hh>
 #include <dune/vtk/vtktypes.hh>
diff --git a/dune/vtk/writers/vtkrectilineargridwriter.hh b/dune/vtk/writers/vtkrectilineargridwriter.hh
index d92f41bf636a1801f17501a1c846f861a6b7576f..8ee4f80b03c3caad1ef01e83369c4eebae9a4d8f 100644
--- a/dune/vtk/writers/vtkrectilineargridwriter.hh
+++ b/dune/vtk/writers/vtkrectilineargridwriter.hh
@@ -4,7 +4,6 @@
 #include <iosfwd>
 #include <map>
 
-#include <dune/vtk/datacollector.hh>
 #include <dune/vtk/filewriter.hh>
 #include <dune/vtk/vtkfunction.hh>
 #include <dune/vtk/vtktypes.hh>
diff --git a/dune/vtk/writers/vtkstructuredgridwriter.hh b/dune/vtk/writers/vtkstructuredgridwriter.hh
index 0a0d47602a667676e07394242026f4ae95707057..e2772b84bd7e60286e9cdb99f7bf0ac0934e57a0 100644
--- a/dune/vtk/writers/vtkstructuredgridwriter.hh
+++ b/dune/vtk/writers/vtkstructuredgridwriter.hh
@@ -4,7 +4,6 @@
 #include <iosfwd>
 #include <map>
 
-#include <dune/vtk/datacollector.hh>
 #include <dune/vtk/filewriter.hh>
 #include <dune/vtk/vtkfunction.hh>
 #include <dune/vtk/vtktypes.hh>
diff --git a/dune/vtk/writers/vtkunstructuredgridwriter.hh b/dune/vtk/writers/vtkunstructuredgridwriter.hh
index 049a495b0d4322ec67d1cc72ab8607b3b5ddda26..2616e0176b6b9a73c92f342a88fc00cb4402b28b 100644
--- a/dune/vtk/writers/vtkunstructuredgridwriter.hh
+++ b/dune/vtk/writers/vtkunstructuredgridwriter.hh
@@ -4,7 +4,6 @@
 #include <iosfwd>
 #include <map>
 
-#include <dune/vtk/datacollector.hh>
 #include <dune/vtk/filewriter.hh>
 #include <dune/vtk/vtkfunction.hh>
 #include <dune/vtk/vtktypes.hh>
@@ -19,7 +18,7 @@ namespace Dune { namespace experimental
    * Requirement:
    * - DataCollector must be a model of \ref DataCollector
    **/
-  template <class GridView, class DataCollector = DefaultDataCollector<GridView>>
+  template <class GridView, class DataCollector = ContinuousDataCollector<GridView>>
   class VtkUnstructuredGridWriter
       : public VtkWriterInterface<GridView, DataCollector>
   {
diff --git a/dune/vtk/writers/vtkwriterinterface.hh b/dune/vtk/writers/vtkwriterinterface.hh
index 26e9202d75901cfa4d87649d3fa1f930c846de7a..aa657965455f0198bf28cde8ebb56b7a55188466 100644
--- a/dune/vtk/writers/vtkwriterinterface.hh
+++ b/dune/vtk/writers/vtkwriterinterface.hh
@@ -6,7 +6,6 @@
 
 #include <dune/common/std/optional.hh>
 
-#include <dune/vtk/datacollector.hh>
 #include <dune/vtk/filewriter.hh>
 #include <dune/vtk/vtkfunction.hh>
 #include <dune/vtk/vtktypes.hh>
diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt
index 96f90d21ea5c9a0f365b93ed92ba4baa3e59230b..0cd13b5fb8fcd2dfd14b95c49b742adbbf971854 100644
--- a/src/CMakeLists.txt
+++ b/src/CMakeLists.txt
@@ -14,10 +14,6 @@ add_executable("datacollector" datacollector.cc)
 target_link_dune_default_libraries("datacollector")
 target_link_libraries("datacollector" dunevtk)
 
-add_executable("quadratic" quadratic.cc)
-target_link_dune_default_libraries("quadratic")
-target_link_libraries("quadratic" dunevtk)
-
 add_executable("structuredgridwriter" structuredgridwriter.cc)
 target_link_dune_default_libraries("structuredgridwriter")
 target_link_libraries("structuredgridwriter" dunevtk)
diff --git a/src/datacollector.cc b/src/datacollector.cc
index 17009534101a06459a658718c9c9ba6ff892e005..fed2887dc80b2ef3af033203afbf6abf77e0bce8 100644
--- a/src/datacollector.cc
+++ b/src/datacollector.cc
@@ -17,100 +17,70 @@
 #include <dune/functions/functionspacebases/interpolate.hh>
 #include <dune/functions/gridfunctions/analyticgridviewfunction.hh>
 #include <dune/functions/gridfunctions/discreteglobalbasisfunction.hh>
-#include <dune/grid/uggrid.hh>
+
 #include <dune/grid/yaspgrid.hh>
 
 #include <dune/vtk/writers/vtkunstructuredgridwriter.hh>
-#include <dune/vtk/vtkreader.hh>
+
+#include <dune/vtk/datacollectors/continuousdatacollector.hh>
 #include <dune/vtk/datacollectors/discontinuousdatacollector.hh>
+#include <dune/vtk/datacollectors/quadraticdatacollector.hh>
 
 using namespace Dune;
 using namespace Dune::experimental;
 using namespace Dune::Functions;
 
-#define GRID_TYPE 2
-static const int dim = 3;
-#if GRID_TYPE == 1
-  using GridType = YaspGrid<dim>;
-#elif GRID_TYPE == 2
-  using GridType = UGGrid<dim>;
-#endif
-using GridViewType = typename GridType::LeafGridView;
-
-void write()
+template <class DataCollector, class GridView, class Fct1, class Fct2>
+void write_dc (std::string prefix, GridView const& gridView, Fct1 const& fct1, Fct2 const& fct2)
 {
-#if GRID_TYPE == 1
-  FieldVector<double,dim> upperRight; upperRight = 1.0;
-  auto numElements = filledArray<dim,int>(4);
-  GridType grid(upperRight,numElements);
-#elif GRID_TYPE == 2
-  FieldVector<double,dim> lowerLeft; lowerLeft = 0.0;
-  FieldVector<double,dim> upperRight; upperRight = 1.0;
-  auto numElements = filledArray<dim,unsigned int>(4);
-  auto gridPtr = StructuredGridFactory<GridType>::createSimplexGrid(lowerLeft, upperRight, numElements);
-  auto& grid = *gridPtr;
-#endif
-
-  GridViewType gridView = grid.leafGridView();
+  VtkUnstructuredGridWriter<GridView, DataCollector> vtkWriter(gridView);
+  vtkWriter.addPointData(fct1, "p1");
+  vtkWriter.addCellData(fct1, "p0");
+  vtkWriter.addPointData(fct2, "q1");
+  vtkWriter.addCellData(fct2, "q0");
+
+  vtkWriter.write(prefix + "_" + std::to_string(GridView::dimension) + "d_ascii.vtu",
+    Vtk::ASCII, Vtk::FLOAT32);
+}
 
+template <class GridView>
+void write (std::string prefix, GridView const& gridView)
+{
   using namespace BasisFactory;
   auto basis = makeBasis(gridView, lagrange<1>());
 
-  std::vector<double> p1function(basis.dimension());
+  FieldVector<double,GridView::dimension> c;
+  if (GridView::dimension > 0) c[0] = 11.0;
+  if (GridView::dimension > 1) c[1] = 7.0;
+  if (GridView::dimension > 2) c[2] = 3.0;
 
-  interpolate(basis, p1function, [](auto const& x) {
-    return 100*x[0] + 10*x[1] + 1*x[2];
-  });
+  std::vector<double> vec(basis.dimension());
+  interpolate(basis, vec, [&c](auto const& x) { return c.dot(x); });
 
   // write discrete global-basis function
-  auto p1FctWrapped = makeDiscreteGlobalBasisFunction<double>(basis, p1function);
+  auto p1Interpol = makeDiscreteGlobalBasisFunction<double>(basis, vec);
+
   // write analytic function
-  auto p1Analytic = makeAnalyticGridViewFunction([](auto const& x) {
-    return std::sin(10*x[0])*std::cos(10*x[1])+std::sin(10*x[2]);
-  }, gridView);
-
-
-  { // Default DataCollector
-    using Writer = VtkUnstructuredGridWriter<GridViewType>;
-    Writer vtkWriter(gridView);
-    vtkWriter.addPointData(p1FctWrapped, "p1");
-    vtkWriter.addCellData(p1FctWrapped, "p0");
-    vtkWriter.addPointData(p1Analytic, "analytic");
-
-    vtkWriter.write("c_ascii_float32.vtu", Vtk::ASCII);
-    vtkWriter.write("c_binary_float32.vtu", Vtk::BINARY);
-    vtkWriter.write("c_compressed_float32.vtu", Vtk::COMPRESSED);
-    vtkWriter.write("c_ascii_float64.vtu", Vtk::ASCII, Vtk::FLOAT64);
-    vtkWriter.write("c_binary_float64.vtu", Vtk::BINARY, Vtk::FLOAT64);
-    vtkWriter.write("c_compressed_float64.vtu", Vtk::COMPRESSED, Vtk::FLOAT64);
-  }
-
-  { // Discontinuous DataCollector
-    using Writer = VtkUnstructuredGridWriter<GridViewType, DiscontinuousDataCollector<GridViewType>>;
-    Writer vtkWriter(gridView);
-    vtkWriter.addPointData(p1FctWrapped, "p1");
-    vtkWriter.addCellData(p1FctWrapped, "p0");
-    vtkWriter.addPointData(p1Analytic, "analytic");
-
-    vtkWriter.write("dc_ascii_float32.vtu", Vtk::ASCII);
-    vtkWriter.write("dc_binary_float32.vtu", Vtk::BINARY);
-    vtkWriter.write("dc_compressed_float32.vtu", Vtk::COMPRESSED);
-    vtkWriter.write("dc_ascii_float64.vtu", Vtk::ASCII, Vtk::FLOAT64);
-    vtkWriter.write("dc_binary_float64.vtu", Vtk::BINARY, Vtk::FLOAT64);
-    vtkWriter.write("dc_compressed_float64.vtu", Vtk::COMPRESSED, Vtk::FLOAT64);
-  }
+  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);
 }
 
+template <int I>
+using int_ = std::integral_constant<int,I>;
+
 int main(int argc, char** argv)
 {
   Dune::MPIHelper::instance(argc, argv);
 
-  write();
+  Hybrid::forEach(std::make_tuple(int_<1>{}, int_<2>{}, int_<3>{}), [](auto dim)
   {
-    auto gridPtr = VtkReader<GridType, ConnectedGridCreator>::read("dc_ascii_float32.vtu");
-    auto& grid = *gridPtr;
-
-    VtkUnstructuredGridWriter<GridViewType> vtkWriter(grid.leafGridView());
-    vtkWriter.write("c_ascii_float32_2.vtu", Vtk::ASCII);
-  }
+    using GridType = YaspGrid<dim.value>;
+    FieldVector<double,dim.value> upperRight; upperRight = 1.0;
+    auto numElements = filledArray<dim.value,int>(4);
+    GridType grid(upperRight, numElements);
+    write("yasp", grid.leafGridView());
+  });
 }
diff --git a/src/legacyvtkwriter.cc b/src/legacyvtkwriter.cc
index 8631e0e74c21342ce9f2bf1576d94dacc8174b9f..52e273574c4a8ffcdae40ac419a9299c84240d83 100644
--- a/src/legacyvtkwriter.cc
+++ b/src/legacyvtkwriter.cc
@@ -27,27 +27,15 @@ using namespace Dune;
 using namespace Dune::experimental;
 using namespace Dune::Functions;
 
-#define GRID_TYPE 1
-
 int main(int argc, char** argv)
 {
   Dune::MPIHelper::instance(argc, argv);
 
   const int dim = 3;
-
-#if GRID_TYPE == 1
   using GridType = YaspGrid<dim>;
   FieldVector<double,dim> upperRight; upperRight = 1.0;
   auto numElements = filledArray<dim,int>(8);
   GridType grid(upperRight,numElements);
-#elif GRID_TYPE == 2
-  using GridType = UGGrid<dim>;
-  FieldVector<double,dim> lowerLeft; lowerLeft = 0.0;
-  FieldVector<double,dim> upperRight; upperRight = 1.0;
-  auto numElements = filledArray<dim,unsigned int>(4);
-  auto gridPtr = StructuredGridFactory<GridType>::createSimplexGrid(lowerLeft, upperRight, numElements);
-  auto& grid = *gridPtr;
-#endif
 
   using GridView = typename GridType::LeafGridView;
   GridView gridView = grid.leafGridView();
@@ -56,7 +44,6 @@ int main(int argc, char** argv)
   auto basis = makeBasis(gridView, lagrange<1>());
 
   std::vector<double> p1function(basis.dimension());
-
   interpolate(basis, p1function, [](auto const& x) {
     return 100*x[0] + 10*x[1] + 1*x[2];
   });
@@ -67,11 +54,5 @@ int main(int argc, char** argv)
   using Writer = VtkUnstructuredGridWriter<GridView>;
   Writer vtkWriter(gridView);
   vtkWriter.addPointData(p1FctWrapped);
-
   vtkWriter.write("test_ascii_float32.vtu", Vtk::ASCII);
-  vtkWriter.write("test_binary_float32.vtu", Vtk::BINARY);
-  vtkWriter.write("test_compressed_float32.vtu", Vtk::COMPRESSED);
-  vtkWriter.write("test_ascii_float64.vtu", Vtk::ASCII, Vtk::FLOAT64);
-  vtkWriter.write("test_binary_float64.vtu", Vtk::BINARY, Vtk::FLOAT64);
-  vtkWriter.write("test_compressed_float64.vtu", Vtk::COMPRESSED, Vtk::FLOAT64);
 }
diff --git a/src/polygongrid.cc b/src/polygongrid.cc
index 233cbae3fc121d70cfd3b3b567c4f16c71a3a486..dcc1ff973f04df0fa76efb63ff3de8f4805089ab 100644
--- a/src/polygongrid.cc
+++ b/src/polygongrid.cc
@@ -60,7 +60,7 @@ int main(int argc, char** argv)
   using Writer = VtkUnstructuredGridWriter<GridView>;
   Writer vtkWriter(gridView);
   auto p1Analytic = makeAnalyticGridViewFunction([](auto const& x) {
-    return std::sin(10*x[0])*std::cos(10*x[1])+std::sin(10*x[2]);
+    return std::sin(10*x[0])*std::cos(10*x[1]);
   }, gridView);
 
   vtkWriter.addPointData(p1Analytic, "p1");
diff --git a/src/quadratic.cc b/src/quadratic.cc
deleted file mode 100644
index 1385434b1e40de5608ca2e8a6a65f009105bcef8..0000000000000000000000000000000000000000
--- a/src/quadratic.cc
+++ /dev/null
@@ -1,85 +0,0 @@
-// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
-// vi: set et ts=4 sw=2 sts=2:
-
-#ifdef HAVE_CONFIG_H
-# include "config.h"
-#endif
-
-#include <iostream>
-#include <vector>
-
-#include <dune/common/parallel/mpihelper.hh> // An initializer of MPI
-#include <dune/common/exceptions.hh> // We use exceptions
-#include <dune/common/filledarray.hh>
-
-#include <dune/functions/functionspacebases/defaultglobalbasis.hh>
-#include <dune/functions/functionspacebases/lagrangebasis.hh>
-#include <dune/functions/functionspacebases/interpolate.hh>
-#include <dune/functions/gridfunctions/analyticgridviewfunction.hh>
-#include <dune/functions/gridfunctions/discreteglobalbasisfunction.hh>
-#include <dune/grid/uggrid.hh>
-#include <dune/grid/yaspgrid.hh>
-
-#include <dune/vtk/writers/vtkunstructuredgridwriter.hh>
-#include <dune/vtk/datacollectors/quadraticdatacollector.hh>
-
-using namespace Dune;
-using namespace Dune::experimental;
-using namespace Dune::Functions;
-
-#define GRID_TYPE 2
-
-int main(int argc, char** argv)
-{
-  Dune::MPIHelper::instance(argc, argv);
-
-  const int dim = 3;
-
-#if GRID_TYPE == 1
-  using GridType = YaspGrid<dim>;
-  FieldVector<double,dim> upperRight; upperRight = 1.0;
-  auto numElements = filledArray<dim,int>(4);
-  GridType grid(upperRight,numElements);
-#elif GRID_TYPE == 2
-  using GridType = UGGrid<dim>;
-  FieldVector<double,dim> lowerLeft; lowerLeft = 0.0;
-  FieldVector<double,dim> upperRight; upperRight = 1.0;
-  auto numElements = filledArray<dim,unsigned int>(4);
-  auto gridPtr = StructuredGridFactory<GridType>::createSimplexGrid(lowerLeft, upperRight, numElements);
-  auto& grid = *gridPtr;
-#endif
-
-  using GridView = typename GridType::LeafGridView;
-  GridView gridView = grid.leafGridView();
-
-  using namespace BasisFactory;
-  auto basis = makeBasis(gridView, lagrange<1>());
-
-  std::vector<double> p1function(basis.dimension());
-
-  interpolate(basis, p1function, [](auto const& x) {
-    return 100*x[0] + 10*x[1] + 1*x[2];
-  });
-
-  // write discrete global-basis function
-  auto p1FctWrapped = makeDiscreteGlobalBasisFunction<double>(basis, p1function);
-
-  using Writer = VtkUnstructuredGridWriter<GridView, QuadraticDataCollector<GridView>>;
-  Writer vtkWriter(gridView);
-  vtkWriter.addPointData(p1FctWrapped, "p1");
-  vtkWriter.addCellData(p1FctWrapped, "p0");
-
-  // write analytic function
-  auto p1Analytic = makeAnalyticGridViewFunction([](auto const& x) {
-    return std::sin(10*x[0])*std::cos(10*x[1])+std::sin(10*x[2]);
-  }, gridView);
-
-  vtkWriter.addPointData(p1Analytic, "analytic");
-
-  vtkWriter.write("p2_ascii_float32.vtu", Vtk::ASCII);
-  vtkWriter.write("p2_binary_float32.vtu", Vtk::BINARY);
-  vtkWriter.write("p2_compressed_float32.vtu", Vtk::COMPRESSED);
-  vtkWriter.write("p2_ascii_float64.vtu", Vtk::ASCII, Vtk::FLOAT64);
-  vtkWriter.write("p2_binary_float64.vtu", Vtk::BINARY, Vtk::FLOAT64);
-  vtkWriter.write("p2_compressed_float64.vtu", Vtk::COMPRESSED, Vtk::FLOAT64);
-}
diff --git a/src/structuredgridwriter.cc b/src/structuredgridwriter.cc
index 535ce1f68151bb61929bf3cbbf5e83e70adca2b8..a331c4b6c31aee8c2279dce0204ff3da0ea98379 100644
--- a/src/structuredgridwriter.cc
+++ b/src/structuredgridwriter.cc
@@ -37,9 +37,12 @@ using int_ = std::integral_constant<int,dim>;
 template <class GridView>
 void write(std::string prefix, GridView const& gridView)
 {
-  auto fct2 = makeAnalyticGridViewFunction([](auto const& x) -> float {
-    return std::sin(10*x[0]) * (x.size() > 1 ? std::cos(10*x[1]) : 1) + (x.size() > 2 ? std::sin(10*x[2]) : 0);
-  }, gridView);
+  FieldVector<double,GridView::dimension> c;
+  if (GridView::dimension > 0) c[0] = 11.0;
+  if (GridView::dimension > 1) c[1] = 7.0;
+  if (GridView::dimension > 2) c[2] = 3.0;
+
+  auto fct2 = makeAnalyticGridViewFunction([&c](auto const& x) -> float { return c.dot(x); }, gridView);
 
   {
     using Writer = VtkImageDataWriter<GridView>;
diff --git a/src/vtkwriter.cc b/src/vtkwriter.cc
index ffe6610e324ec4b242e52c904618fdd65d7d8ad3..001d0d3562847bb28b836dd86b02b43429392a14 100644
--- a/src/vtkwriter.cc
+++ b/src/vtkwriter.cc
@@ -20,65 +20,83 @@
 #include <dune/grid/uggrid.hh>
 #include <dune/grid/yaspgrid.hh>
 
-#include <dune/vtk/writers/vtkunstructuredgridwriter.hh>
+#include <dune/vtk/vtkwriter.hh>
 
 using namespace Dune;
 using namespace Dune::experimental;
 using namespace Dune::Functions;
 
-#define GRID_TYPE 1
-
-int main(int argc, char** argv)
+using TestCases = std::set<std::tuple<std::string,Vtk::FormatTypes,Vtk::DataTypes>>;
+static TestCases test_cases = {
+  {"ascii32", Vtk::ASCII, Vtk::FLOAT32},
+  {"bin32", Vtk::BINARY, Vtk::FLOAT32},
+  {"zlib32", Vtk::COMPRESSED, Vtk::FLOAT32},
+  {"ascii64", Vtk::ASCII, Vtk::FLOAT64},
+  {"bin64", Vtk::BINARY, Vtk::FLOAT64},
+  {"zlib64", Vtk::COMPRESSED, Vtk::FLOAT64},
+};
+
+template <class GridView>
+void write (std::string prefix, GridView const& gridView)
 {
-  Dune::MPIHelper::instance(argc, argv);
-
-  const int dim = 3;
-
-#if GRID_TYPE == 1
-  using GridType = YaspGrid<dim>;
-  FieldVector<double,dim> upperRight; upperRight = 1.0;
-  auto numElements = filledArray<dim,int>(8);
-  GridType grid(upperRight,numElements);
-#elif GRID_TYPE == 2
-  using GridType = UGGrid<dim>;
-  FieldVector<double,dim> lowerLeft; lowerLeft = 0.0;
-  FieldVector<double,dim> upperRight; upperRight = 1.0;
-  auto numElements = filledArray<dim,unsigned int>(4);
-  auto gridPtr = StructuredGridFactory<GridType>::createSimplexGrid(lowerLeft, upperRight, numElements);
-  auto& grid = *gridPtr;
-#endif
-
-  using GridView = typename GridType::LeafGridView;
-  GridView gridView = grid.leafGridView();
-
   using namespace BasisFactory;
   auto basis = makeBasis(gridView, lagrange<1>());
 
-  std::vector<double> p1function(basis.dimension());
+  FieldVector<double,GridView::dimension> c;
+  if (GridView::dimension > 0) c[0] = 11.0;
+  if (GridView::dimension > 1) c[1] = 7.0;
+  if (GridView::dimension > 2) c[2] = 3.0;
 
-  interpolate(basis, p1function, [](auto const& x) {
-    return 100*x[0] + 10*x[1] + 1*x[2];
-  });
+  std::vector<double> vec(basis.dimension());
+  interpolate(basis, vec, [&c](auto const& x) { return c.dot(x); });
 
   // write discrete global-basis function
-  auto p1FctWrapped = makeDiscreteGlobalBasisFunction<double>(basis, p1function);
-
-  using Writer = VtkUnstructuredGridWriter<GridView>;
-  Writer vtkWriter(gridView);
-  vtkWriter.addPointData(p1FctWrapped, "p1");
-  vtkWriter.addCellData(p1FctWrapped, "p0");
+  auto p1Interpol = makeDiscreteGlobalBasisFunction<double>(basis, vec);
 
   // write analytic function
-  auto p1Analytic = makeAnalyticGridViewFunction([](auto const& x) {
-    return std::sin(10*x[0])*std::cos(10*x[1])+std::sin(10*x[2]);
-  }, gridView);
-
-  vtkWriter.addPointData(p1Analytic, "analytic");
-
-  vtkWriter.write("test_ascii_float32.vtu", Vtk::ASCII);
-  vtkWriter.write("test_binary_float32.vtu", Vtk::BINARY);
-  vtkWriter.write("test_compressed_float32.vtu", Vtk::COMPRESSED);
-  vtkWriter.write("test_ascii_float64.vtu", Vtk::ASCII, Vtk::FLOAT64);
-  vtkWriter.write("test_binary_float64.vtu", Vtk::BINARY, Vtk::FLOAT64);
-  vtkWriter.write("test_compressed_float64.vtu", Vtk::COMPRESSED, Vtk::FLOAT64);
+  auto p1Analytic = makeAnalyticGridViewFunction([&c](auto const& x) { return c.dot(x); }, gridView);
+
+  VtkWriter<GridView> vtkWriter(gridView);
+  vtkWriter.addPointData(p1Interpol, "p1");
+  vtkWriter.addCellData(p1Interpol, "p0");
+  vtkWriter.addPointData(p1Analytic, "q1");
+  vtkWriter.addCellData(p1Analytic, "q0");
+  for (auto const& test_case : test_cases) {
+    vtkWriter.write(prefix + "_" + std::to_string(GridView::dimension) + "d_" + std::get<0>(test_case) + ".vtu",
+      std::get<1>(test_case), std::get<2>(test_case));
+  }
 }
+
+template <int I>
+using int_ = std::integral_constant<int,I>;
+
+int main (int argc, char** argv)
+{
+  Dune::MPIHelper::instance(argc, argv);
+
+#ifdef HAVE_UG
+  // Test VtkWriter for UGGrid
+  Hybrid::forEach(std::make_tuple(int_<2>{}, int_<3>{}), [](auto dim)
+  {
+    using GridType = UGGrid<dim.value>;
+    {
+      FieldVector<double,dim.value> lowerLeft; lowerLeft = 0.0;
+      FieldVector<double,dim.value> upperRight; upperRight = 1.0;
+      auto numElements = filledArray<dim.value,unsigned int>(4);
+      auto gridPtr = StructuredGridFactory<GridType>::createSimplexGrid(lowerLeft, upperRight, numElements);
+
+      write("ug", gridPtr->leafGridView());
+    }
+  });
+#endif
+
+  // Test VtkWriter for YaspGrid
+  Hybrid::forEach(std::make_tuple(int_<1>{}, int_<2>{}, int_<3>{}), [](auto dim)
+  {
+    using GridType = YaspGrid<dim.value>;
+    FieldVector<double,dim.value> upperRight; upperRight = 1.0;
+    auto numElements = filledArray<dim.value,int>(4);
+    GridType grid(upperRight, numElements);
+    write("yasp", grid.leafGridView());
+  });
+}
\ No newline at end of file