%Cholesky factorization of a symmetric positive definite matrix A

A = [5,7,6,5;7,10,8,7;6,8,10,9;5,7,9,10]

[n,n] = size(A);      %A is nxn (A must be square)
L = zeros(n,n);

for k=1:n
    for j=1:k-1
        L(k,k) = L(k,k)-L(k,j)^2;
    end
    L(k,k) = sqrt(L(k,k)+A(k,k));
    for i=k+1:n
        for j=1:k-1
            L(i,k) = L(i,k)-L(i,j)*L(k,j);
        end
        L(i,k) = (L(i,k)+A(i,k))/L(k,k);
    end
end

disp('the L matrix is...')
L;
