How do you write MATLAB code for LU factorization when U is the unit matrix instead of L. The upper triangular matrix will have the diagonal of 1s instead of the lower triangular matrix.
2 Answers
You are looking for Crout's Method for LU decomposition.
The Wikipedia article has the following code
function [L, U] = LUdecompCrout(A) [R, C] = size(A); for i = 1:R L(i, 1) = A(i, 1); U(i, i) = 1; end for j = 2:R U(1, j) = A(1, j) / L(1, 1); end for i = 2:R for j = 2:i L(i, j) = A(i, j) - L(i, 1:j - 1) * U(1:j - 1, j); end for j = i + 1:R U(i, j) = (A(i, j) - L(i, 1:i - 1) * U(1:i - 1, j)) / L(i,i); end end end
6 Comments
Ben Grossmann
Interesting! I'm surprised that Matlab doesn't have this built in
John Alexiou
@BenGrossmann - why would it, for all practical purposes Doolittle method is sufficient. Only for academic reasons, it would be included. But the big question is why to choose one over the other.
Ben Grossmann
Right... I suppose there isn't much difference in performance. Really, I'm surprised that the LU method doesn't have something like a
unit_upper option, but I guess the extra transposing isn't too much of a bother.John Alexiou
@BenGrossmann - that would be a really nice option. Better yet to decompose into
A = L D U where both L and U are unit diagonals, and D is only diagonal, so you can choose your flavor afterward, and help in symbolic evaluations.Ben Grossmann
Yes, that is a sensible workaround, but it would be nice to also not have to deal with the D term "by hand". In fact, it looks like there's an option to get
[L,U,P,Q,D] = lu(S) so that L*(P'*D*Q')*U for permutation matrices P,Q, which is a sort of maximally flexible approach |
Here's one way to do this. If M is your matrix,
[P,Q] = lu(M');
L = Q'; U = P';
Alternatively, assuming that M has an LU decomposition, we could do the following:
[L,U,P,Q,D] = lu(M);
L = L*P'*D*Q';
This should yield a lower triangular L and upper triangular U where U has 1's on the diagonal whenever M has an LU decomposition.
1 Comment
John Alexiou
I was about to write the exact same answer since
A' = U' * L' and the Doolittle method makes the L factors with unity in the diagonals.
lu factorization doolittle's vs.and scrolling way down.