From 2bbd7ca0369cb394638057f792e70014d6fe2513 Mon Sep 17 00:00:00 2001
From: podlesny <podlesny@zedat.fu-berlin.de>
Date: Sat, 6 Mar 2021 13:27:55 +0100
Subject: [PATCH] rescale DirectionalRestriction to direction.two_norm = 1.0
---
.../spatial-solving/tnnmg/functional.hh | 25 ++++++++++++++++---
1 file changed, 21 insertions(+), 4 deletions(-)
diff --git a/dune/tectonic/spatial-solving/tnnmg/functional.hh b/dune/tectonic/spatial-solving/tnnmg/functional.hh
index 2e76d0aa..15674fb3 100755
--- a/dune/tectonic/spatial-solving/tnnmg/functional.hh
+++ b/dune/tectonic/spatial-solving/tnnmg/functional.hh
@@ -162,11 +162,12 @@ class DirectionalRestriction :
using Matrix = typename Base::Matrix;
using Vector = typename Base::Vector;
- DirectionalRestriction(const GlobalMatrix& matrix, const GlobalVector& linearTerm, const Nonlinearity& phi, const GlobalLowerObstacle& lower, const GlobalLowerObstacle& upper, const GlobalVector& origin, const GlobalVector& direction) :
+ DirectionalRestriction(const GlobalMatrix& matrix, const GlobalVector& linearTerm, const Nonlinearity& phi, const GlobalLowerObstacle& lower, const GlobalLowerObstacle& upper, const GlobalVector& origin, const GlobalVector& direction, const double scaling) :
Base(matrix, linearTerm, lower, upper, origin, direction),
origin_(origin),
direction_(direction),
- phi_(phi)
+ phi_(phi),
+ scaling_(scaling)
{
phi_.directionalDomain(origin_, direction_, domain_);
@@ -207,12 +208,17 @@ class DirectionalRestriction :
return direction_;
}
+ double scaling() const {
+ return scaling_;
+ }
protected:
const GlobalVector& origin_;
- const GlobalVector& direction_;
+ GlobalVector direction_;
const Nonlinearity& phi_;
Interval domain_;
+
+ const double scaling_;
};
@@ -474,7 +480,18 @@ class Functional : public Dune::TNNMG::BoxConstrainedQuadraticFunctional<M, V, L
friend auto directionalRestriction(const Functional& f, const Vector& origin, const Vector& direction)
-> DirectionalRestriction<Matrix, Vector, Nonlinearity, LowerObstacle, UpperObstacle, Range>
{
- return DirectionalRestriction<Matrix, Vector, Nonlinearity, LowerObstacle, UpperObstacle, Range>(f.quadraticPart(), f.linearPart(), f.phi(), f.lowerObstacle(), f.upperObstacle(), origin, direction);
+ // rescale direction for numerical stability
+ auto dir = direction;
+ auto dirNorm = dir.two_norm();
+ if (dirNorm > 0.0)
+ dir /= dirNorm;
+ else {
+ dirNorm = 0.0;
+ dir = 0.0;
+ }
+
+
+ return DirectionalRestriction<Matrix, Vector, Nonlinearity, LowerObstacle, UpperObstacle, Range>(f.quadraticPart(), f.linearPart(), f.phi(), f.lowerObstacle(), f.upperObstacle(), origin, dir, dirNorm);
}
friend auto shift(const Functional& f, const Vector& origin)
--
GitLab