# What happens when I use a conjugate gradient solver with a symmetric positive semi-definite matrix?

I have a symmetric positive semi-definite matrix, i.e., a laplacian and wonder what may happen when I use a CG solver, that is an algorithm for positive definite matrices.

• What happens when the matrix is only positive semi-definite?

• Can I avoid problems by somehow substracting the eigenvector for the eigenvalue 0?

For the laplacian I know the eigenvector(s) for eigenvalue 0, because they are the indicator vectors for the connected components of the domain.

I know there are some numeric tricks like adding a small value to the main diagonal, but I guess this will either make the solution worse when the value is too large or the algorithm will still be instable when the value is too small.

• If a matrix has an eigenvalue $0$, then $\text{det}(A)=0$ and we are dealing with a singular matrix. Solving a system of linear equations with a singular matrix poses problems of its own, independent of whether you use CG or ILU to solve it. May 28, 2020 at 6:54
• @ThijsSteel Oh, you're right. So for most applications my question should probably be what boundary conditions make a system solvable. But I was still interested in what will actually happen when trying to use a CG solver for such a system. Another follow up question would probably be what kind of solutions the "more general" solvers provide in that case.
– allo
May 28, 2020 at 11:02

You first have to make sure that your system is solvable: this happens iff the right-hand side $$b$$ is orthogonal to the kernel of $$A$$. If $$A$$ has a dimension-1 kernel spanned by $$v$$, you need to have $$v^* b=0$$. If that is not the case, then go back to the modelling stage and ask yourself if what you did makes sense.

Convergence of iterative methods for singular systems is tricky in general, but for positive-semidefinite matrices and CG all works nicely: CG will converge on your singular system without any modifications. More precisely, it will converge to the unique solution orthogonal to the kernel ($$v^Tx=0$$), which is also the minimum-norm one.

For a proof, see https://arxiv.org/abs/1809.00793 . Intuitively, what happens is that, after a suitable orthogonal change of basis, you are solving the system $$\begin{bmatrix} A_{11} & 0\\ 0 & 0 \end{bmatrix} \begin{bmatrix} x_1\\ x_2 \end{bmatrix} = \begin{bmatrix} b_1\\0 \end{bmatrix},$$ with $$A_{11}$$ positive definite. In this form, it is not complicated to see that all the iterates that CG produces have a zero in their second block: $$(x_k)_2=0$$, for each $$k$$. In addition, the iterates coincide with the ones that you'd get by applying CG to the (nonsingular) linear system $$A_{11}x_1=b_1$$.

Numerically, you may want to re-orthogonalize $$x_k \leftarrow x_k - vv^Tx_k$$ every few steps to make sure that your iterates stay orthogonal to the kernel.

• If it's not orthogonal to the kernel, the residual of $b-Ax$ simply converges to the component of $b$ spanned by the kernel, right? Apr 10 at 5:58
• @syockit This does not sound immediate to prove to me at the moment (but I might be just missing things, I did not think thoroughly about this issue). If $b_2$ is not zero, all iterates $x_k$ of CG have a nonzero second block $(x_k)_2$, and this block may change the scalar products and affect the convergence behavior of CG. Apr 10 at 7:58

To at least partially answer the question:

You have to make sure the system is solvable in the first place. As mentined in the accepted answer, "the system $$Ax=b$$ is solvable iff $$b \bot Ker(A)$$", but that was never obvious to me.

The obvious condition of solvability is $$b \in Im(A)$$. Next, it holds, that $$Im(A)$$ and $$Ker(A^\top)$$ are each other's orthogonal complements (for any matrix $$A$$). And since $$A$$ is symmetric in our case, $$Ker(A^\top) = Ker(A)$$.