Leo, thanks for your reply.
Thanks for your reply too professor. :)
I was able to successfully compute the LU decomposition that way.
Regards,
Nuno
Jaime Villate wrote:
> 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
>
>
>
>
>
>