From 3ea1c8e913307dce98e3f206974fd386958b01e7 Mon Sep 17 00:00:00 2001
From: Elias Pipping <elias.pipping@fu-berlin.de>
Date: Mon, 19 Dec 2011 14:53:26 +0100
Subject: [PATCH] Rescale for numerical stability

---
 dune/tectonic/myblockproblem.hh | 10 +++++++---
 src/TODO.org                    |  1 -
 2 files changed, 7 insertions(+), 4 deletions(-)

diff --git a/dune/tectonic/myblockproblem.hh b/dune/tectonic/myblockproblem.hh
index 6e43c2b1..06b05444 100644
--- a/dune/tectonic/myblockproblem.hh
+++ b/dune/tectonic/myblockproblem.hh
@@ -74,7 +74,11 @@ template <class MyConvexProblemTypeTEMPLATE> class MyBlockProblem {
 
   double computeDampingParameter(VectorType const &u,
                                  VectorType const &projected_v) const {
-    VectorType const v = projected_v;
+    VectorType v = projected_v;
+
+    double const vnorm = v.two_norm();
+
+    v /= vnorm; // Rescale for numerical stability
 
     VectorType tmp = problem.f;
     problem.A.mmv(u, tmp);
@@ -94,8 +98,8 @@ template <class MyConvexProblemTypeTEMPLATE> class MyBlockProblem {
         localA, localb, problem.phi, u, v);
 
     int bisectionsteps = 0;
-    Bisection bisection(0.0, 1.0, 1e-12, true, 0);            // TODO
-    return bisection.minimize(psi, 1.0, 0.0, bisectionsteps); // TODO
+    Bisection bisection(0.0, 1.0, 1e-12, true, 0);                      // TODO
+    return bisection.minimize(psi, vnorm, 0.0, bisectionsteps) / vnorm; // TODO
   }
 
   void assembleTruncate(VectorType const &u, Linearization &linearization,
diff --git a/src/TODO.org b/src/TODO.org
index 5e66a10c..60d0487b 100644
--- a/src/TODO.org
+++ b/src/TODO.org
@@ -1,3 +1,2 @@
-* check if v should be normalised by default
 * use nested iteration to obtain better iteration to start with
 * fix up octave bindings (low-pri)
-- 
GitLab