# Null Space Projection for Singular Systems

Let $A$ be a general symmetric matrix operator and $P$ be the unique orthogonal projection onto $\operatorname{Range}(A) = \operatorname{Null}(A)^\perp$.

Analytically, the system $$Ax = Pb$$ should have a solution for any vector $b$. If we further stipulate that this solution lie in $\operatorname{Range}(A)$, the solution should be unique. (And in the positive definite case, this is precisely the solution that the conjugate gradient algorithm will find in exact arithmetic.)

The above is a common tactic when solving the Poisson equation with either Neumann or periodic boundary conditions; if $A$ is taken to be the Neumann or periodic Laplacian, then $$P = I - \frac{1}{N}(1_N)(1_N)^T,$$ where $1_N = [1, \cdots, 1]^T$, and we can solve $Ax = Pb$ using a conjugate gradient approach.

But what if the projection operator $P$ doesn't have a simple analytic expression? Can it be calculated numerically from $A$? Is there a way to solve this type of singular problem in the case that $A$ is a black-box symmetric operator?

I realize one approach would be to solve the normal equations $A^2x = Ab$, which avoids using a projection operator entirely, but this also severely worsens the conditioning of the system, so I'm interested to know if there's another way.

Jed has already answered your question, implicitly assuming that your matrix is large and sparse. If this should not be the case, you can do an eigenvalue decomposition and write $A=X\Lambda X^T$ where the columns of $X$ correspond to the eigenvectors and the entries of the diagonal matrix $\Lambda$ are the eigenvalues. (Note that this decomposition is not easily computable if $A$ is large, as Jed mentioned.) Then, the null space consists of exactly those vectors spanned by the eigenvectors corresponding to zero (or very small) eigenvalues, i.e., the projector onto the null space can be written as $$1_N = X \tilde\Lambda X^T$$ where $(\tilde\Lambda)_{ii} = 0$ if $|\Lambda_{ii}|\ge \varepsilon$ and $(\tilde\Lambda)_{ii} = \Lambda_{ii}$ if $|\Lambda_{ii}|<\varepsilon$.
• @GeoffOxberry: Yes, you are correct. I meant to add that (the above is the conceptual view, not how it's implemented) but forgot. That said, if $A$ is a symmetric matrix, isn't the SVD the same as the eigenvalue decomposition, and isn't QR just one way of computing it? Jun 2, 2013 at 17:30
• Yes. All of the above also apply to the general case, although of course the SVD and eigendecomposition are not necessarily equal if $A$ is not Hermitian. Although SVD can be calculated using QR (and vice versa), I've been told the SVD is the more accurate algorithm. Jun 2, 2013 at 22:23
• To follow up, I looked this up, and the SVD is only the same as the eigenvalue decomposition if $A$ is Hermitian and positive semidefinite. If $A$ is Hermitian, but not positive semidefinite, it still has orthogonal eigenvectors, which can be used to calculate a projector. Jun 4, 2013 at 5:15