function [L,R] = lr(A) % calculates an LR-decomposition of a matrix A without % pivoting, i.e. returns a lower and a upper triangular % matrix L and R, respectively, such that A = L*R. n = size(A,1); for k = 1:n A(k+1:n,k) = A(k+1:n,k) / A(k,k); A(k+1:n,k+1:n) = A(k+1:n,k+1:n) - A(k+1:n,k) * A(k,k+1:n); end R = triu(A); L = tril(A,-1) + diag(ones(n,1)); end