function ShowChol()
% Illustrates how the Cholesky factorization can be computed and how
% to use the Matlab chol function
format short
clc
% Generate a small symmetric poistive definite matrix...
n = 5;
X = randn(n,n);
A = X'*X;
% Compute the lower triangular Cholesky factor...
G = MyChol(A) % A = G*G' G lower triangular
L = chol(A,'lower') % A = L*L' L lower triangular
R = chol(A) % A = R'*R R upper triangular
% Note, to solve Ax = b we would proceed as follows...
% x = G'\(G\b)
% x = L'\(L\b)
% x = R\(R'\b)

function G = MyChol(A)
% A is nxn symmetric and positive definite.
% G is lower triangular so A = GG'
[n,n] = size(A);
G = zeros(n,n);
v = zeros(n,1);
for k=1:n
% Compute the kth column of G assuming that you know G(:,1:k-1).

% Compare kth columns in A = GG':
% A(:,k) = G(:,1:k)*G(k,1:k)'
% = G(:,1:k-1)*G(k,1:k-1)' + G(:,k)*G(k,k)

% Thus, if v = A(:,k) - G(:,1:k-1)*G(k,1:k-1)' (computable) then
% v = G(:,k)*G(k,k).
% This shows that we can compute G(:,k) once we have v.
% But G is lower triangular so all we care about is v(k:n)...

v(k:n) = A(k:n,k) - G(k:n,1:k-1)*G(k,1:k-1)';

% Since v(k:n) = G(k:n,k)*G(k,k), it follows that

G(k:n,k) = v(k:n)/sqrt(v(k));
end

% n^3/3 flops = half of LU.