Sparse Eigensolvers

Sparse eigensolvers compute a small spectral window without storing all eigenvectors or reducing a large matrix densely. Their central operation is usually $v\mapsto Hv$; sparsity matters through the cost of this action, not through a guarantee of fast convergence. A request for the lowest energies is distinct from a request for eigenvalues of smallest magnitude.

Helpful background. Matrix diagonalization defines residuals and the sensitivity of nearly degenerate eigenspaces.

Spectral targets and operator actions

For Hermitian $H$, extremal methods approximate eigenpairs from a trial subspace $Q$ by diagonalizing $Q^\dagger HQ$. The columns of $Q$ should remain orthonormal. Lanczos builds this space through repeated Hamiltonian actions, while Davidson methods use preconditioned residual directions.

A matrix-free implementation should pass linearity and adjoint tests before it is given to either method. For arbitrary test vectors $u,v$, compare $u^\dagger Hv$ with $(Hu)^\dagger v$. Check each symmetry sector separately. A callback that changes its approximation tolerance between actions may no longer represent one fixed Hermitian operator.

Memory includes the search basis, often $O(nm)$ for $m$ stored vectors, in addition to sparse matrix storage. Restarting reduces this cost but can discard useful directions. Clustered eigenvalues commonly require a larger retained subspace or a block method.

These Hermitian algorithms are not appropriate merely because an operator is sparse. Absorbing potentials, effective decay Hamiltonians, and Liouvillians can be genuinely non-Hermitian. Arnoldi iteration then supplies the general Krylov projection

$$AQ_m=Q_m B_m+h_{m+1,m}q_{m+1}e_m^\dagger,$$

where $B_m=Q_m^\dagger A Q_m$ is upper Hessenberg. Orthogonalization generally involves all preceding basis vectors, instead of the exact-arithmetic three-term Lanczos recurrence. For a nonnormal operator, small residuals alone do not give the Hermitian spectral-distance bound: left/right eigenvector conditioning or pseudospectral information may also matter. Preserve the physical non-Hermitian part rather than symmetrizing it to fit a Hermitian solver.

Shift-invert targeting of interior levels

To find eigenvalues near a real shift $\sigma$, use

$$B=(H-\sigma I)^{-1},\qquad Bu_j=\mu_j u_j,\qquad \mu_j=\frac1{E_j-\sigma}.$$

Levels nearest $\sigma$ become large in magnitude. The inverse is an operation: factor $H-\sigma I$ once and solve linear systems, or use controlled iterative solves. Do not form a dense inverse. Sparse factorization can create substantial fill-in, so a matrix with few nonzero entries can still exceed memory during shift-invert.

A shift exactly on an eigenvalue makes the system singular. A nearby shift improves spectral targeting but worsens conditioning of the solves. The outer algorithm must account for inner errors, particularly if it relies on self-adjointness. Recompute every accepted residual with the original Hamiltonian,

$$r_j=Hx_j-\theta_jx_j, \qquad \theta_j=x_j^\dagger Hx_j.$$

The transformed residual alone can be misleading. If $Bx=\mu x+s$, then $Hx-(\sigma+1/\mu)x=-(H-\sigma I)s/\mu$. Thus the size of $s$ has to be interpreted through a further operator and scale factor.

A target-selection example

Consider $H=\operatorname{diag}(-10,-1,0.2,3)$. The algebraically smallest energy is $-10$, whereas the smallest-magnitude energy is $0.2$. With $\sigma=0$, the inverse eigenvalues are $-0.1,-1,5,1/3$; the largest magnitude identifies $0.2$. With $\sigma=-0.8$, the level $-1$ becomes the closest target.

This example is deliberately diagonal: it isolates target semantics from iteration. A second test should rotate it by a random unitary matrix while preserving eigenvalues. The solver must recover the same spectrum, although iteration counts can depend on the start vector and preconditioner.

Convergence and missing states

Accept a block only when original residuals meet the energy scale of interest and $X^\dagger X$ is close to $I$. Repeat with a different starting vector or larger search space when completeness of a narrow window matters. A start vector exactly orthogonal to a symmetry sector never finds states in that sector through polynomial actions alone.

A converged Ritz value does not establish that every lower eigenvalue was found. Check known symmetry multiplicities, compare level counts with small dense cases, and search overlapping spectral windows. Separate algebraic residual error from basis truncation and physical finite-size error; tightening the solver cannot remove either of the latter.

Track the observables used in the final claim as well as energies: a stable energy can coexist with a less stable eigenvector. For Hermitian problems, the Rayleigh-quotient residual equals the energy standard deviation. The same residual interpretation and cluster-projector checks apply to dense and sparse solvers.

Exercises

For the example spectrum, choose a shift that targets $3$ with a largest-magnitude inverse solver and explain why choosing $\sigma=3$ is unsuitable.

Solution

A shift such as $2.8$ makes $1/(3-2.8)=5$ dominate the other inverse magnitudes. At $\sigma=3$ the shifted matrix is singular, so its inverse does not exist. A small offset balances selectivity against conditioning.

References

• Bai, Z., J. Demmel, J. Dongarra, A. Ruhe, and H. van der Vorst, eds., Templates for the Solution of Algebraic Eigenvalue Problems: A Practical Guide, SIAM, 2000, online text — spectral transformations and iterative eigenproblems.
• Lehoucq, R. B., D. C. Sorensen, and C. Yang, ARPACK Users’ Guide: Solution of Large-Scale Eigenvalue Problems with Implicitly Restarted Arnoldi Methods, SIAM, 1998 — restart and target-selection conventions.
• Saad, Y., Numerical Methods for Large Eigenvalue Problems, revised ed., SIAM, 2011, author’s text — Krylov projection, Arnoldi iteration, and sparse eigenvalue algorithms.