Fix evil bug in output code

The computation of element centers was wrong (it was in global
coordinates rather than in local ones).  This did not matter
for quasiconvexity-test.cc (because there all shape function
derivatives are constant, and it doesn't matter where you
evaluate them).  But it lead to wrong results in
quasiconvexity-test-micromorphic.cc.

src/quasiconvexity-test-micromorphic.cc

@@ -76,7 +76,7 @@ void writeDisplacement(const Basis& basis, const Vector& deformation, const Fiel
     // store values of shape functions
     std::vector<Dune::FieldVector<double,dim> > shapeFunctionValues(localFiniteElement.size());
 
-    auto p = element.geometry().center();
+    auto p = referenceElement(element).position(0,0);  // Element center in local coordinates
 
     // get gradients of shape functions
     localFiniteElement.localBasis().evaluateFunction(p, shapeFunctionValues);

src/quasiconvexity-test.cc

@@ -139,7 +139,7 @@ void writeDisplacement(const Basis& basis, const Vector& periodicX, const FieldM
     std::vector<Dune::FieldMatrix<double,1,dim> > referenceGradients(localFiniteElement.size());
     std::vector<Dune::FieldVector<double,dim> > gradients(localFiniteElement.size());
 
-    auto p = element.geometry().center();
+    auto p = referenceElement(element).position(0,0);  // Element center in local coordinates
 
     const auto jacobianInverseTransposed = element.geometry().jacobianInverseTransposed(p);