function [L U]=LUCrout(a)
% Factorización LU de matriz a por Crout
%
 n=size(a,1); L=zeros(n); U=eye(n);
 for k=1:n
    i=k:n;   L(i,k)=a(i,k)-L(i,1:k-1)*U(1:k-1,k);
    i=k+1:n; U(k,i)=(a(k,i)-L(k,1:k-1)*U(1:k-1,i))/L(k,k);
 end
end
