MathOverflow will be down for maintenance for approximately 3 hours, starting Monday evening (06/24/2013) at approximately 9:00 PM Eastern time (UTC-4).
0

I'm having problems understanding why the conjugate gradient method breaks down for singular matrices. I've read a good introduction to intuitively understanding the CG method through visualizing the quadratic form, but I'm not really understanding mathematically why this seems to be. Does anyone have a good mathematical explanation or any sources that explain why conjugate gradient iteration doesn't work for singular matrices?

Thanks in advance.

flag
In what sense would you want it to work? Typically, a system won't have any solution at all if the coefficient matrix is singular. – Darsh Ranjan Dec 5 2009 at 22:04
My question was geared more towards any mathematical explanation of why singular matrices don't work. My understanding, from the previously referenced article (pg 5), is that singular matrices have multiple solutions rather than a single solution. – john Dec 5 2009 at 22:30

4 Answers

2

If the matrix is singular then then the images of some vectors may be zero and that could mean some of the steps may end up with division by zero.

link|flag
3

Let my reply here to your comment:

My question was geared more towards any mathematical explanation of why singular matrices don't work. My understanding, from the previously referenced article (pg 5), is that singular matrices have multiple solutions rather than a single solution.

Unfortunately, that's not so. If $A$ is a singular matrix, then for most $b$, the equation $Ax=b$ has no solution. It has a solution if and only if $b$ is in the column space of $A$ (or orthogonal to the kernel of $A$, if $A$ is symmetric), in which case there are infinitely many solutions. This is basic linear algebra.

Let's say $A$ is symmetric positive semidefinite, so we can at least entertain the thought of conjugate gradients. Conjugate gradients can be interpreted as maximizing the quadratic function $$f(x) = b^Tx-\frac12 x^TAx.$$ If $A$ is singular, then this function typically has no maximum value, as I will demonstrate. Since $A$ is singular, there are vectors $v$ such that $Av = 0$. If there is any such vector $v$ that is not orthogonal to $b$, consider substituting $x=tv$ into $f(x)$: $$f(tv) = t(b^Tv) - \frac12t^2(v^TAv) = t(b^Tv) = ct,$$ where $c$ is the nonzero real number $b^Tv$. Obviously, this is unbounded as a function of $t$, so we see that there is no maximum value if $b$ is not orthogonal to the kernel of $A$.

link|flag
1

The conjugate gradient method becomes unstable when the matrix A is singular, i.e., if you compute the output of the k-th iteration for a small change in the initial data, large deviations will occur. A similar large deviation occurs due the computer roundoff error. Even for non singular matrices A that are close to to be singular, for example when one of the eigenvalues of A is smaller by orders of magnitudes than the rest, instability problems can occur: The convergence rate becomes slower, and even convergence can be lost due to roundoff errors.

link|flag
2

A search for "conjugate gradient singular matrix" took me to this question. While the answer is obviously given by the responses, the question can be refined: Can CG still give a working algorithm if the matrix is singular, but behaves as a symmetric positive definite form on a (large) subspace?

A standard example is given by the finite element discretization of the Neumann problem on a simply connected domain. The constant functions are both the kernel and the cokernel of the Laplacian. On functions with vanishing mean, the Laplacian is still a positive definite symmetric operator, and we would like to leverage this structure.

This is non-trivial and best our numerical method is derived from a fully analytic setting, because this might provide us the convergence analysis as well. --- This appraoch is for example elaborated in

On the Finite Element Solution of the Pure Neumann Problem
Pavel Bochev and R. B. Lehoucq
SIAM Review
Vol. 47, No. 1 (Mar., 2005), pp. 50-66
Published by: Society for Industrial and Applied Mathematics
Article Stable URL: http://www.jstor.org/stable/20453601

Apart from this canon standard example for the Laplacian, system matrices with larger kernel appear in numerical methods for the de-Rham-complex, in particular if the domain is topologically non-trivial (Finite Element Exterior Calculus, Discrete Exterior Calculus). Singular system solves are still no standard material for education in computational science. As far as I may dare to give an estimate, there is still much room for a better theory building.

link|flag

Your Answer

Get an OpenID
or

Not the answer you're looking for? Browse other questions tagged or ask your own question.