function [ L, U ] = LU_decomp( A )

[n,m] = size(A);

L = diag(ones(n,1));
for k = 1:n-1
    for i = k+1:n
        L(i,k)= A(i,k)/A(k,k);
        for j = 1:m
            A(i,j) = A(i,j) - s * A(k,j);
        end
    end
end
U = A;

end

