Unexpected behavior of backslash for singular matrix

[Matthias Gobbert]

Dear Colleagues,

The following small linear system has a singular system matrix, as Octave 
confirms readily by det(A) giving an exact 0, but when attempting to solve 
the linear system by the obvious sequence of commands

A = [1 1 0 1; 2 1 -1 1; 4 -1 -2 2; 3 -1 -1 2]
b = [2; 1; 0; -3]
x = A \ b

results in the output below.

Octave gives a warning, but then outputs a 'solution' vector x that has no 
"inf" or "NaN" in it. I find this behavior unexpected. Naturally, 
mathematically one should not proceed after the warning, or one should not 
trust the result after seeing the warning, but if one insists anyway and 
since the behavior is really inside of the backslash operator so that I 
cannot stop the process manually, should x not involve entries of "inf" or 
"NaN"?

I thought that Octave calls LAPACK routines, and these to my knowledge 
would execute the arithmetic no matter what and thus automatically return 
"inf" or "NaN" or at minimum very large or very small numbers. Yet, here, 
Octave gives quite reasonably scaled entries in x that belie the fact that 
the result is completely invalid, see the computation of A*x and A*x-b at 
the end below.

Can you explain what Octave does internally to the backslash operator that 
might explain this situation?

Thanks in advance,

Matthias


> > A = [1 1 0 1; 2 1 -1 1; 4 -1 -2 2; 3 -1 -1 2]
A =
   1   1   0   1
   2   1  -1   1
   4  -1  -2   2
   3  -1  -1   2

> > b = [2; 1; 0; -3]
b =
   2
   1
   0
  -3

> > x = A \ b
warning: matrix singular to machine precision
x =
   0.11111
   1.66667
  -0.88889
  -0.77778

> > A*x
ans =

   1.00000
   2.00000
  -1.00000
  -2.00000

> > A*x - b
ans =
  -1.00000
   1.00000
  -1.00000
   1.00000