# [Dune] [#778] MatrixHelper::detAAT is unstable

Dune flyspray at dune-project.org
Thu Apr 22 08:36:24 CEST 2010

```THIS IS AN AUTOMATED MESSAGE, DO NOT REPLY.

FS#778 - MatrixHelper::detAAT is unstable
User who did this - Martin Nolte (nolte)

----------
I think I don't get your exact point. The determinant in your example is very nearly zero:

gnuplot> print 0.099999999999999867 * (-0.0099999999999998979) - (-0.010000000000002118) * 0.099999999999999867
2.22044604925031e-16

I have to admit, though, that the above formula should be numerically stable in this context.

Now, detAAT first (A * A^T) and then its determinant, i.e.,  the determinant of (computed by gnuplot)
AAT = 0.0101
AAT = 0.0101
AAT = 0.0101
AAT = 0.0101
The determinant of this matrix is obviously zero.

Moreover, the matrix has one eigenvalue pretty close to 0 and one eigenvalue around 0.01, i.e., it is quite ill-conditioned. I wonder if the inverse can be computed numerically (with double precision)?

Since this MatrixHelper shall only be used for the Jacobians of geometries, the following questions arise:
(a) is there an algorithm that is numerically more stable (and works in the nonsquare case)?
(b) do we really need to consider geometries where the jacobian is this ill-conditioned?

If someone really must use such ill-conditioned geometries, I think he is well-advised to use some higher precision floating point type for the geometry's ctype.
----------