diff --git a/dune/microstructure/CorrectorComputer.hh b/dune/microstructure/CorrectorComputer.hh
index c364e83f78c9c34050268a721d72aca487661253..c11727750567d796725c90ac33c43f817a0437ba 100644
--- a/dune/microstructure/CorrectorComputer.hh
+++ b/dune/microstructure/CorrectorComputer.hh
@@ -241,7 +241,12 @@ public:
   //                                   const localFunction2& lambda
   //                                   )
 
+  /** \brief Compute the stiffness matrix for one element
+   *
+   * \param quadRule The quadrature rule to be used
+   */
   void computeElementStiffnessMatrix(const typename Basis::LocalView& localView,
+                                     const Dune::QuadratureRule<double,dim>& quadRule,
                                     ElementMatrixCT& elementMatrix
                                     )
   {
@@ -274,20 +279,8 @@ public:
   //   std::cout << "localView.size(): " << localView.size() << std::endl;
   //   std::cout << "nSf : " << nSf << std::endl;
 
-  //   int orderQR = 2*(dim*localFiniteElement.localBasis().order()-1)+5;  // TEST
-    int orderQR = 2*(dim*localFiniteElement.localBasis().order()-1);
-  //   int orderQR = 0;
-  //   int orderQR = 1;
-  //   int orderQR = 2;
-  //   int orderQR = 3;
-    const auto& quad = Dune::QuadratureRules<double,dim>::rule(element.type(), orderQR);
-    
-  //   double elementContribution = 0.0;
-    
-  //   std::cout << "Print QuadratureOrder:" << orderQR << std::endl;  //TEST`
-
     int QPcounter= 0;
-    for (const auto& quadPoint : quad)
+    for (const auto& quadPoint : quadRule)
     {
       QPcounter++;
       const auto& quadPos = quadPoint.position();
@@ -659,9 +652,15 @@ public:
       const auto localPhiOffset = localView.size();
   //     Dune::Timer Time;
       //std::cout << "localPhiOffset : " << localPhiOffset << std::endl;
+
+      // Determine a suitable quadrature rule
+      const auto& localFiniteElement = localView.tree().child(0).finiteElement();
+      int orderQR = 2*(dim*localFiniteElement.localBasis().order()-1);
+      const auto& quadRule = Dune::QuadratureRules<double,dim>::rule(element.type(), orderQR);
+
       ElementMatrixCT elementMatrix;
       // computeElementStiffnessMatrix(localView, elementMatrix, muLocal, lambdaLocal);
-      computeElementStiffnessMatrix(localView, elementMatrix);
+      computeElementStiffnessMatrix(localView, quadRule, elementMatrix);
       
       
   //     std::cout << "local assembly time:" << Time.elapsed() << std::endl;