Commit 9698e497 authored by Jö Fahlke's avatar Jö Fahlke

[!262] [!221] [TLIME] Relax convergence criterion.

Merge branch 'cherry-pick-9da2de98' into 'releases/2.6'

ref:core/dune-istl Merge branch 'fix-tlime-residual-limit' into 'master'

The previous convergence limit was originally determined experimentally as
1e-11. This worked for many blas implementations and architectures. However,
when used with openblas on skylake, apparently the residual norm would not go
below ~1e-10, so convergence was never achieved. In fact, even on non-skylake
the residual norm would go above 1e-11 again after briefly dipping below, if
iterating further.

We believe that this is due to openblas selecting \-- at runtime \-- some
skylake specific algorithm leading to a different ordering of operations, in
turn leading to differences in numerical cancellation. We have however not
verified this conclusively, nor have we identified precisely which blas
algorithm is causing this.

This patch raises the convergence limit to
sqrt(numeric_limits<field_type>::epsilon()). This limit has no theoretical
justification \-- it was selected because it usually works as a convergence
limit for other (completely unrelated) algorithms, and because it works for
both Skylake and other architectures (AMD Epyc) in this particular case.

Developed together with Sebastian Westerheide.

Fixes: [#48].

See merge request [!221]

(cherry picked from commit 9da2de98)

0bb1e69d \[TLIME\] Relax convergence criterion.

Closes: [#61] (on 2.6).

See merge request [!262]

  [#48]: gitlab.dune-project.org/NoneNone/issues/48
  [!221]: gitlab.dune-project.org/core/dune-istl/merge_requests/221
  [#61]: gitlab.dune-project.org/NoneNone/issues/61
  [!262]: gitlab.dune-project.org/core/dune-istl/merge_requests/262
parents 8f33d649 ae21be67
Pipeline #15294 passed with stage
in 12 minutes and 21 seconds
......@@ -5,7 +5,7 @@
#include <cmath> // provides std::abs and std::sqrt
#include <cassert> // provides assert
#include <limits>
#include <iostream> // provides std::cout, std::endl
#include <dune/common/exceptions.hh> // provides DUNE_THROW(...), Dune::Exception
......@@ -191,7 +191,7 @@ protected:
// 7) get smallest magnitude eigenvalue via TLIME iteration
// (assume that m_ is symmetric)
{
const Real epsilon = 1e-11;
const Real epsilon = std::sqrt(std::numeric_limits<Real>::epsilon());
x = 1.0;
// 7.1) perform TLIME iteration for smallest magnitude
// eigenvalue
......@@ -322,7 +322,7 @@ protected:
// 9) get smallest singular value as square root of the smallest
// magnitude eigenvalue of m^t*m via TLIME iteration
{
const Real epsilon = 1e-11;
const Real epsilon = std::sqrt(std::numeric_limits<Real>::epsilon());
x = 1.0;
// 9.1) perform TLIME iteration for smallest magnitude
// eigenvalue of m^t*m
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment