diff --git a/dune/tectonic/circularconvexfunction.hh b/dune/tectonic/circularconvexfunction.hh deleted file mode 100644 index 65f79b5d65971fdfa84e868d3975a8ea7657236f..0000000000000000000000000000000000000000 --- a/dune/tectonic/circularconvexfunction.hh +++ /dev/null @@ -1,71 +0,0 @@ -#ifndef CIRCULAR_CONVEX_FUNCTION_HH -#define CIRCULAR_CONVEX_FUNCTION_HH - -#include <cmath> -#include <dune/fufem/interval.hh> - -namespace Dune { -template <class MatrixType, class VectorType> class CircularConvexFunction { -public: - CircularConvexFunction(MatrixType const &A, VectorType const &b, - VectorType const &x, VectorType const &dir) - : A(A), - b(b), - x(x), - dir(dir), - xnorm(x.two_norm()), - dnorm(dir.two_norm()) {} - - double quadraticPart() const { return 0; } - - double linearPart() const { return 0; } - - void subDiff(double m, Interval<double> &D) const { - VectorType x; - cartesian(m, x); - VectorType t; - tangent(m, t); - - VectorType tmp; - A.mv(x, tmp); // Ax - tmp -= b; // Ax - b - D[0] = D[1] = tmp * t; // <Ax - b,t> - } - - void domain(Interval<double> &domain) const { - domain[0] = -2 * M_PI; - domain[1] = 2 * M_PI; - } - - void cartesian(double m, VectorType &y) const { - y = 0; - y.axpy(std::cos(m), x); - y.axpy(std::sin(m) * xnorm / dnorm, dir); - } - -private: - MatrixType const &A; - VectorType const &b; - VectorType const &x; - VectorType const &dir; - - double const dnorm; - double const xnorm; - - /* If x and d were normalised, we would have - - cartesian(a) = cos(a) * x + sin(a) * d and - tangent(a) = -sin(a) * x + cos(a) * d. - - Since we x and d are not normalised and the return of - cartesian() is fixed, we scale the tangent. - */ - void tangent(double m, VectorType &y) const { - y = 0; - y.axpy(-std::sin(m) * dnorm, x); - y.axpy(std::cos(m) * xnorm, dir); - } -}; -} - -#endif diff --git a/dune/tectonic/minimisation.hh b/dune/tectonic/minimisation.hh index 83488b0f3fd02c9985160c3a5751dafb30f9be79..9bae48be8e2f0189db270f91c27ce1f9077842ed 100644 --- a/dune/tectonic/minimisation.hh +++ b/dune/tectonic/minimisation.hh @@ -8,7 +8,6 @@ #include <dune/fufem/interval.hh> #include <dune/tnnmg/problem-classes/bisection.hh> -#include "circularconvexfunction.hh" #include "mydirectionalconvexfunction.hh" namespace Dune { @@ -70,31 +69,6 @@ void descentMinimisation(Functional const &J, x.axpy(stepsize, descDir); } -template <class Functional> -void tangentialMinimisation(Functional const &J, - typename Functional::SmallVector &x, - typename Functional::SmallVector const &descDir, - Bisection const &bisection) { - using SmallMatrix = typename Functional::SmallMatrix; - using SmallVector = typename Functional::SmallVector; - - // We completely ignore the nonlinearity here -- when restricted - // to a circle, it just enters as a constant! - CircularConvexFunction<SmallMatrix, SmallVector> const JRest(J.A, J.b, x, - descDir); - - int count; - double const stepsize = bisection.minimize(JRest, 0.0, 1.0, count); - dverb << "Number of iterations in the bisection method: " << count - << std::endl; - ; - - // Since x is used in the computation of the rhs, do not write to it directly - SmallVector tmp; - JRest.cartesian(stepsize, tmp); - x = tmp; -} - template <class Functional> void minimise(Functional const &J, typename Functional::SmallVector &x, size_t steps, Bisection const &bisection) { @@ -102,19 +76,12 @@ void minimise(Functional const &J, typename Functional::SmallVector &x, for (size_t step = 0; step < steps; ++step) { SmallVector descDir; - bool const linesearchp = J.descentDirection(x, descDir); + J.descentDirection(x, descDir); if (descDir.two_norm() < 1e-14) // TODO: Make controllable return; - if (linesearchp) { - descentMinimisation(J, x, descDir, bisection); - } else { - Bisection slowBisection(bisection); - slowBisection.setFastQuadratic(false); - - tangentialMinimisation(J, x, descDir, slowBisection); - } + descentMinimisation(J, x, descDir, bisection); } } }