wrong LU decomposition

[Thomas Weber]

Am Dienstag, den 05.02.2008, 02:19 +0100 schrieb David Bateman:
> miro simko wrote:
> > --------
> > Bug report for Octave 2.1.73 configured for i486-pc-linux-gnu
> > 
> > Description:
> > -----------
> > 
> > octave computes wrong LU decomposition for given matrix A
> > correct decomposition A=L*U:
> > (1   2   3   4)   (1  0  0  0)   (1  2  3  4)
> > (2   9  12  15) = (2  1  0  0) * (0  5  6  7)
> > (3  26  41  49)   (3  4  1  0)   (0  0  8  9)
> > (5  40 107 135)   (5  6  7  1)   (0  0  0 10)
> > 
> > octave gives:
> > L= (0.20000   0.85714  -0.25000   1.00000)
> >    (0.40000   1.00000   0.00000   0.00000)
> >    (0.60000  -0.28571   1.00000   0.00000)
> >    (1.00000   0.00000   0.00000   0.00000)
> > 
> > U=(5.00000    40.00000   107.00000   135.00000)
> >   (0.00000    -7.00000   -30.80000   -39.00000)
> >   (0.00000     0.00000   -32.00000   -43.14286)
> >   (0.00000     0.00000     0.00000    -0.35714)
> > 
> > 
> > Repeat-By:
> > ---------
> > L = [1,0,0,0; 2,1,0,0; 3,4,1,0; 5,6,7,1];
> > U = [1,2,3,4; 0,5,6,7; 0,0,8,9; 0,0,0,10];
> > A = L*U;
> > %A = [1,2,3,4; 2,9,12,15; 3,26,41,49; 5,40,107,135];
> > [L2,U2] = lu(A);
> 
> norm (L2*U2 - A)
> 
> ans =  4.4409e-16
> 
> Nothing here, move along..

I guess lu()'s documentation should be improved. I had to read some
other description to understand why L2 is not a lower diagonal matrix
(and that the returned L matrix changes if you ask for P).

Added to my TODO list.

	Thomas