Factor LU de matriz-- en Matlab

function [L,U]= factorlu (A)
n=length(A);
U=A;
L=eye(n);
for p=1:n
    for q = p+1:n
        L(q,p)=U(q,p)/U(p,p);
        U(q,p:n)=U(q,p:n)-U(p,p:n)*L(q,p);
    end
end

ejecutando la funcion

>> A = [1 2 2; 4 4 2; 2 6 4];
>> [L,U]=factorlu(A)

L =

    1.0000         0         0
    4.0000    1.0000         0
    2.0000   -0.5000    1.0000


U =

     1     2     2
     0    -4    -6
     0     0    -3