Let's say we have a block matrix $ M =\left( \begin{array}{ccc}
A & B\\
B^{*} & C \end{array} \right)$ where $M$ is positive definite. ($A$ and $C$ are also positive definite.)
There is a formula for carrying out block Cholesky decomposition. See Wikipedia: Block LU decomposition. Summarising we have the following result.
The matrix $M = LU$ can be decomposed in an algebraic manner into
$$L =
\begin{pmatrix}
A^{\frac{1}{2}} & 0 \\
B^{*} A^{-\frac{*}{2}} & Q^{\frac{1}{2}}
\end{pmatrix}$$where $\begin{matrix}
Q = C – B^{*} A^{-1} B
\end{matrix}$ ($*$ indicates transpose in this case)
Now lets say we have already carried out the Cholesky decomposition for A, and C. So we have already calculated $A^{1/2}$, and $C^{1/2}$ (It is therefore straightforward to calculate the inverses $A^{-1/2}$, and $C^{-1/2}$ using forward substitution).
Rewriting $Q$ in terms of these quantities we now have
$$Q = Q^{1/2}Q^{*/2} = C^{1/2} C^{*/2} – (B^{*} A^{-*/2})(A^{-1/2} B)$ = $(C^{1/2} + B^{*}A^{-*/2})(C^{1/2} – B^{*}A^{-*/2})^{*}.$$
My question is this: Given this set up is it possible to algebraically calculate $Q^{1/2}$ without having to apply Cholesky decomposition to $Q$. Or, in other words, can I use $C^{1/2}$ to help me in the calculation of $Q^{1/2}$.
Best Answer
If A,C are fixed, and B is variable but nice (low-rank), then you want what is called "Cholesky update". If A,B,C are fixed, then probably you should not be picky about how the blocking is done, and you want to use a standard "block Cholesky". I have not found a clear answer for A,C fixed, B variable and not nice (I can ask around, so let me know if that really is your case).
Cholesky update
Rank one updates, chol(A) to chol(A+xx*), are easy and safe. Rank one "downdates", chol(A) to chol(A-xx*), are easy but require a little care: stable algorithms are given in Stewart's Matrix Algorithms Vol 1, Algorithm 4.3.8, p. 347. Chapter 12.5 of Golub–Van Loan has some similar stuff, and Cholesky down-dating in 12.5.4. This function has been widely implemented, and the cholupdate command in matlab dates back to 1979 code from LINPACK.
[0,B*;B,0] is a sum of rank one matrices, and so by updating and downdating those rank one guys, you could probably get what you want, and it might even be faster than chol(Q). However, it can be a lot better to update more ranks at a time.
Apparently this is a common request in machine learning, and M. Seeger wrote a technical report on this problem of low rank updates to a Cholesky factorization, and mentions several common pitfalls, especially as regards to actually doing it with existing software.
A more scholarly (and older) treatment is in section 3 of this article version of Ch. 12.5 of GvL:
Gill, P. E.; Golub, G. H.; Murray, W.; Saunders, M. A. "Methods for modifying matrix factorizations." Math. Comp. 28 (1974), 505–535. MR343558 DOI:10.2307/2005923
Davis and Hager in MR1824053 note that algorithm C1 can be used for a reasonably efficient, multiple rank, single pass, update of a dense matrix (and go on to describe sparse techniques).
Note that these mostly do not take advantage of the block structure of [0,B*;B,0], so you might find something better that is more specialized.
Block Cholesky
Blocking the Cholesky decomposition is often done for an arbitrary (symmetric positive definite) matrix. I didn't immediately find a textbook treatment, but the description of the algorithm used in PLAPACK is simple and standard.
In their algorithm they do not use the factorization of C, just of A. That allows them to reduce the problem of chol([A,B*;B,C]) to just chol(A) and chol(Q). The point of the algorithm is that you do not choose A and C to have the same size. You choose A to fit nicely in cache, and do your work at a higher BLAS level. In other words, A is a real block, and C is just leftover garbage you'll need to sweep up next.
In particular, C is discarded and replaced by Q during the algorithm, but chol(Q) is also computed by decomposing Q itself into a block matrix. This means that the algorithm is discarding any information you had about C, so if C is fixed while B varies, this would be quite wasteful.