Implement method derivativeOfQuaternionToMatrix

This is necessary to compute H1 errors in matrix space.

dune/gfe/rotation.hh

@@ -892,6 +892,27 @@ public:
 
     }
 
+    /** \brief Derivative of the map from unit quaternions to orthogonal matrices
+     */
+    static Tensor3<T,3,3,4> derivativeOfQuaternionToMatrix(const Dune::FieldVector<T,4>& p)
+    {
+      Tensor3<T,3,3,4> result;
+
+      result[0][0] = { 2*p[0], -2*p[1], -2*p[2],  2*p[3]};
+      result[0][1] = { 2*p[1],  2*p[0], -2*p[3], -2*p[2]};
+      result[0][2] = { 2*p[2],  2*p[3],  2*p[0],  2*p[1]};
+
+      result[1][0] = { 2*p[1],  2*p[0],  2*p[3],  2*p[2]};
+      result[1][1] = {-2*p[0],  2*p[1], -2*p[2],  2*p[3]};
+      result[1][2] = {-2*p[3],  2*p[2],  2*p[1], -2*p[0]};
+
+      result[2][0] = { 2*p[2], -2*p[3],  2*p[0], -2*p[1]};
+      result[2][1] = { 2*p[3],  2*p[2],  2*p[1],  2*p[0]};
+      result[2][2] = {-2*p[0], -2*p[1],  2*p[2],  2*p[3]};
+
+      return result;
+    }
+
     /** \brief Set rotation from orthogonal matrix
 
     We tacitly assume that the matrix really is orthogonal */