diff --git a/dumux/io/vtkfunction.hh b/dumux/io/vtkfunction.hh
index 6da7340f7d0cf273c6e75e3ed4be27f1179f8961..7a1c88044f8c28635d77e98ac6aa41531137a668 100644
--- a/dumux/io/vtkfunction.hh
+++ b/dumux/io/vtkfunction.hh
@@ -36,6 +36,24 @@
 
 namespace Dumux::Vtk {
 
+namespace Impl {
+
+template<class Field>
+double accessEntry(const Field& f, int mycomp, int i)
+{
+    if constexpr (IsIndexable<std::decay_t<decltype(std::declval<Field>()[0])>>{})
+    {
+        if constexpr (IsIndexable<std::decay_t<decltype(std::declval<Field>()[0][0])>>{})
+            DUNE_THROW(Dune::InvalidStateException, "Invalid field type");
+        else
+            return f[i][mycomp];
+    }
+    else
+        return f[i];
+}
+
+} // end namespace Impl
+
 /*!
  * \ingroup InputOutput
  * \brief a VTK function that supports both scalar and vector values for each element
@@ -61,7 +79,7 @@ public:
 
     //! evaluate
     virtual double evaluate(int mycomp, const Element& e, const Dune::FieldVector<ctype, dim>&) const
-    { return accessChooser_(mycomp, mapper_.index(e), IsIndexable<decltype(field_[0])>()); }
+    { return Impl::accessEntry(field_, mycomp, mapper_.index(e)); }
 
 #if DUNE_VERSION_GTE(DUNE_GRID, 2, 7)
     //! get output precision for the field
@@ -84,20 +102,6 @@ public:
     }
 
 private:
-
-    //! access for vectorial fields
-    double accessChooser_(int mycomp, int i, std::true_type) const
-    { return vectorFieldAccess_(mycomp, i, IsIndexable<decltype(field_[0][0])>()); }
-
-    //! access for scalar fields
-    double accessChooser_(int mycomp, int i, std::false_type) const { return field_[i]; }
-
-    //! access to permissive vectorial fields
-    double vectorFieldAccess_(int mycomp, int i, std::false_type) const { return field_[i][mycomp]; }
-
-    //! if the field is indexable more than two times, throw error
-    double vectorFieldAccess_(int mycomp, int i, std::true_type) const { DUNE_THROW(Dune::InvalidStateException, "Invalid field type"); }
-
     const F& field_;
     const std::string name_;
     int nComps_;
@@ -136,7 +140,7 @@ public:
 
         std::vector<Dune::FieldVector<ctype, 1>> cornerValues(nVertices);
         for (unsigned i = 0; i < nVertices; ++i)
-            cornerValues[i] = accessChooser_(mycomp, mapper_.subIndex(e, i, dim), IsIndexable<decltype(field_[0])>());
+            cornerValues[i] = Impl::accessEntry(field_, mycomp, mapper_.subIndex(e, i, dim));
 
         // (Ab)use the MultiLinearGeometry class to do multi-linear interpolation between scalars
         const Dune::MultiLinearGeometry<ctype, dim, 1> interpolation(e.type(), std::move(cornerValues));
@@ -164,20 +168,6 @@ public:
                                       << name << " (" << field.size() << ") and mapper (" << mapper.size() << ")");
     }
 private:
-
-    //! access for vectorial fields
-    double accessChooser_(int mycomp, int i, std::true_type) const
-    { return vectorFieldAccess_(mycomp, i, IsIndexable<decltype(field_[0][0])>()); }
-
-    //! access for scalar fields
-    double accessChooser_(int mycomp, int i, std::false_type) const { return field_[i]; }
-
-    //! access to permissive vectorial fields
-    double vectorFieldAccess_(int mycomp, int i, std::false_type) const { return field_[i][mycomp]; }
-
-    //! if the field is indexable more than two times, throw error
-    double vectorFieldAccess_(int mycomp, int i, std::true_type) const { DUNE_THROW(Dune::InvalidStateException, "Invalid field type"); }
-
     const F& field_;
     const std::string name_;
     int nComps_;
@@ -220,7 +210,7 @@ public:
 
         std::vector<Dune::FieldVector<ctype, 1>> cornerValues(nVertices);
         for (unsigned i = 0; i < nVertices; ++i)
-            cornerValues[i] = accessChooser_(mycomp, mapper_.index(e), i, IsIndexable<decltype(field_[0])>());
+            cornerValues[i] = accessEntry_(mycomp, mapper_.index(e), i);
 
         // (Ab)use the MultiLinearGeometry class to do multi-linear interpolation between scalars
         const Dune::MultiLinearGeometry<ctype, dim, 1> interpolation(e.type(), std::move(cornerValues));
@@ -247,22 +237,19 @@ public:
                                       << name << " (" << field.size() << ") and mapper (" << mapper.size() << ")");
     }
 private:
-
     //! access to the field
-    double accessChooser_(int mycomp, int eIdx, int cornerIdx, std::true_type) const
-    { return fieldAccess_(mycomp, eIdx, cornerIdx, IsIndexable<decltype(field_[0][0])>()); }
-
-    //! fields have to be indexable at least twice
-    double accessChooser_(int mycomp, int eIdx, int cornerIdx, std::false_type) const
-    { DUNE_THROW(Dune::InvalidStateException, "Invalid field type"); }
-
-    //! scalar field access
-    double fieldAccess_(int mycomp, int eIdx, int cornerIdx, std::false_type) const
-    { return field_[eIdx][cornerIdx]; }
-
-    //! vectorial field access
-    double fieldAccess_(int mycomp, int eIdx, int cornerIdx, std::true_type) const
-    { return field_[eIdx][cornerIdx][mycomp]; }
+    double accessEntry_(int mycomp, int eIdx, int cornerIdx) const
+    {
+        if constexpr (IsIndexable<decltype(std::declval<F>()[0])>{})
+        {
+            if constexpr (IsIndexable<decltype(std::declval<F>()[0][0])>{})
+                return field_[eIdx][cornerIdx][mycomp];
+            else
+                return field_[eIdx][cornerIdx];
+        }
+        else
+            DUNE_THROW(Dune::InvalidStateException, "Invalid field type");
+    }
 
     const F& field_;
     const std::string name_;