Skip to content

#778 MatrixHelper::detAAT is unstable

Metadata

Property Value
Reported by Oliver Sander (oliver.sander@tu-dresden.de)
Reported at Apr 21, 2010 20:53
Type Bug Report
Version Git (pre2.4) [autotools]
Operating System Unspecified / All
Last edited by Christian Engwer (christi@conan.iwr.uni-heidelberg.de)
Last edited at Nov 1, 2011 17:10

Description

For certain almost singular matrices, for example

A[0][0] =  0.099999999999999867;
A[0][1] = -0.010000000000002118;
A[1][0] =  0.099999999999999867;
A[1][1] = -0.0099999999999998979;

the method detAAT in genericgeometries/matrix.hh aborts with

testmatrix: ../matrix.hh:243: static void Dune::GenericGeometry::MatrixHelper<Traits>::cholesky_L(const typename Traits::Matrix<n, n>::type&, typename Traits::Matrix<n, n>::type&) [with int n = 2, Traits = Dune::GenericGeometry::DuneCoordTraits<double>]: Assertion `x > FieldType( 0 )' failed.

rather than computing something close to zero.

This is not an artificial example. I found it when integrating over remote intersections between two unrelated grids.

A test case has been added in the file genericgeometries/test/testmatrix.cc.