Cholesky factorisation (thing)

The Cholesky Factorisation of a matrix A is given by a matrix L such that LL^T=A, where L is lower triangular.

This factorisation is a special case of LU factorisation which is possible if and only if the original matrix was symmetric positive definite: that is, A=A^T (the symmetric part) and for any vector x other than the zero vector, x^TAx is strictly greater than zero (the positive definite part).

Once such a factorisation has been found, the original system Ax=b can be solved as for LU factorisation, with L^T taking the place of U. Thus the only difference is making use of the symmetry to simplify the derivation.

Deriving a Cholesky Factorisation
As with LU factorisation, the matrix L is built up by a series of column vectors which are obtained by an iterative process using the matrix A and previously determined columns. There is no longer any need to track U, since this will simply be L^T; however, additional computation is introduced by the use of square roots.
So to start the process, define L₁ as the first column of A divided through by sqrt(a_1,1).
Construct a new matrix A² by A² = A -L₁L₁^T (again, this is not a scalar product, else the difference makes no sense). Obtain L₂ from the second column of A² divided by sqrt(a²_2,2).
Proceed iteratively, finding L_i by Aⁱ_i/sqrt(aⁱ_i,i), until all N columns have been found. Then check LL^T gives A.

Possible problems
If any of the diagonal entries of A is zero or negative, then dividing by the square root won't be possible (working in the reals). However, by the iff relation between being symmetric positive definite and having a Cholesky factorisation, this should not occur. In fact, this can be used as a test for A being SPD, rather than the potentially awkward x^TAx > 0: attempt to find the Cholesky factors and if the process ever fails, you didn't have an SPD matrix.

---

A sample Cholesky factorisation
Let A be the matrix

 2 -1  0
-1  2 -1
 0 -1  2

Then L₁ is A₁/sqrt(a_1,1)= [2, -1, 0]^T/sqrt(2). It's easiest to carry the square root outside the vector, as it'll be squared back up again when we construct A²:

              2 -1  0    2  -1  0    0  0  0
A - L1L1T =   -1 2 -1 - -1  0.5  0 = 0 1.5 -1
              0 -1  2    0   0  0    0  -1  2

The second column of A² divided by sqrt(a²_2,2)=sqrt(3/2) gives us L₂= [0 3/2 -1]^T/sqrt(3/2). Repeat again to get A³ which is just

0 0 0
0 0 0
0 0 4/3

Which makes L₃=[0 0 4/3]^T/sqrt(4/3)= [0 0 2/√3]^T. Combining these columns we get for L:

2/√2    0        0
-1/√2  √(3/2)   0
0     -√(2/3)  2/√3