Does Lanczos have trouble with large matrix elements?

Question

I have a large, yet very sparse, matrix that I'd like to diagonalize. Both my own Lanczos implementation and the ARPACK that's built in with scipy fail to converge properly, though. I know that my Lanczos code itself is correct, and the ARPACK implementation should be correct too.

One thing that I note about my matrix is that the entries in the $n$-th row grow like $\mathcal{O}(n^2)$. Could that contribute to the problem? Especially since - for physical reasons - I expect the eigenvector associated with the smallest (algebraic) eigenvalue to have little weight at these indices.

EDIT: If I start with a totally not random starting vector that just has weight 1 in the first element and is zero everywhere else, then both my own Lanczos and ARPACK converge to the correct (LAPACK dense matrix diagonalization) solution. In a general setting where I wouldn't already have some physical intuition what eigenvector I'm looking for, what can one do to combat large matrix elements (both diagonal and off-diagonal)?

Answer 1

It sounds like for your general case, you should be doing some kind of column and row scaling, effectively left and right preconditioning with diagonal scaling matrices. You would computing these scaling matrices to try to get all matrix elements on the same order of magnitude.

Example Let's work with an example $2\times 2$ matrix $$ \begin{bmatrix} 2 & 1 \\ 10^{10} & 2\times 10^{10} \end{bmatrix}$$ It has eigenvectors (approximately) $[5\times 10^{-11},1]^T$ and $[-2,1]^T$. Suppose we scale it as $$ \begin{bmatrix}1 & 0 \\ 0 & 10^{-5} \end{bmatrix}\begin{bmatrix} 2 & 1 \\ 10^{10} & 2\times 10^{10} \end{bmatrix} \begin{bmatrix}1 & 0 \\ 0 & 10^{5} \end{bmatrix}$$ Obviously the eigenvalues don't change, but the eigenvectors are now $[5\times 10^{-6},1]^T$ and $[-200000,1]$. Notice that the elements of the first eigenvector are within 6 orders of magnitude of each other, instead of 11. Of course, the second eigenvector has worse scaling, but no worse than the first. This scaling should improve things overall since the accuracy of the small elements of the eigenvectors depends on how they small they are relative to the largest elements.

Edit for symmetry In the case that $A$ is symmetric, then the original eigenvalue problem is $Ax=\lambda x$. If we apply symmetric scaling matrices $S$ on both sides to balance $A$, then $$ SASy = \lambda S^2 y$$ for $y=S^{-1}x$. This preserves symmetry at the cost of having to solve a generalized eigenvalue problem (but the $B$ matrix is simply diagonal, so this may not be too bad).

Note that $S$ should be chosen to have nice round elements, like exact powers of 2. This makes it so that $S$ can be applied exactly without rounding error, and $S^2$ can be stored accurately so long as it is within the floating point range.