diff --git a/examples/ellipt.cc b/examples/ellipt.cc
index ddaad29b138ebdfd9a486c6f5f910851025bc160..6baf47bfb8cbba1adf2a4ffc10a5cb46f742a29f 100644
--- a/examples/ellipt.cc
+++ b/examples/ellipt.cc
@@ -32,6 +32,8 @@ int main(int argc, char** argv)
   auto opForce = makeOperator(tag::test{}, [](auto const& x) { return -1.0; }, 0);
   prob.addVectorOperator(opForce, _0);
 
+  auto opForce2 = makeOperator(tag::test{}, [](auto const& x) { return -2.0; }, 0);
+  prob.addVectorOperator(BoundaryType{0}, opForce2, _0);
 
   // set boundary condition
   auto predicate = [](auto const& x){ return x[0] < 1.e-8 || x[1] < 1.e-8; }; // define boundary
diff --git a/src/amdis/ContextGeometry.hpp b/src/amdis/ContextGeometry.hpp
index 0d5659bc0efad5389900ec6902f1d3212adce531..6307603e484c7d19668e360b1a659a6e3809cab0 100644
--- a/src/amdis/ContextGeometry.hpp
+++ b/src/amdis/ContextGeometry.hpp
@@ -151,7 +151,7 @@ namespace AMDiS
     Geometry const* geometry_;
 
     // The localGeometry may be constructed only if needed
-    Dune::Std::optional<LocalGeometry> localGeometry_;
+    mutable Dune::Std::optional<LocalGeometry> localGeometry_;
   };
 
 } // end namespace AMDiS
diff --git a/src/amdis/GridFunctionOperator.hpp b/src/amdis/GridFunctionOperator.hpp
index 61c7e019f1e91ca9c707bbb08110613bb9e58d94..1cda525f18a074f413632cf229f2fb20a88ec4e3 100644
--- a/src/amdis/GridFunctionOperator.hpp
+++ b/src/amdis/GridFunctionOperator.hpp
@@ -24,31 +24,29 @@ namespace AMDiS
    *
    * The class is specialized, by deriving from it, in \ref GridFunctionOperator.
    *
+   * \tparam LocalContext The Element or Intersection type
    * \tparam GridFunction The GridFunction, a LocalFunction is created from, and
    *                      that is evaluated at quadrature points.
-   * \tparam QuadratureCreator A functor that provides a \ref Dune::QuadratureRule.
    *
    * **Requirements:**
+   * - `LocalContext` models either Entity (of codim 0) or Intersection
    * - `GridFunction` models the \ref Concepts::GridFunction
-   * - `QuadratureCreator` models \ref Concepts::Callable<Dune::GeometryType, LocalFunction, F>
-   *   where `F` is a functor of the signature `int(int)` that calculates the
-   *   degree of the (bi)linear-form. The argument passed to `F` is the polynomial
-   *   order of the GridFunction.
    **/
-  template <class Derived, class GridFunction>
+  template <class Derived, class LocalContext, class GridFunction>
   class GridFunctionOperatorBase
