I need to write a program to solve matrix equations Ax=b where A is an nxn matrix, and b is a vector with n entries using LU decomposition. Unfortunately I'm not allowed to use any prewritten codes in Matlab. I am having problems with the first part of my code where i decompose the matrix in to an upper and lower matrix.
I have written the following code, and cannot work out why it is giving me all zeros on the diagonal in my lower matrix. Any improvements would be greatly appreciated.
function[L R]=LR2(A)
%Decomposition of Matrix AA: A = L R
z=size(A,1);
L=zeros(z,z);
R=zeros(z,z);
for i=1:z
% Finding L
for k=1:i-1
L(i,k)=A(i,k);
for j=1:k-1
L(i,k)= L(i,k)-L(i,j)*R(j,k);
end
L(i,k) = L(i,k)/R(k,k);
end
% Finding R
for k=i:z
R(i,k) = A(i,k);
for j=1:i-1
R(i,k)= R(i,k)-L(i,j)*R(j,k);
end
end
end
R
L
end
I know that i could simply assign all diagonal components of L to be 1, but would like to understand what the problem is with my program!
I am also wondering how to change this program to include pivoting. I understand I need to say that if a diagonal element is equal to zero something needs to be changed. How would I go about this? Thanks in advance!