Give the DR tensor variable the correct dimensions

The old version kept was too big, and kept a useless row of zeros around.
Removing this row should lead to a tiny increase in efficiency.  More importantly,
it gets us closer to the CosseratEnergyLocalStiffness class again, which
is scheduled to replace MixedCosseratEnergy eventually.

dune/gfe/mixedcosseratenergy.hh

@@ -90,7 +90,7 @@ public:  // for testing
      */
     static void computeDR(const Rotation<field_type,3>& value,
                           const Dune::FieldMatrix<field_type,4,gridDim>& derivative,
-                          Tensor3<field_type,3,3,3>& DR)
+                          Tensor3<field_type,3,3,gridDim>& DR)
     {
         // The LocalGFEFunction class gives us the derivatives of the orientation variable,
         // but as a map into quaternion space.  To obtain matrix coordinates we use the
@@ -210,7 +210,7 @@ public:
                 + (mu_*lambda_)/(2*mu_ + lambda_) * 0.5 * ((detU-1)*(detU-1) + (1.0/detU -1)*(1.0/detU -1));
     }
 
-    RT curvatureEnergy(const Tensor3<field_type,3,3,3>& DR) const
+    RT curvatureEnergy(const Tensor3<field_type,3,3,gridDim>& DR) const
     {
 #ifdef DONT_USE_CURL
         return mu_ * std::pow(L_c_ * L_c_ * DR.frobenius_norm2(),q_/2.0);
@@ -219,18 +219,14 @@ public:
 #endif
     }
 
-    RT bendingEnergy(const Dune::FieldMatrix<field_type,dim,dim>& R, const Tensor3<field_type,3,3,3>& DR) const
+    RT bendingEnergy(const Dune::FieldMatrix<field_type,dim,dim>& R, const Tensor3<field_type,3,3,gridDim>& DR) const
     {
         // left-multiply the derivative of the third director (in DR[][2][]) with R^T
-        Dune::FieldMatrix<field_type,3,3> RT_DR3;
+        Dune::FieldMatrix<field_type,3,3> RT_DR3(0);
         for (int i=0; i<3; i++)
-            for (int j=0; j<3; j++) {
-                RT_DR3[i][j] = 0;
+            for (int j=0; j<gridDim; j++)
                 for (int k=0; k<3; k++)
                     RT_DR3[i][j] += R[k][i] * DR[k][2][j];
-            }
-
-
 
         return mu_ * sym(RT_DR3).frobenius_norm2()
                + mu_c_ * skew(RT_DR3).frobenius_norm2()
@@ -335,7 +331,7 @@ energy(const typename Basis::LocalView& localView,
         //  Note: we need it in matrix coordinates
         //////////////////////////////////////////////////////////
 
-        Tensor3<field_type,3,3,3> DR;
+        Tensor3<field_type,3,3,gridDim> DR;
         computeDR(orientationValue, orientationDerivative, DR);
 
         // Add the local energy density