-      : public LocalOperator<Derived>
+      : public LocalOperator<Derived, LocalContext>
   {
     using LocalFunction = decltype(localFunction(std::declval<GridFunction>()));
-    using Element = typename GridFunction::EntitySet::Element;
+
+    using Element = typename Impl::ContextTypes<LocalContext>::Entity;
     using Geometry = typename Element::Geometry;
 
     using LocalCoordinate = typename GridFunction::EntitySet::LocalCoordinate;
 
-    using QuadFactory = QuadratureFactory<typename Geometry::ctype, Element::dimension, LocalFunction>;
+    using QuadFactory = QuadratureFactory<typename Geometry::ctype, LocalContext::mydimension, LocalFunction>;
 
   public:
-    /// \brief Constructor. Stores a copy of `expr`.
+    /// \brief Constructor. Stores a copy of `gridFct`.
     /**
      * An expression operator takes an expression, following the interface of
      * \ref ExpressionBase, and stores a copy. Additionally, it gets the
@@ -85,10 +83,8 @@ namespace AMDiS
     void bind_impl(Element const& element, Geometry const& geometry)
     {
       assert( bool(quadFactory_) );
-      if (!this->bound_) {
-        localFct_.bind(element);
-        quadFactory_->bind(localFct_);
-      }
+      localFct_.bind(element);
+      quadFactory_->bind(localFct_);
     }
 
     /// Unbinds operator from element.
@@ -100,7 +96,7 @@ namespace AMDiS
     template <class PreQuadFactory>
     void setQuadFactory(PreQuadFactory const& pre)
     {
-      quadFactory_ = makeQuadratureFactoryPtr<typename Geometry::ctype, Element::dimension, LocalFunction>(pre);
+      quadFactory_ = makeQuadratureFactoryPtr<typename Geometry::ctype, LocalContext::mydimension, LocalFunction>(pre);
     }
 
   protected:
@@ -131,13 +127,16 @@ namespace AMDiS
     /// Assign each element type a quadrature rule
     std::shared_ptr<QuadFactory> quadFactory_ = std::make_shared<QuadFactory>();
 
-    int termOrder_ = 0;    //< the derivative order of this operator
+    /// the derivative order of this operator (in {0, 1, 2})
+    int termOrder_ = 0;
   };
 
 
+  /// \brief The transposed operator, implemented in term of its transposed by
+  /// calling \ref getElementMatrix with inverted arguments.
   template <class Derived, class Transposed>
   class GridFunctionOperatorTransposed
-      : public LocalOperator<Derived>
+      : public LocalOperator<Derived, typename Transposed::LocalContext>
   {
     template <class T, class... Args>
     using Constructable = decltype( new T(std::declval<Args>()...) );
@@ -152,12 +151,12 @@ namespace AMDiS
     template <class Element, class Geometry>
     void bind_impl(Element const& element, Geometry const& geometry)
     {
-      transposedOp_.bind(element, geometry);
+      transposedOp_.bind_impl(element, geometry);
     }
 
     void unbind_impl()
     {
-      transposedOp_.unbind();
+      transposedOp_.unbind_impl();
     }
 
     template <class PreQuadFactory>
@@ -181,17 +180,6 @@ namespace AMDiS
     Transposed transposedOp_;
   };
 
-
-  /// A default implementation of an GridFunctionOperator if no specialization is available.
-  /**
-   * An operator must implement either \ref getElementVector, or
-   * \ref getElementMatrix, if it is a vector or matrix operator,
-   * respectively.
-   **/
-  template <class Tag, class GridFct>
-  class GridFunctionOperator;
-
-
   template <class Tag, class PreGridFct, class PreQuadFactory>
   struct PreGridFunctionOperator
   {
@@ -200,7 +188,6 @@ namespace AMDiS
     PreQuadFactory quadFactory;
   };
 
-
   /// Store tag and expression to create a \ref GridFunctionOperator
   template <class Tag, class PreGridFct, class... QuadratureArgs>
   auto makeOperator(Tag tag, PreGridFct const& expr, QuadratureArgs&&... args)
@@ -213,24 +200,41 @@ namespace AMDiS
 
 #ifndef DOXYGEN
 
+  /// The base-template for GridFunctionOperators
+  /**
+   * An operator can specialize this class, by deriving from \ref GridFunctionOperatorBase.
+   * With the generic function \ref makeLocalOperator, an instance is created. To
+   * distinguisch different GridFunction operators, a tag can be provided that has no
+   * other effect.
+   *
+   * \tparam Tag  An Identifier for this GridFunctionOperator
+   * \tparam LocalContext  An Element or Intersection the operator is evaluated on
+   * \tparam GridFct  A GridFunction evaluated in local coordinates on the bound element
+   **/
+  template <class Tag, class LocalContext, class GridFct>
+  class GridFunctionOperator
+      : public GridFunctionOperatorBase<GridFunctionOperator<Tag, LocalContext, GridFct>, LocalContext, GridFct>
+  {};
+
+
   /// Generate an \ref GridFunctionOperator from a PreOperator (tag, expr).
   /// @{
-  template <class Tag, class... Args, class GridView>
+  template <class LocalContext, class Tag, class... Args, class GridView>
   auto makeLocalOperator(PreGridFunctionOperator<Tag, Args...> const& op, GridView const& gridView)
   {
     auto gridFct = makeGridFunction(op.expr, gridView);
-    using GridFctOp = GridFunctionOperator<Tag, decltype(gridFct)>;
+    using GridFctOp = GridFunctionOperator<Tag, LocalContext, decltype(gridFct)>;
     GridFctOp localOperator{op.tag, gridFct};
     localOperator.setQuadFactory(op.quadFactory);
     return localOperator;
   }
 
-  template <class Tag, class... Args, class GridView>
+  template <class LocalContext, class Tag, class... Args, class GridView>
   auto makeLocalOperator(std::reference_wrapper<PreGridFunctionOperator<Tag, Args...>> op, GridView const& gridView)
   {
     PreGridFunctionOperator<Tag, Args...> const& op_ref = op;
     auto gridFct = makeGridFunction(std::ref(op_ref.expr), gridView);
-    using GridFctOp = GridFunctionOperator<Tag, decltype(gridFct)>;
+    using GridFctOp = GridFunctionOperator<Tag, LocalContext, decltype(gridFct)>;
     GridFctOp localOperator{op_ref.tag, gridFct};
     localOperator.setQuadFactory(op_ref.quadFactory);
     return localOperator;
@@ -240,22 +244,22 @@ namespace AMDiS
 
   /// Generate a shared_ptr to \ref GridFunctionOperator from a PreOperator (tag, expr).
   /// @{
-  template <class Tag, class... Args, class GridView>
+  template <class LocalContext, class Tag, class... Args, class GridView>
   auto makeLocalOperatorPtr(PreGridFunctionOperator<Tag, Args...> const& op, GridView const& gridView)
   {
     auto gridFct = makeGridFunction(op.expr, gridView);
-    using GridFctOp = GridFunctionOperator<Tag, decltype(gridFct)>;
+    using GridFctOp = GridFunctionOperator<Tag, LocalContext, decltype(gridFct)>;
     auto localOperator = std::make_shared<GridFctOp>(op.tag, gridFct);
     localOperator->setQuadFactory(op.quadFactory);
     return localOperator;
   }
 
-  template <class Tag, class... Args, class GridView>
+  template <class LocalContext, class Tag, class... Args, class GridView>
   auto makeLocalOperatorPtr(std::reference_wrapper<PreGridFunctionOperator<Tag, Args...>> op, GridView const& gridView)
   {
     PreGridFunctionOperator<Tag, Args...> const& op_ref = op;
     auto gridFct = makeGridFunction(std::ref(op_ref.expr), gridView);
-    using GridFctOp = GridFunctionOperator<Tag, decltype(gridFct)>;
+    using GridFctOp = GridFunctionOperator<Tag, LocalContext, decltype(gridFct)>;
     auto localOperator = std::make_shared<GridFctOp>(op_ref.tag, gridFct);
     localOperator->setQuadFactory(op_ref.quadFactory);
     return localOperator;
diff --git a/src/amdis/LocalAssembler.hpp b/src/amdis/LocalAssembler.hpp
index 85a1db5a3a684e38036232ef576562082b5580f2..1d47532dfd273c401c1c5182483ba0e19ecf9ba2 100644
--- a/src/amdis/LocalAssembler.hpp
+++ b/src/amdis/LocalAssembler.hpp
@@ -65,14 +65,13 @@ namespace AMDiS
 
     /// Implementation of \ref LocalAssemblerBase::assemble
     /**
-     * Stores geometry and localGeometry and a \ref ContextGeometry and calls
+     * Stores geometry and localGeometry and calls
      * \ref calculateElementVector or \ref calculateElementMatrix on the
      * vector or matrix operator, respectively.
      **/
-    virtual void assemble(
-        LocalContext const& localContext,
-        Nodes const&... nodes,
-        ElementMatrixVector& elementMatrixVector) final
+    virtual void assemble(LocalContext const& localContext,
+                          Nodes const&... nodes,
+                          ElementMatrixVector& elementMatrixVector) final
     {
       ContextGeometry<LocalContext> context{localContext, getElement(), getGeometry()};
       assembleImpl(context, nodes..., elementMatrixVector);
@@ -122,7 +121,8 @@ namespace AMDiS
 
   private:
 
-    std::unique_ptr<Operator> storage_;     //< the stored operator, implementing \ref GridFunctionOperatorBase
+    /// the stored operator, implementing \ref GridFunctionOperatorBase
+    std::unique_ptr<Operator> storage_;
     Operator& op_;
 
     Element const* element_ = nullptr;
diff --git a/src/amdis/LocalOperator.hpp b/src/amdis/LocalOperator.hpp
index 5f4ac34cdbdb34b2dbd3249f4e0cc8858eb3b101..87d2f3d38bca5cf2c01fbaa2664b479c0bee0379 100644
--- a/src/amdis/LocalOperator.hpp
+++ b/src/amdis/LocalOperator.hpp
@@ -16,12 +16,18 @@ namespace AMDiS
 
   /// \brief The main implementation of an operator to be used in a \ref LocalAssembler.
   /**
-   * An Operator that takes the element it is assembled on.
+   * The CRTP Base class for local operators.
+   *
+   * \tparam Derived  The class that derives from this base class
+   * \tparam LocalContextType  The type of the element or intersection the operator is evaluated on
    **/
-  template <class Derived>
+  template <class Derived, class LocalContextType>
   class LocalOperator
   {
   public:
+    using LocalContext = LocalContextType;
+    using Element = typename Impl::ContextTypes<LocalContext>::Entity;
+    using Geometry = typename Element::Geometry;
 
     /// Initialize the local operator on the current gridView
     template <class GridView>
@@ -38,7 +44,6 @@ namespace AMDiS
      *
      * By default, it binds the \ref localFct_ to the `element`.
      **/
-    template <class Element, class Geometry>
     void bind(Element const& element, Geometry const& geometry)
     {
       if (!bound_) {
@@ -105,7 +110,7 @@ namespace AMDiS
     // default implementation. Can be overridden in the derived classes
     void unbind_impl() {}
 
-
+    // default implementation called by \ref calculateElementMatrix
     template <class Context, class RowNode, class ColNode, class ElementMatrix>
     void getElementMatrix(Context const& context,
                           RowNode const& rowNode,
@@ -115,8 +120,9 @@ namespace AMDiS
       error_exit("Needs to be implemented by derived class!");
     }
 
+    // default implementation called by \ref calculateElementVector
     template <class Context, class Node, class ElementVector>
-    void getElementVector(Context const& contextGeometry,
+    void getElementVector(Context const& context,
                           Node const& node,
                           ElementVector& elementVector)
     {
@@ -181,15 +187,15 @@ namespace AMDiS
 
 
   /// Generate an \ref GridFunctionOperator from a PreOperator.
-  template <class Derived, class GridView>
-  auto makeLocalOperator(LocalOperator<Derived> const& localOp, GridView const& /*gridView*/)
+  template <class Derived, class LocalContext, class GridView>
+  auto makeLocalOperator(LocalOperator<Derived, LocalContext> const& localOp, GridView const& /*gridView*/)
   {
     return localOp.derived();
   }
 
   /// Generate a shared_ptr to \ref GridFunctionOperator from a PreOperator.
-  template <class Derived, class GridView>
-  auto makeLocalOperatorPtr(LocalOperator<Derived> const& localOp, GridView const& /*gridView*/)
+  template <class Derived, class LocalContext, class GridView>
+  auto makeLocalOperatorPtr(LocalOperator<Derived, LocalContext> const& localOp, GridView const& /*gridView*/)
   {
     return std::make_shared<Derived>(localOp.derived());
   }
diff --git a/src/amdis/ProblemStat.inc.hpp b/src/amdis/ProblemStat.inc.hpp
index 37d72c9472bfcf4a7536291307b9f6ac7fae6a9e..4f08455f00d4d64a853b81d0f33fe518c3cbd767 100644
--- a/src/amdis/ProblemStat.inc.hpp
+++ b/src/amdis/ProblemStat.inc.hpp
@@ -176,7 +176,7 @@ void ProblemStat<Traits>::addMatrixOperator(
   auto i = child(localView_->tree(), makeTreePath(row));
   auto j = child(localView_->tree(), makeTreePath(col));
 
-  auto op = makeLocalOperator(preOp, globalBasis_->gridView());
+  auto op = makeLocalOperator<Element>(preOp, globalBasis_->gridView());
   auto localAssembler = makeLocalAssemblerPtr<Element>(std::move(op), i, j);
 
   matrixOperators_[i][j].element.push_back({localAssembler, factor, estFactor});
@@ -201,7 +201,7 @@ void ProblemStat<Traits>::addMatrixOperator(
   auto j = child(localView_->tree(), makeTreePath(col));
   using Intersection = typename GridView::Intersection;
 
-  auto op = makeLocalOperator(preOp, globalBasis_->gridView());
+  auto op = makeLocalOperator<Intersection>(preOp, globalBasis_->gridView());
   auto localAssembler = makeLocalAssemblerPtr<Intersection>(std::move(op), i, j);
 
   matrixOperators_[i][j].boundary.push_back({localAssembler, factor, estFactor, b});
@@ -221,7 +221,7 @@ void ProblemStat<Traits>::addVectorOperator(
 
   auto i = child(localView_->tree(), makeTreePath(path));
 
-  auto op = makeLocalOperator(preOp, globalBasis_->gridView());
+  auto op = makeLocalOperator<Element>(preOp, globalBasis_->gridView());
   auto localAssembler = makeLocalAssemblerPtr<Element>(std::move(op), i);
 
   rhsOperators_[i].element.push_back({localAssembler, factor, estFactor});
@@ -243,7 +243,7 @@ void ProblemStat<Traits>::addVectorOperator(
   auto i = child(localView_->tree(), makeTreePath(path));
   using Intersection = typename GridView::Intersection;
 
-  auto op = makeLocalOperator(preOp, globalBasis_->gridView());
+  auto op = makeLocalOperator<Intersection>(preOp, globalBasis_->gridView());
   auto localAssembler = makeLocalAssemblerPtr<Intersection>(std::move(op), i);
 
   rhsOperators_[i].boundary.push_back({localAssembler, factor, estFactor, b});
diff --git a/src/amdis/assembler/ConvectionDiffusionOperator.hpp b/src/amdis/assembler/ConvectionDiffusionOperator.hpp
index 78477656f758cac2bf71a5176af80dc14aca0fa3..a2f938e7a718a1c353d2052097b8beefb0f0c232 100644
--- a/src/amdis/assembler/ConvectionDiffusionOperator.hpp
+++ b/src/amdis/assembler/ConvectionDiffusionOperator.hpp
@@ -17,9 +17,10 @@ namespace AMDiS
   /// convection-diffusion operator
   /// <A*grad(u),grad(v)> - <b*u, grad(v)> + <c*u, v> = <f, v> (conserving) or
   /// <A*grad(u),grad(v)> + <b*grad(u), v> + <c*u, v> = <f, v> (non conserving)
-  template <class GridFctA, class GridFctB, class GridFctC, class GridFctF, bool conserving = true>
+  template <class LocalContext, class GridFctA, class GridFctB, class GridFctC, class GridFctF, bool conserving = true>
   class ConvectionDiffusionOperator
-      : public LocalOperatorBase<ConvectionDiffusionOperator<GridFctA, GridFctB, GridFctC, GridFctF, conserving>>
+      : public LocalOperator<ConvectionDiffusionOperator<LocalContext, GridFctA, GridFctB, GridFctC, GridFctF, conserving>,
+                             LocalContext>
   {
     using A_range_type = typename GridFctA::Range;
     static_assert( Category::Scalar<A_range_type> || Category::Matrix<A_range_type>,
@@ -36,8 +37,7 @@ namespace AMDiS
 
   public:
     ConvectionDiffusionOperator(GridFctA const& gridFctA, GridFctB const& gridFctB,
-                                GridFctC const& gridFctC, GridFctF const& gridFctF,
-                                bool_t<conserving> = {})
+                                GridFctC const& gridFctC, GridFctF const& gridFctF)
       : gridFctA_(gridFctA)
       , gridFctB_(gridFctB)
       , gridFctC_(gridFctC)
@@ -193,13 +193,37 @@ namespace AMDiS
     GridFctF gridFctF;
   };
 
+  template <class PreGridFctA, class PreGridFctB, class PreGridFctC, class PreGridFctF, bool conserving>
+  struct PreConvectionDiffusionOperator
+  {
+    PreGridFctA gridFctA;
+    PreGridFctB gridFctB;
+    PreGridFctC gridFctC;
+    PreGridFctF gridFctF;
+  };
 
-  template <class GridFctA, class GridFctB, class GridFctC, class GridFctF, bool conserving = true>
-  auto convectionDiffusion(GridFctA const& gridFctA, GridFctB const& gridFctB,
-                           GridFctC const& gridFctC, GridFctF const& gridFctF,
+  template <class PreGridFctA, class PreGridFctB, class PreGridFctC, class PreGridFctF, bool conserving = true>
+  auto convectionDiffusion(PreGridFctA const& gridFctA, PreGridFctB const& gridFctB,
+                           PreGridFctC const& gridFctC, PreGridFctF const& gridFctF,
                            bool_t<conserving> = {})
   {
-    return ConvectionDiffusionOperator<GridFctA, GridFctB, GridFctC, GridFctF, conserving>{gridFctA, gridFctB, gridFctC, gridFctF};
+    using Pre = PreConvectionDiffusionOperator<PreGridFctA, PreGridFctB, PreGridFctC, PreGridFctF, conserving>;
+    return Pre{gridFctA, gridFctB, gridFctC, gridFctF};
+  }
+
+  template <class LocalContext, class... GrdFcts, bool conserving, class GridView>
+  auto makeLocalOperator(PreConvectionDiffusionOperator<GridFcts..., conserving> const& pre, GridView const& gridView)
+  {
+    auto gridFctA = makeGridFunction(pre.gridFctA, gridView);
+    auto gridFctB = makeGridFunction(pre.gridFctB, gridView);
+    auto gridFctC = makeGridFunction(pre.gridFctC, gridView);
+    auto gridFctF = makeGridFunction(pre.gridFctF, gridView);
+
+    using GridFctOp = ConvectionDiffusionOperator<LocalContext,
+      decltype(gridFctA), decltype(gridFctB), decltype(gridFctC), decltype(gridFctF), conserving>;
+
+    GridFctOp localOperator{gridFctA, gridFctB, gridFctC, gridFctF};
+    return localOperator;
   }
 
   /** @} **/
diff --git a/src/amdis/assembler/FirstOrderDivTestvecTrial.hpp b/src/amdis/assembler/FirstOrderDivTestvecTrial.hpp
index f4524fd0457525b489e43b8045fc9bf33cd4cd60..dc49905dec8c8242c1a80952be2bb022275675b4 100644
--- a/src/amdis/assembler/FirstOrderDivTestvecTrial.hpp
+++ b/src/amdis/assembler/FirstOrderDivTestvecTrial.hpp
@@ -18,13 +18,13 @@ namespace AMDiS
 
 
   /// first-order operator \f$ \langle\nabla\cdot\Psi, c\,\phi\rangle \f$
-  template <class GridFct>
-  class GridFunctionOperator<tag::divtestvec_trial, GridFct>
-      : public GridFunctionOperatorTransposed<GridFunctionOperator<tag::divtestvec_trial, GridFct>,
-                                              GridFunctionOperator<tag::test_divtrialvec, GridFct>>
+  template <class LocalContext, class GridFct>
+  class GridFunctionOperator<tag::divtestvec_trial, LocalContext, GridFct>
+      : public GridFunctionOperatorTransposed<GridFunctionOperator<tag::divtestvec_trial, LocalContext, GridFct>,
+                                              GridFunctionOperator<tag::test_divtrialvec, LocalContext, GridFct>>
   {
     using Self = GridFunctionOperator;
-    using Transposed = GridFunctionOperator<tag::test_divtrialvec, GridFct>;
+    using Transposed = GridFunctionOperator<tag::test_divtrialvec, LocalContext, GridFct>;
     using Super = GridFunctionOperatorTransposed<Self, Transposed>;
 
   public:
diff --git a/src/amdis/assembler/FirstOrderGradTestTrial.hpp b/src/amdis/assembler/FirstOrderGradTestTrial.hpp
index eee21e9d8e3e5cc57e1d48ab5ca1629978b2f91c..6808103e4b6013135f7c278caed611dea2776dcc 100644
--- a/src/amdis/assembler/FirstOrderGradTestTrial.hpp
+++ b/src/amdis/assembler/FirstOrderGradTestTrial.hpp
@@ -18,13 +18,13 @@ namespace AMDiS
 
 
   /// first-order operator \f$ \langle\nabla\psi, \mathbf{b}\,\phi\rangle \f$
-  template <class GridFct>
-  class GridFunctionOperator<tag::gradtest_trial, GridFct>
-      : public GridFunctionOperatorTransposed<GridFunctionOperator<tag::gradtest_trial, GridFct>,
-                                              GridFunctionOperator<tag::test_gradtrial, GridFct>>
+  template <class LocalContext, class GridFct>
+  class GridFunctionOperator<tag::gradtest_trial, LocalContext, GridFct>
+      : public GridFunctionOperatorTransposed<GridFunctionOperator<tag::gradtest_trial, LocalContext, GridFct>,
+                                              GridFunctionOperator<tag::test_gradtrial, LocalContext, GridFct>>
   {
     using Self = GridFunctionOperator;
-    using Transposed = GridFunctionOperator<tag::test_gradtrial, GridFct>;
+    using Transposed = GridFunctionOperator<tag::test_gradtrial, LocalContext, GridFct>;
     using Super = GridFunctionOperatorTransposed<Self, Transposed>;
 
   public:
diff --git a/src/amdis/assembler/FirstOrderGradTestTrialvec.hpp b/src/amdis/assembler/FirstOrderGradTestTrialvec.hpp
index e1cf8d08a0f4c50c2f504a46bcb61840dd2e0db7..83c6c9e86c8eb9dceef7638ca56d3f612eae309e 100644
--- a/src/amdis/assembler/FirstOrderGradTestTrialvec.hpp
+++ b/src/amdis/assembler/FirstOrderGradTestTrialvec.hpp
@@ -18,13 +18,13 @@ namespace AMDiS
 
 
   /// first-order operator \f$ \langle\nabla\psi, c\,\Phi\rangle \f$
-  template <class GridFct>
-  class GridFunctionOperator<tag::gradtest_trialvec, GridFct>
-      : public GridFunctionOperatorTransposed<GridFunctionOperator<tag::gradtest_trialvec, GridFct>,
-                                              GridFunctionOperator<tag::testvec_gradtrial, GridFct>>
+  template <class LocalContext, class GridFct>
+  class GridFunctionOperator<tag::gradtest_trialvec, LocalContext, GridFct>
+      : public GridFunctionOperatorTransposed<GridFunctionOperator<tag::gradtest_trialvec, LocalContext, GridFct>,
+                                              GridFunctionOperator<tag::testvec_gradtrial, LocalContext, GridFct>>
   {
     using Self = GridFunctionOperator;
-    using Transposed = GridFunctionOperator<tag::testvec_gradtrial, GridFct>;
+    using Transposed = GridFunctionOperator<tag::testvec_gradtrial, LocalContext, GridFct>;
     using Super = GridFunctionOperatorTransposed<Self, Transposed>;
 
   public:
diff --git a/src/amdis/assembler/FirstOrderPartialTestTrial.hpp b/src/amdis/assembler/FirstOrderPartialTestTrial.hpp
index 14f9a892476197566a0aaf138648e0b19b32ae1a..0c9435f8e40c0d80328ab42a914fd1fe6f84273b 100644
--- a/src/amdis/assembler/FirstOrderPartialTestTrial.hpp
+++ b/src/amdis/assembler/FirstOrderPartialTestTrial.hpp
@@ -21,13 +21,13 @@ namespace AMDiS
 
 
   /// first-order operator \f$ \langle\partial_i\psi, c\,\phi\rangle \f$
-  template <class GridFct>
-  class GridFunctionOperator<tag::partialtest_trial, GridFct>
-      : public GridFunctionOperatorTransposed<GridFunctionOperator<tag::partialtest_trial, GridFct>,
-                                              GridFunctionOperator<tag::test_partialtrial, GridFct>>
+  template <class LocalContext, class GridFct>
+  class GridFunctionOperator<tag::partialtest_trial, LocalContext, GridFct>
+      : public GridFunctionOperatorTransposed<GridFunctionOperator<tag::partialtest_trial, LocalContext, GridFct>,
+                                              GridFunctionOperator<tag::test_partialtrial, LocalContext, GridFct>>
   {
     using Self = GridFunctionOperator;
-    using Transposed = GridFunctionOperator<tag::test_partialtrial, GridFct>;
+    using Transposed = GridFunctionOperator<tag::test_partialtrial, LocalContext, GridFct>;
     using Super = GridFunctionOperatorTransposed<Self, Transposed>;
 
   public:
diff --git a/src/amdis/assembler/FirstOrderTestDivTrialvec.hpp b/src/amdis/assembler/FirstOrderTestDivTrialvec.hpp
index 8c067d84ef35d96619a42d3ffca15717395ab73e..5cdbeab30c818ca90454eeffdba526943da80ad5 100644
--- a/src/amdis/assembler/FirstOrderTestDivTrialvec.hpp
+++ b/src/amdis/assembler/FirstOrderTestDivTrialvec.hpp
@@ -20,12 +20,13 @@ namespace AMDiS
 
 
   /// first-order operator \f$ \langle\psi, c\,\nabla\cdot\Phi\rangle \f$
-  template <class GridFct>
-  class GridFunctionOperator<tag::test_divtrialvec, GridFct>
-      : public GridFunctionOperatorBase<GridFunctionOperator<tag::test_divtrialvec, GridFct>, GridFct>
+  template <class LocalContext, class GridFct>
+  class GridFunctionOperator<tag::test_divtrialvec, LocalContext, GridFct>
+      : public GridFunctionOperatorBase<GridFunctionOperator<tag::test_divtrialvec, LocalContext, GridFct>,
+                                        LocalContext, GridFct>
   {
     using Self = GridFunctionOperator;
-    using Super = GridFunctionOperatorBase<Self, GridFct>;
+    using Super = GridFunctionOperatorBase<Self, LocalContext, GridFct>;
 
     static_assert( Category::Scalar<typename GridFct::Range>, "Expression must be of scalar type." );
 
diff --git a/src/amdis/assembler/FirstOrderTestGradTrial.hpp b/src/amdis/assembler/FirstOrderTestGradTrial.hpp
index 017f69341836fed3038dc8074d0d0d3869f117fa..bccf890d72ed0d0289b9df7755692603e6a1e0e3 100644
--- a/src/amdis/assembler/FirstOrderTestGradTrial.hpp
+++ b/src/amdis/assembler/FirstOrderTestGradTrial.hpp
@@ -20,12 +20,13 @@ namespace AMDiS
 
 
   /// first-order operator \f$ \langle\psi, \mathbf{b}\cdot\nabla\phi\rangle \f$
-  template <class GridFct>
-  class GridFunctionOperator<tag::test_gradtrial, GridFct>
-      : public GridFunctionOperatorBase<GridFunctionOperator<tag::test_gradtrial, GridFct>, GridFct>
+  template <class LocalContext, class GridFct>
+  class GridFunctionOperator<tag::test_gradtrial, LocalContext, GridFct>
+      : public GridFunctionOperatorBase<GridFunctionOperator<tag::test_gradtrial, LocalContext, GridFct>,
+                                        LocalContext, GridFct>
   {
     using Self = GridFunctionOperator;
-    using Super = GridFunctionOperatorBase<Self, GridFct>;
+    using Super = GridFunctionOperatorBase<Self, LocalContext, GridFct>;
 
     static_assert( Category::Vector<typename GridFct::Range>, "Expression must be of vector type." );
 
diff --git a/src/amdis/assembler/FirstOrderTestPartialTrial.hpp b/src/amdis/assembler/FirstOrderTestPartialTrial.hpp
index 6595956af6261ec4a5cacff34235ba507d643a36..d69dcb68b094c22ae90b48a6fc41590e93015f21 100644
--- a/src/amdis/assembler/FirstOrderTestPartialTrial.hpp
+++ b/src/amdis/assembler/FirstOrderTestPartialTrial.hpp
@@ -23,12 +23,13 @@ namespace AMDiS
 
 
   /// first-order operator \f$ \langle\psi, c\,\partial_i\phi\rangle \f$
-  template <class GridFct>
-  class GridFunctionOperator<tag::test_partialtrial, GridFct>
-      : public GridFunctionOperatorBase<GridFunctionOperator<tag::test_partialtrial, GridFct>, GridFct>
+  template <class LocalContext, class GridFct>
+  class GridFunctionOperator<tag::test_partialtrial, LocalContext, GridFct>
+      : public GridFunctionOperatorBase<GridFunctionOperator<tag::test_partialtrial, LocalContext, GridFct>,
+                                        LocalContext, GridFct>
   {
     using Self = GridFunctionOperator;
-    using Super = GridFunctionOperatorBase<Self, GridFct>;
+    using Super = GridFunctionOperatorBase<Self, LocalContext, GridFct>;
 
     static_assert( Category::Scalar<typename GridFct::Range>, "Expression must be of scalar type." );
 
diff --git a/src/amdis/assembler/FirstOrderTestvecGradTrial.hpp b/src/amdis/assembler/FirstOrderTestvecGradTrial.hpp
index 729edec4d6965d8458147e27627e435a6b4139d3..d06a44dbc79a2538e0933246edf81a30c1e80f98 100644
--- a/src/amdis/assembler/FirstOrderTestvecGradTrial.hpp
+++ b/src/amdis/assembler/FirstOrderTestvecGradTrial.hpp
@@ -20,12 +20,13 @@ namespace AMDiS
 
 
   /// first-order operator \f$ \langle\Psi, c\,\nabla\phi\rangle \f$
-  template <class GridFct>
-  class GridFunctionOperator<tag::testvec_gradtrial, GridFct>
-      : public GridFunctionOperatorBase<GridFunctionOperator<tag::testvec_gradtrial, GridFct>, GridFct>
+  template <class LocalContext, class GridFct>
+  class GridFunctionOperator<tag::testvec_gradtrial, LocalContext, GridFct>
+      : public GridFunctionOperatorBase<GridFunctionOperator<tag::testvec_gradtrial, LocalContext, GridFct>,
+                                        LocalContext, GridFct>
   {
     using Self = GridFunctionOperator;
-    using Super = GridFunctionOperatorBase<Self, GridFct>;
+    using Super = GridFunctionOperatorBase<Self, LocalContext, GridFct>;
 
     static_assert( Category::Scalar<typename GridFct::Range>, "Expression must be of scalar type." );
 
diff --git a/src/amdis/assembler/SecondOrderDivTestvecDivTrialvec.hpp b/src/amdis/assembler/SecondOrderDivTestvecDivTrialvec.hpp
index 350bc58fc39cb091cd9ad8cd28e893bda4ee318c..0e841d53cfef6760ae18741161538680b1922dcf 100644
--- a/src/amdis/assembler/SecondOrderDivTestvecDivTrialvec.hpp
+++ b/src/amdis/assembler/SecondOrderDivTestvecDivTrialvec.hpp
@@ -20,12 +20,13 @@ namespace AMDiS
 
 
   /// second-order operator \f$ \langle\nabla\cdot\Psi, c\,\nabla\cdot\Phi\rangle \f$
-  template <class GridFct>
-  class GridFunctionOperator<tag::divtestvec_divtrialvec, GridFct>
-      : public GridFunctionOperatorBase<GridFunctionOperator<tag::divtestvec_divtrialvec, GridFct>, GridFct>
+  template <class LocalContext, class GridFct>
+  class GridFunctionOperator<tag::divtestvec_divtrialvec, LocalContext, GridFct>
+      : public GridFunctionOperatorBase<GridFunctionOperator<tag::divtestvec_divtrialvec, LocalContext, GridFct>,
+                                        LocalContext, GridFct>
   {
     using Self = GridFunctionOperator;
-    using Super = GridFunctionOperatorBase<Self, GridFct>;
+    using Super = GridFunctionOperatorBase<Self, LocalContext, GridFct>;
 
     static_assert( Category::Scalar<typename GridFct::Range>, "Expression must be of scalar type." );
 
diff --git a/src/amdis/assembler/SecondOrderGradTestGradTrial.hpp b/src/amdis/assembler/SecondOrderGradTestGradTrial.hpp
index 30c46b1b55ef6704fee4e3629b838b9008b9d7f5..d8f6ea1ddda6a0418390ee8e290e1c8f56953066 100644
--- a/src/amdis/assembler/SecondOrderGradTestGradTrial.hpp
+++ b/src/amdis/assembler/SecondOrderGradTestGradTrial.hpp
@@ -20,12 +20,13 @@ namespace AMDiS
 
 
   /// second-order operator \f$ \langle\nabla\psi, c\,\nabla\phi\rangle \f$, or \f$ \langle\nabla\psi, A\,\nabla\phi\rangle \f$
-  template <class GridFct>
-  class GridFunctionOperator<tag::gradtest_gradtrial, GridFct>
-      : public GridFunctionOperatorBase<GridFunctionOperator<tag::gradtest_gradtrial, GridFct>, GridFct>
+  template <class LocalContext, class GridFct>
+  class GridFunctionOperator<tag::gradtest_gradtrial, LocalContext, GridFct>
+      : public GridFunctionOperatorBase<GridFunctionOperator<tag::gradtest_gradtrial, LocalContext, GridFct>,
+                                        LocalContext, GridFct>
   {
     using Self = GridFunctionOperator;
-    using Super = GridFunctionOperatorBase<Self, GridFct>;
+    using Super = GridFunctionOperatorBase<Self, LocalContext, GridFct>;
 
     using expr_value_type = typename GridFct::Range;
     static_assert( Category::Scalar<expr_value_type> || Category::Matrix<expr_value_type>,
diff --git a/src/amdis/assembler/SecondOrderPartialTestPartialTrial.hpp b/src/amdis/assembler/SecondOrderPartialTestPartialTrial.hpp
index bf964bbcbf86776d744be86887042af8c7ed6281..6539acf8ce1ec785c1a2cc689884f5e9fb8a2295 100644
--- a/src/amdis/assembler/SecondOrderPartialTestPartialTrial.hpp
+++ b/src/amdis/assembler/SecondOrderPartialTestPartialTrial.hpp
@@ -24,12 +24,13 @@ namespace AMDiS
 
 
   /// second-order operator \f$ \langle\partial_i\psi, c\,\partial_j\phi\rangle \f$
-  template <class GridFct>
-  class GridFunctionOperator<tag::partialtest_partialtrial, GridFct>
-      : public GridFunctionOperatorBase<GridFunctionOperator<tag::partialtest_partialtrial, GridFct>, GridFct>
+  template <class LocalContext, class GridFct>
+  class GridFunctionOperator<tag::partialtest_partialtrial, LocalContext, GridFct>
+      : public GridFunctionOperatorBase<GridFunctionOperator<tag::partialtest_partialtrial, LocalContext, GridFct>,
+                                        LocalContext, GridFct>
   {
     using Self = GridFunctionOperator;
-    using Super = GridFunctionOperatorBase<Self, GridFct>;
+    using Super = GridFunctionOperatorBase<Self, LocalContext, GridFct>;
 
     static_assert( Category::Scalar<typename GridFct::Range>, "Expression must be of scalar type." );
 
diff --git a/src/amdis/assembler/StokesOperator.hpp b/src/amdis/assembler/StokesOperator.hpp
index 8b6aa6a5b751dbaa6829b6c21c05d5ac4391846e..b9dd503273666303b696d2b6a01effdcf4220366 100644
--- a/src/amdis/assembler/StokesOperator.hpp
+++ b/src/amdis/assembler/StokesOperator.hpp
@@ -20,12 +20,13 @@ namespace AMDiS
 
 
   /// Stokes operator \f$ \langle\nabla\mathbf{v}, \nu \nabla\mathbf{u}\rangle + \langle\nabla\cdot\mathbf{v}, p\rangle + \langle q, \langle\nabla\cdot\mathbf{u}\rangle \f$
-  template <class ViscosityExpr>
-  class GridFunctionOperator<tag::stokes, ViscosityExpr>
-      : public GridFunctionOperatorBase<GridFunctionOperator<tag::stokes, ViscosityExpr>, ViscosityExpr>
+  template <class LocalContext, class ViscosityExpr>
+  class GridFunctionOperator<tag::stokes, LocalContext, ViscosityExpr>
+      : public GridFunctionOperatorBase<GridFunctionOperator<tag::stokes, LocalContext, ViscosityExpr>,
+                                        LocalContext, ViscosityExpr>
   {
     using Self = GridFunctionOperator;
-    using Super = GridFunctionOperatorBase<Self, ViscosityExpr>;
+    using Super = GridFunctionOperatorBase<Self, LocalContext, ViscosityExpr>;
 
     static_assert( Category::Scalar<typename ViscosityExpr::Range>, "Viscosity must be of scalar type." );
 
diff --git a/src/amdis/assembler/ZeroOrderTest.hpp b/src/amdis/assembler/ZeroOrderTest.hpp
index 39dbdb1f7d62c0e51fcd9b23dad923ba903c7197..5639f0c57c26e05bf7f862c823ba9d9d29bc22ef 100644
--- a/src/amdis/assembler/ZeroOrderTest.hpp
+++ b/src/amdis/assembler/ZeroOrderTest.hpp
@@ -19,12 +19,13 @@ namespace AMDiS
 
 
   /// zero-order vector-operator \f$ (c\, \psi) \f$
-  template <class GridFct>
-  class GridFunctionOperator<tag::test, GridFct>
-      : public GridFunctionOperatorBase<GridFunctionOperator<tag::test, GridFct>, GridFct>
+  template <class LocalContext, class GridFct>
+  class GridFunctionOperator<tag::test, LocalContext, GridFct>
+      : public GridFunctionOperatorBase<GridFunctionOperator<tag::test, LocalContext, GridFct>,
+                                        LocalContext, GridFct>
   {
     using Self = GridFunctionOperator;
-    using Super = GridFunctionOperatorBase<Self, GridFct>;
+    using Super = GridFunctionOperatorBase<Self, LocalContext, GridFct>;
 
     static_assert( Category::Scalar<typename GridFct::Range>, "Expression must be of scalar type." );
 
diff --git a/src/amdis/assembler/ZeroOrderTestTrial.hpp b/src/amdis/assembler/ZeroOrderTestTrial.hpp
index 2315ea92b8d0fcc886e68f768f0ccda719e36c9f..47b0e4d3676e9165604c1c7ef475f07d5f86e52c 100644
--- a/src/amdis/assembler/ZeroOrderTestTrial.hpp
+++ b/src/amdis/assembler/ZeroOrderTestTrial.hpp
@@ -19,12 +19,13 @@ namespace AMDiS
 
 
   /// zero-order operator \f$ \langle\psi, c\,\phi\rangle \f$
-  template <class GridFct>
-  class GridFunctionOperator<tag::test_trial, GridFct>
-      : public GridFunctionOperatorBase<GridFunctionOperator<tag::test_trial, GridFct>, GridFct>
+  template <class LocalContext, class GridFct>
+  class GridFunctionOperator<tag::test_trial, LocalContext, GridFct>
+      : public GridFunctionOperatorBase<GridFunctionOperator<tag::test_trial, LocalContext, GridFct>,
+                                        LocalContext, GridFct>
   {
     using Self = GridFunctionOperator;
-    using Super = GridFunctionOperatorBase<Self, GridFct>;
+    using Super = GridFunctionOperatorBase<Self, LocalContext, GridFct>;
 
     static_assert( Category::Scalar<typename GridFct::Range>, "Expression must be of scalar type." );
 
diff --git a/src/amdis/assembler/ZeroOrderTestTrialvec.hpp b/src/amdis/assembler/ZeroOrderTestTrialvec.hpp
index a4dbec12faf19910b049af8a46c3aa802e44511e..352cacfc0b58455043df285da8222924caba2127 100644
--- a/src/amdis/assembler/ZeroOrderTestTrialvec.hpp
+++ b/src/amdis/assembler/ZeroOrderTestTrialvec.hpp
@@ -19,12 +19,13 @@ namespace AMDiS
 
 
   /// zero-order operator \f$ \langle\psi, \mathbf{b}\cdot\Phi\rangle \f$
-  template <class GridFct>
-  class GridFunctionOperator<tag::test_trialvec, GridFct>
-      : public GridFunctionOperatorBase<GridFunctionOperator<tag::test_trialvec, GridFct>, GridFct>
+  template <class LocalContext, class GridFct>
+  class GridFunctionOperator<tag::test_trialvec, LocalContext, GridFct>
+      : public GridFunctionOperatorBase<GridFunctionOperator<tag::test_trialvec, LocalContext, GridFct>,
+                                        LocalContext, GridFct>
   {
     using Self = GridFunctionOperator;
-    using Super = GridFunctionOperatorBase<Self, GridFct>;
+    using Super = GridFunctionOperatorBase<Self, LocalContext, GridFct>;
 
     static_assert( Category::Vector<typename GridFct::Range>, "Expression must be of vector type." );
 
diff --git a/src/amdis/assembler/ZeroOrderTestvec.hpp b/src/amdis/assembler/ZeroOrderTestvec.hpp
index 68a3fe3f9922ef51fce5e2250633d3952658bc7b..0a6c876c5d5f5a8fc6c26654e4e461fe50cb1ddd 100644
--- a/src/amdis/assembler/ZeroOrderTestvec.hpp
+++ b/src/amdis/assembler/ZeroOrderTestvec.hpp
@@ -19,12 +19,13 @@ namespace AMDiS
 
 
   /// zero-order vector-operator \f$ (\mathbf{b}\cdot\Psi) \f$
-  template <class GridFct>
-  class GridFunctionOperator<tag::testvec, GridFct>
-      : public GridFunctionOperatorBase<GridFunctionOperator<tag::testvec, GridFct>, GridFct>
+  template <class LocalContext, class GridFct>
+  class GridFunctionOperator<tag::testvec, LocalContext, GridFct>
+      : public GridFunctionOperatorBase<GridFunctionOperator<tag::testvec, LocalContext, GridFct>,
+                                        LocalContext, GridFct>
   {
     using Self = GridFunctionOperator;
-    using Super = GridFunctionOperatorBase<Self, GridFct>;
+    using Super = GridFunctionOperatorBase<Self, LocalContext, GridFct>;
 
     static_assert( Category::Vector<typename GridFct::Range>, "Expression must be of vector type." );
 
diff --git a/src/amdis/assembler/ZeroOrderTestvecTrial.hpp b/src/amdis/assembler/ZeroOrderTestvecTrial.hpp
index ad67fa578977a89dd2b5a4a60977feb6526eab4f..f37c6aa46cc7bd46e7de3e1983c46b0c22c0b1b4 100644
--- a/src/amdis/assembler/ZeroOrderTestvecTrial.hpp
+++ b/src/amdis/assembler/ZeroOrderTestvecTrial.hpp
@@ -18,13 +18,13 @@ namespace AMDiS
 
 
   /// zero-order operator \f$ \langle\Psi, \mathbf{b}\,\phi\rangle \f$
-  template <class GridFct>
-  class GridFunctionOperator<tag::testvec_trial, GridFct>
-      : public GridFunctionOperatorTransposed<GridFunctionOperator<tag::testvec_trial, GridFct>,
-                                              GridFunctionOperator<tag::test_trialvec, GridFct>>
+  template <class LocalContext, class GridFct>
+  class GridFunctionOperator<tag::testvec_trial, LocalContext, GridFct>
+      : public GridFunctionOperatorTransposed<GridFunctionOperator<tag::testvec_trial, LocalContext, GridFct>,
+                                              GridFunctionOperator<tag::test_trialvec, LocalContext, GridFct>>
   {
     using Self = GridFunctionOperator;
-    using Transposed = GridFunctionOperator<tag::test_trialvec, GridFct>;
+    using Transposed = GridFunctionOperator<tag::test_trialvec, LocalContext, GridFct>;
     using Super = GridFunctionOperatorTransposed<Self, Transposed>;
 
   public:
diff --git a/src/amdis/assembler/ZeroOrderTestvecTrialvec.hpp b/src/amdis/assembler/ZeroOrderTestvecTrialvec.hpp
index fefccfbe0427f8351076a6e8f5375b0d519404bd..f8028777533a6387a55f8b1ef12cca1183f03b9b 100644
--- a/src/amdis/assembler/ZeroOrderTestvecTrialvec.hpp
+++ b/src/amdis/assembler/ZeroOrderTestvecTrialvec.hpp
@@ -19,12 +19,13 @@ namespace AMDiS
 
 
   /// zero-order operator \f$ \langle\Psi, c\,\Phi\rangle \f$, or \f$ \langle\Psi, A\,\Phi\rangle \f$
-  template <class GridFct>
-  class GridFunctionOperator<tag::testvec_trialvec, GridFct>
-      : public GridFunctionOperatorBase<GridFunctionOperator<tag::testvec_trialvec, GridFct>, GridFct>
+  template <class LocalContext, class GridFct>
+  class GridFunctionOperator<tag::testvec_trialvec, LocalContext, GridFct>
+      : public GridFunctionOperatorBase<GridFunctionOperator<tag::testvec_trialvec, LocalContext, GridFct>,
+                                        LocalContext, GridFct>
   {
     using Self = GridFunctionOperator;
-    using Super = GridFunctionOperatorBase<Self, GridFct>;
+    using Super = GridFunctionOperatorBase<Self, LocalContext, GridFct>;
 
     using expr_value_type = typename GridFct::Range;
     static_assert( Category::Scalar<expr_value_type> || Category::Matrix<expr_value_type>,
diff --git a/src/amdis/gridfunctions/AnalyticGridFunction.hpp b/src/amdis/gridfunctions/AnalyticGridFunction.hpp
index 0faaa28dba1fef64e0ae68923e6c167625d39ec2..055bfff1f0e2e4c08b288239f462699c65a82181 100644
--- a/src/amdis/gridfunctions/AnalyticGridFunction.hpp
+++ b/src/amdis/gridfunctions/AnalyticGridFunction.hpp
@@ -42,6 +42,7 @@ namespace AMDiS
 
     Range operator()(Domain const& local) const
     {
+      assert(geometry_.has_value());
       return fct_(geometry_.value().global(local));
     }
 
diff --git a/src/amdis/utility/Tuple.hpp b/src/amdis/utility/Tuple.hpp
index 68e47efcf9b23ae2c4eddcac93493ac92f8744b0..e6b946deb38fe90699f7495589d7cb916171f23d 100644
--- a/src/amdis/utility/Tuple.hpp
+++ b/src/amdis/utility/Tuple.hpp
@@ -29,7 +29,7 @@ namespace AMDiS
 
       // create value through a map of arguments
       template <class Map, class... Args>
-      TupleEntry(tag::mapped_arg, Map const& map, Args&&... args)
+      TupleEntry(tag::mapped_arg, Map&& map, Args&&... args)
         : value_(map(std::forward<Args>(args)...))
       {}
 
@@ -52,7 +52,7 @@ namespace AMDiS
     {
       TupleImpl() = default;
 
-      /// constructor from sequence of types, disabling this for copy and move constructor
+      // constructor from sequence of types, disabling this for copy and move constructor
       template <class... OtherTypes,
         std::enable_if_t<(sizeof...(OtherTypes) == sizeof...(Types)), int> = 0,
         std::enable_if_t<none_of_v<IsTupleImpl<std::decay_t<OtherTypes>>::value...>, int> = 0>
@@ -60,10 +60,11 @@ namespace AMDiS
         : TupleEntry<Indices, Types>{std::forward<OtherTypes>(other)}...
       {}
 
+      // constructor using a map to build the elements
       template <class Map, class... OtherTypes,
         std::enable_if_t<(sizeof...(OtherTypes) == sizeof...(Types)), int> = 0>
-      TupleImpl(tag::mapped_arg, Map const& map, OtherTypes&&... other)
-        : TupleEntry<Indices, Types>{tag::mapped_arg{}, map, std::forward<OtherTypes>(other)}...
+      TupleImpl(tag::mapped_arg, Map&& map, OtherTypes&&... other)
+        : TupleEntry<Indices, Types>{tag::mapped_arg{}, std::forward<Map>(map), std::forward<OtherTypes>(other)}...
       {}
     };
 
@@ -77,38 +78,39 @@ namespace AMDiS
     using Super = Impl::TupleImpl<std::index_sequence_for<Types...>, Types...>;
 
   public:
+    // some default constructors
     Tuple() = default;
     Tuple(Tuple const&) = default;
     Tuple(Tuple&&) = default;
     Tuple& operator=(Tuple const&) = default;
     Tuple& operator=(Tuple&&) = default;
 
-    // construct from types different from Types, but convertible
+    /// Construct from types different from Types, but convertible
     template <class... OtherTypes,
       std::enable_if_t<(sizeof...(OtherTypes) == sizeof...(Types)), int> = 0>
     explicit Tuple(OtherTypes&&... other)
       : Super{std::forward<OtherTypes>(other)...}
     {}
 
-    // constructor from tuple of possibly other types
+    /// Constructor from tuple of possibly other types
     template <class... OtherTypes,
       std::enable_if_t<(sizeof...(OtherTypes) == sizeof...(Types)), int> = 0>
     explicit Tuple(Tuple<OtherTypes...> const& that)
       : Super{that}
     {}
 
-    // constructor from tuple of possibly other types
+    /// Constructor from tuple of possibly other types
     template <class... OtherTypes,
       std::enable_if_t<(sizeof...(OtherTypes) == sizeof...(Types)), int> = 0>
     explicit Tuple(Tuple<OtherTypes...>&& that)
       : Super{std::move(that)}
     {}
 
-    // constructor using a map to build the elements
+    /// Constructor using a map to build the elements
     template <class Map, class... OtherTypes,
       std::enable_if_t<(sizeof...(OtherTypes) == sizeof...(Types)), int> = 0>
-    Tuple(tag::mapped_arg, Map const& map, OtherTypes&&... other)
-      : Super{tag::mapped_arg{}, map, std::forward<OtherTypes>(other)...}
+    Tuple(tag::mapped_arg, Map&& map, OtherTypes&&... other)
+      : Super{tag::mapped_arg{}, std::forward<Map>(map), std::forward<OtherTypes>(other)...}
     {}
 
     /// Get the I'th tuple element by reference
@@ -126,6 +128,12 @@ namespace AMDiS
       Impl::TupleEntry<I, TupleElement_t<I, AMDiS::Types<Types...>>> const& base = *this;
       return base.value_;
     }
+
+    /// Return the number of elements in this tuple
+    constexpr std::size_t size() const
+    {
+      return sizeof...(Types);
+    }
   };
 
 
@@ -233,8 +241,7 @@ namespace Dune
 {
   namespace Std
   {
-#if ! DUNE_HAVE_CXX_APPLY && ! DUNE_HAVE_CXX_EXPERIMENTAL_APPLY
-  using AMDiS::Impl::apply;
-#endif
+    using AMDiS::Impl::apply;
+
   } // end namespace Std
 } // end namespace Dune