diff --git a/dune/tectonic/minimisation.hh b/dune/tectonic/minimisation.hh index 389a0a1498576f2e7fd97efd66954c4b0f31420c..71edb09b22115375c13e9f6076d4767b2497f56c 100644 --- a/dune/tectonic/minimisation.hh +++ b/dune/tectonic/minimisation.hh @@ -12,22 +12,16 @@ #include "mydirectionalconvexfunction.hh" template <class Functional> -void descentMinimisation(Functional const &J, - typename Functional::LocalVector &x, - typename Functional::LocalVector const &v, - Bisection const &bisection) { - using LocalVector = typename Functional::LocalVector; - using LocalNonlinearity = typename Functional::Nonlinearity; - - MyDirectionalConvexFunction<LocalNonlinearity> const JRest( +double lineSearch(Functional const &J, + typename Functional::LocalVector const &x, + typename Functional::LocalVector const &v, + Bisection const &bisection) { + MyDirectionalConvexFunction<typename Functional::Nonlinearity> const JRest( computeDirectionalA(J.A, v), computeDirectionalb(J.A, J.b, x, v), *J.phi, x, v); - // }}} int count; - double const stepsize = bisection.minimize(JRest, 0.0, 0.0, count); - - Arithmetic::addProduct(x, stepsize, v); + return bisection.minimize(JRest, 0.0, 0.0, count); } template <class Functional> @@ -38,12 +32,14 @@ void minimise(Functional const &J, typename Functional::LocalVector &x, for (size_t step = 0; step < steps; ++step) { LocalVector v; J.gradient(x, v); - if (v.two_norm() == 0.0) - return; + double const vnorm = v.two_norm(); + if (vnorm <= 0.0) + return; v *= -1; - descentMinimisation(J, x, v, bisection); + double const alpha = lineSearch(J, x, v, bisection); + Arithmetic::addProduct(x, alpha, v); } } #endif