From d96d93a81ed03e0539bbb5b9ccc9f23c4ae34dfd Mon Sep 17 00:00:00 2001
From: Elias Pipping <elias.pipping@fu-berlin.de>
Date: Wed, 3 Jul 2013 13:01:32 +0200
Subject: [PATCH] [Cleanup] Helpers for MyDirectionalConvexFunctional
---
dune/tectonic/minimisation.hh | 16 ++------------
dune/tectonic/myblockproblem.hh | 18 +++------------
dune/tectonic/mydirectionalconvexfunction.hh | 23 ++++++++++++++++++++
3 files changed, 28 insertions(+), 29 deletions(-)
diff --git a/dune/tectonic/minimisation.hh b/dune/tectonic/minimisation.hh
index 597ba23d..9ed960ea 100644
--- a/dune/tectonic/minimisation.hh
+++ b/dune/tectonic/minimisation.hh
@@ -7,7 +7,6 @@
#include <dune/fufem/arithmetic.hh>
#include <dune/fufem/interval.hh>
-#include <dune/solvers/computeenergy.hh>
#include <dune/tnnmg/problem-classes/bisection.hh>
#include "mydirectionalconvexfunction.hh"
@@ -21,20 +20,9 @@ void descentMinimisation(Functional const &J,
using SmallVector = typename Functional::SmallVector;
using LocalNonlinearityType = typename Functional::NonlinearityType;
- // {{{ Construct a restriction of J to the line x + t * v
-
- /* We have
-
- 1/2 <A(x+tv),x+tv>-<b,x+tv> = 1/2 <Av,v> t^2 - <b-Ax,v> t + <1/2 Ax-b,x>
-
- since A is symmetric.
- */
- SmallVector tmp = J.b;
- Arithmetic::addProduct(tmp, -1.0, J.A, x);
- double const JRestb = tmp * v;
-
MyDirectionalConvexFunction<LocalNonlinearityType> const JRest(
- 2.0 * computeEnergy(J.A, v), JRestb, *J.phi, x, v);
+ computeDirectionalA(J.A, v), computeDirectionalb(J.A, J.b, x, v), *J.phi,
+ x, v);
// }}}
int count;
diff --git a/dune/tectonic/myblockproblem.hh b/dune/tectonic/myblockproblem.hh
index cb88e8fd..f1c8eae9 100644
--- a/dune/tectonic/myblockproblem.hh
+++ b/dune/tectonic/myblockproblem.hh
@@ -105,22 +105,10 @@ template <class MyConvexProblemTypeTEMPLATE> class MyBlockProblem {
Arithmetic::addProduct(tmp, -1.0, problem.A, u);
double const localb = tmp * v;
- problem.A.mv(v, tmp);
- double const localA = tmp * v;
-
- /*
- 1/2 <A(u + hv),u + hv> - <b, u + hv>
- = 1/2 <Av,v> h^2 - <b - Au, v> h + const.
-
- localA = <Av,v>
- localb = <b - Au, v>
- */
-
MyDirectionalConvexFunction<
- Dune::GlobalNonlinearity<MatrixType, VectorType>> const psi(localA,
- localb,
- problem.phi,
- u, v);
+ Dune::GlobalNonlinearity<MatrixType, VectorType>> const
+ psi(computeDirectionalA(problem.A, v),
+ computeDirectionalb(problem.A, problem.f, u, v), problem.phi, u, v);
Interval<double> D;
psi.subDiff(0, D);
diff --git a/dune/tectonic/mydirectionalconvexfunction.hh b/dune/tectonic/mydirectionalconvexfunction.hh
index 154aa7a5..0b2016d4 100644
--- a/dune/tectonic/mydirectionalconvexfunction.hh
+++ b/dune/tectonic/mydirectionalconvexfunction.hh
@@ -9,6 +9,29 @@
#include <dune/fufem/interval.hh>
+/*
+ 1/2 <A(u + hv),u + hv> - <b, u + hv>
+ = 1/2 <Av,v> h^2 - <b - Au, v> h + const.
+
+ localA = <Av,v>
+ localb = <b - Au, v>
+*/
+
+template <class MatrixType, class VectorType>
+double computeDirectionalA(MatrixType const &A, VectorType const &v) {
+ VectorType tmp(v.size());
+ A.mv(v, tmp);
+ return tmp * v;
+}
+
+template <class MatrixType, class VectorType>
+double computeDirectionalb(MatrixType const &A, VectorType const &b,
+ VectorType const &u, VectorType const &v) {
+ VectorType tmp = b;
+ Arithmetic::addProduct(tmp, -1.0, A, u);
+ return tmp * v;
+}
+
template <class NonlinearityType> class MyDirectionalConvexFunction {
public:
using VectorType = typename NonlinearityType::VectorType;
--
GitLab