applied_numerics_uulm/w3/solveLR.m

30 lines
514 B
Matlab

function x_ret = solveLR(A, b)
len = length(b)
z = zeros(len,1)
% solve Lz = b
for i = 1:len
z(i) = b(i);
for j = 1:(i-1)
z(i) = z(i) - z(j) * A(i, j);
end
% no division necessary as L is unipotent
end
x = zeros(len,1)
% solve Rx = z
for i = len:-1:1
x(i) = z(i);
for j = len:-1:(i+1)
x(i) = x(i) - x(j) * A(i, j);
end
x(i) = x(i) / A(i, i)
end
x_ret = x
end