5
$\begingroup$

Suppose I have a symmetric $N \times N$ matrix A which has a one-dimensional Nullspace $N$. A is positive definite on $N^\bot$. In my case $N$ is the space of constant vectors (i.e. generated by the all-one vector).

I have to solve the problem $Ax = b$, with $b \in R(A)$ which has infinitely many solutions. I am looking for the minimum norm solution. The matrix $A$ is very large and sparse, direct methods aren't really an option. The rank-deficient least squares algorithms I have seen also appear to be prohibitive.

I was solving a non-singular version of this problem with Conjugate Gradient. Is there anyway I can modify the algorithm to solve this particular problem?

EDIT: Boiling the problem down to the bone the question is if $A$ is positive semi-definite with exactly one 0 eigenvalue, does CG work?

Thanks.

$\endgroup$
1
7
$\begingroup$

I'd suggest you to shift away the singularity: solve $(A+ee^T)y=b$ instead and then orthogonalize $y$ with respect to $e$ to get $x$. $A+ee^T$ is not sparse but you can compute matrix-vector products cheaply, and that's all you need for CG.

EDIT: forgot to define it, $e$ is the vector of all ones

$\endgroup$
1
  • $\begingroup$ This is exactly what I did. With the help of a preconditioner it converges sufficiently fast. In hindsight it was the obvious route to take. Thanks for pointing it out to me. $\endgroup$ – RadonNikodym Aug 16 '10 at 8:21
3
$\begingroup$

Plain CG, without the regularization, will converge to the pseudo-inverse solution provided you start the iteration with $y=0$. The regularized (preconditioned) solution will be different and the operation count will depend on the sparsity, so it's good to maintain this if you can. I guess partly it must depend on whether you can tolerate a solution that has a nonzero component in the null space.

$\endgroup$
0
$\begingroup$

You can add a linear constraint, stipulating that your solution should be orthogonal to the null space, assuming it is known to you ahead of time.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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