Add descentAtZero()

dune/tectonic/ellipticenergy.hh

@@ -31,27 +31,32 @@ template <int dim> class EllipticEnergy {
     return y * v + (*phi)(v); // <1/2 Av - b,v> + H(|v|)
   }
 
+  void descentAtZero(SmallVector &ret) const {
+    SmallVector const zero(0);
+    // If there is a direction of descent, this is it
+    SmallVector d;
+    smoothGradient(zero, d);
+    d *= -1;
+
+    Interval<double> D;
+    phi->directionalSubDiff(zero, d, D);
+    double const nonlinearDecline = D[1];
+    double const smoothDecline = -(d * d);
+    double const combinedDecline = smoothDecline + nonlinearDecline;
+
+    if (combinedDecline < 0) {
+      ret = d;
+    } else {
+      ret = 0;
+    }
+  }
+
   // returns false if the direction is tangential
   bool descentDirection(SmallVector const &x, SmallVector &ret) const {
     // Check the squared norm rather than each component because
     // complementaryProjection() divides by it
     if (x.two_norm2() == 0.0) {
-      // If there is a direction of descent, this is it
-      SmallVector d;
-      smoothGradient(x, d);
-      d *= -1;
-
-      Interval<double> D;
-      phi->directionalSubDiff(x, d, D);
-      double const nonlinearDecline = D[1];
-      double const smoothDecline = -(d * d);
-      double const combinedDecline = smoothDecline + nonlinearDecline;
-
-      if (combinedDecline < 0) {
-        ret = d;
-      } else {
-        ret = 0;
-      }
+      descentAtZero(ret);
       return true;
     }