From 6588b185a8b76041afc6c9ce521c53a33b6f5ff1 Mon Sep 17 00:00:00 2001 From: Elias Pipping <elias.pipping@fu-berlin.de> Date: Fri, 9 Sep 2011 18:19:03 +0200 Subject: [PATCH] Use bisection class from dune-tnnmg --- src/bisection-example-new.cc | 65 +++++++----------------------------- 1 file changed, 12 insertions(+), 53 deletions(-) diff --git a/src/bisection-example-new.cc b/src/bisection-example-new.cc index 6df06d37..c44375fe 100644 --- a/src/bisection-example-new.cc +++ b/src/bisection-example-new.cc @@ -49,7 +49,7 @@ class SampleFunctional { return (x * dir > 0) ? PlusGrad(x) * dir : MinusGrad(x) * dir; } - SmallVector minimise(const SmallVector x, unsigned int iterations) const { + SmallVector minimise(const SmallVector x) const { SmallVector descDir = ModifiedGradient(x); if (descDir == SmallVector(0.0)) @@ -75,57 +75,14 @@ class SampleFunctional { MyDirectionalConvexFunctionType; MyDirectionalConvexFunctionType rest(rest_A, rest_b, phi, x, descDir); - Interval<double> D; - - // { Debugging - Dune::dverb << "Starting at x with J(x) = " << operator()(x) << std::endl; - SmallVector debug_normalised_dir = descDir; - debug_normalised_dir /= debug_normalised_dir.two_norm(); - Dune::dverb << "Minimizing in direction w with dJ(x,w) = " - << directionalDerivative(x, debug_normalised_dir) << std::endl; - // } - - double l = 0; - double r = 1; - while (true) { - rest.subDiff(r, D); - if (D[1] >= 0) - break; - - l = r; - r *= 2; - Dune::dverb << "Widened interval!" << std::endl; - } - Dune::dverb << "Interval now [" << l << "," << r << "]" << std::endl; - - { // Debugging - SmallVector tmpl = x; - tmpl.axpy(l, descDir); - SmallVector tmpr = x; - tmpr.axpy(r, descDir); - assert(directionalDerivative(tmpl, debug_normalised_dir) < 0); - assert(directionalDerivative(tmpr, debug_normalised_dir) > 0); - } + Bisection bisection; // FIXME: default values + int count; + // FIXME: does x_old = 1 make any sense?! + double const m = bisection.minimize(rest, 0.0, 1.0, count); + Dune::dverb << "Number of iterations in the bisection method: " << count + << std::endl; + ; - double m = l / 2 + r / 2; - for (unsigned int count = 0; count < iterations; ++count) { - { // Debugging - Dune::dverb << "now at m = " << m << std::endl; - SmallVector tmp = x; - tmp.axpy(m, descDir); - Dune::dverb << "Value of J here: " << operator()(tmp) << std::endl; - } - - rest.subDiff(m, D); - if (D[1] < 0) { - l = m; - m = (m + r) / 2; - } else if (D[1] > 0) { - r = m; - m = (l + m) / 2; - } else - break; - } SmallVector middle = descDir; middle *= m; return middle; @@ -220,8 +177,9 @@ void testSampleFunction() { SampleFunctional::SmallVector correction; for (int i = 1; i <= 6; ++i) { - correction = J.minimise(start, 20); + correction = J.minimise(start); start += correction; + std::cout << "New value: J(...) = " << J(start) << std::endl; } std::cout << "Arrived at J(...) = " << J(start) << std::endl; } @@ -251,8 +209,9 @@ void testTrivialFunction() { SampleFunctional::SmallVector correction; for (int i = 1; i <= 5; ++i) { - correction = J.minimise(start, 20); + correction = J.minimise(start); start += correction; + std::cout << "New value: J(...) = " << J(start) << std::endl; } std::cout << "Arrived at J(...) = " << J(start) << std::endl; } -- GitLab