# [Maxima] Khaletsky method (LU decomposition)

Jaime Villate villate at fe.up.pt
Wed Feb 10 03:20:38 CST 2010

```Hi Nuno,
it is nice to read you in this list :)

On Tue, 9 Feb 2010, Nuno Santos wrote:
>  A: matrix( [1, 2, 3], [4, 5, 6], [7, 8, 9] );
>  lu_factor(A)\$
>  get_lu_factors(%);
>
>  and that would give me L and U in %[2], and [%3], respectively.
>
>  However, the LU decomposition, but I am looking for it's:
>      L - lower triangular matrix (not lower triangular matrix with unitary
>  diagonal, which is what get_lu_factors returns)
>    U - as upper triangular matrix with unitary diagonal. (not upper triangular
>  without unitary diagonal, which is what get_lu_factors returns)

What I usually did was to construct a diagonal matrix D with its
diagonal equal to that of U and then rewrite LU = (L.D).(inverse(D).U)

L2=L.D will no longer have ones in the diagonal and U2=inverse(D).U will
have ones in the diagonal. That is done in Maxima in this way:

[P,L,U]: get_lu_factors(lu_factor(A))\$
D:genmatrix(lambda ([i,j], if (i=j) then U[i,j] else 0), 3, 3)\$
L2: L.D;
U2: invert(D).U;

However, in your particular example that method fails because U has a zero
in the diagonal, so inverse(D) does not exist.
A more clever way of doing it is what Leo has proposed:

On Tue, 2010-02-09 at 15:18 +0000, Leo Butler wrote:
>
> I think if you want to change this behaviour, you should compute the
> factorization of transpose(A), then transpose L & U.

(%i69) A: matrix( [1, 2, 3], [4, 5, 6], [7, 8, 9] )\$
(%i70) [P,L,U]: get_lu_factors(lu_factor(transpose(A)))\$
(%i71) L2: transpose(U);
(%o71) matrix([1,0,0],[4,-3,0],[7,-6,0])
(%i72) U2: transpose(L);
(%o72) matrix([1,2,3],[0,1,2],[0,0,1])
(%i73) L2.U2;
(%o73) matrix([1,2,3],[4,5,6],[7,8,9])

Regards,
Jaime Villate

```