# Solving extended positive definite system with Cholesky factorization

Suppose we have a system of equations:

$$\begin{bmatrix} \boldsymbol{A} & \boldsymbol{1} \\ \boldsymbol{1}^\top & 0 \end{bmatrix} \begin{bmatrix} \boldsymbol{x} \\ x_0 \end{bmatrix} = \begin{bmatrix} \boldsymbol{b} \\ 1 \end{bmatrix}$$

where $\boldsymbol{A} \succ 0$ is symmetric positive definitive. How would you exploit the Cholesky factorization of $\boldsymbol{A} = \boldsymbol{L}\boldsymbol{L}^\top$ to efficiently solve this system for various right-hand-side $\boldsymbol{b}$?

If there is no efficient trick with the Cholesky factorization, what would you recommend to guarantee decent numerical stability and speed?

We have $$\mathbf{A}\mathbf{x}+x_0\mathbf{1}=\mathbf{b}, \quad \mathbf{1}^T\mathbf{x}=1.$$ Eliminate $\mathbf{x}$ from the first equation and substitute to the second: $$\mathbf{x}=\mathbf{A}^{-1}(\mathbf{b}-x_0\mathbf{1}) \quad\implies\quad \mathbf{1}^T\mathbf{A}^{-1}(\mathbf{b}-x_0\mathbf{1})=1.$$ This gives $$x_0=\frac{\mathbf{1}^T\mathbf{A}^{-1}\mathbf{b}-1}{\mathbf{1}^T\mathbf{A}^{-1}\mathbf{1}}.$$
So you can first solve using Cholesky the systems $$\mathbf{y}=\mathbf{A}^{-1}\mathbf{b}, \quad \mathbf{f}=\mathbf{A}^{-1}\mathbf{1},$$ and then form $x_0$ and $\mathbf{x}$ as $$x_0=\frac{\mathbf{1}^T\mathbf{y}-1}{\mathbf{1}^T\mathbf{f}}, \quad \mathbf{x}=\mathbf{y}-x_0\mathbf{f}.$$