Numerical problems in matrix factorization

[fabio.schoen]

I am using MTL (a c++ interface to lapack), so I am not sure this is a Lapack or a MTL problem. In any case I ask this user group hoping in some hints.

Assume I try to solve the linear system Ax=b with
A=
[1,1,2 ],
[2,2,2],
[2,2,1]
and b:
[15],
[15],
[15]

obviously I get an error (A is singular and with this b there is no solution). In fact mtl method getsv, which interfaces to lapack, returns:

> Factorization failed with INFO = 2


However if I change A to
[1,1,2 ],
[2.1,2.1,2.1],
[2,2,1 ]
I get a solution!:

> A:
> 3x3
> [
> [1,1,2],
> [2.1,2.1,2.1],
> [2,2,1]
> ]
> b:
> 3x1
> [
> [15],
> [15],
> [15]
> ]
> x:
> 3x1
> [
> [-2.16004e+17],
> [2.16004e+17],
> [5]
> ]
> A:
> 3x3
> [
> [2.1,2.1,2.1],
> [0.952381,2.64545e-17,-1],
> [0.47619,0.5,1.5]
> ]


Thank you for any help,
fabio
This is clearly an error. Did anybody noticed it before?

[sven]

I don't know which LAPACK routine MTL interfaces to, but if for example it is DGESV then what you see is possible. DGESV only tests for an exact zero in the factorization, but in your second example a(2,2) is very small but not exactly zero.

It is a good idea to also obtain an estimate of the condition number of A to see whether or not the matrix is close to singularity. (For LAPACK routine DGECON returns an estimate of the (reciprocal of the) condition number.)

Best wishes,

Sven Hammarling.