# Preconditioner for finding the smallest eigenpairs of a large, but structured, matrix

I'm trying to find the eigenvector corresponding to the second smallest eigenvalue of a large $(4,000,000 \times 4,000,000)$ matrix $L$. $L$ is a graph Laplacian, with the following structure: $L = D - SBS^\intercal$, where $D$ is a diagonal matrix and $S$ and $B$ are a large sparse matrices. $S$ has dimension $4,000,000 \times 10,000,000$, but only about $40,000,000$ non-zero entries. $B$ has similar scale and sparsity. So I can rapidly perform matrix-vector multiplication: $Lv = Dv - S(B(S^\intercal v))$.

Currently I'm using Scipy (which calls an Arnoldi algorithm implemented in ARPACK) to find the smallest eigenvalues and corresponding eigenvectors of $L$. Rather than directly finding the smallest eigenpairs of $L$, I find the largest eigenpairs of $M^{-1}$, where $M = L + cI$. (Adding $cI$ and inverting doesn't change the eigenvectors.) To compute $M^{-1}v$, where $v$ is an arbitrary vector, I use the conjugate gradient algorithm.

But this approach is slow: the Arnoldi algorithm often takes hundreds of iterations to converge, and each Arnoldi iteration requires a few hundred conjugate gradient (CG) iterations. Since I'm solving hundreds of CG problems with the same design matrix $M$, I was thinking there could be much efficiency gained by preconditioning. But what preconditioner would be suitable for this problem? I've read that provably good sparse approximations of Laplacians exists ("ultrasparsifiers"), but can they be shifted and inverted, without destroying their sparsity, to form good approximations of $M^{-1}$? Or are there techniques for preconditioning that exploit the structure of $M$?

Thanks much.

• Have you looked at scicomp.stackexchange.com ? They might have more expertise. – S. Carnahan Sep 27 '13 at 3:43
• I got a promising answer from there: scicomp.stackexchange.com/questions/8621 – Jeff Sep 30 '13 at 2:20
• Well, you can apply the Arnoldi iteration to $L^{-1}$ where each solve with $L$ you approximate by an inner iteration with some efficient preconditioning, e.g., a multigrid one as provided in the linked answer. The AMG package PyAMG (pyamg.org) is built on SciPy and NumPy. – Algebraic Pavel Oct 13 '13 at 22:22