Add some functions for debugging

dune/tectonic/myblockproblem.hh

@@ -15,6 +15,22 @@
 #include "mydirectionalconvexfunction.hh"
 #include "samplefunctional.hh"
 
+/* Just for debugging */
+template <int dim, class VectorType, class MatrixType>
+double computeEnergy(
+    MatrixType const &A, VectorType const &b,
+    Dune::GlobalNonlinearity<dim, VectorType, MatrixType> const &phi,
+    VectorType const &x) {
+  double ret;
+  VectorType tmp(x.size());
+  A.mv(x, tmp);          //            Ax
+  ret = 0.5 * (tmp * x); // ret = 1/2 <Ax,x>
+
+  ret -= b * x;  // ret = 1/2 <Ax,x> - <b,x>
+  ret += phi(x); // ret = 1/2 <Ax,x> - <b,x> + phi(x)
+  return ret;
+}
+
 /** \brief Base class for problems where each block can be solved with a
  * modified gradient method */
 template <class MyConvexProblemTypeTEMPLATE> class MyBlockProblem {

dune/tectonic/mydirectionalconvexfunction.hh

@@ -17,6 +17,13 @@ template <class NonlinearityType> class MyDirectionalConvexFunction {
     phi_.directionalDomain(u_, v_, dom_);
   }
 
+  /* Just for debugging */
+  double operator()(double x) const {
+    VectorType tmp = v_;
+    tmp *= x;
+    return (0.5 * A * x * x) - (b * x) + phi_(tmp);
+  }
+
   double quadraticPart() const { return A; }
 
   double linearPart() const { return b; }