Davidson Methods

Davidson methods enlarge an eigenvector search space with an approximate correction to the current residual. Their advantage appears when an inexpensive model of $H-\theta I$ captures the difficult directions, as often happens for determinant-basis electronic Hamiltonians. Sparsity alone does not supply that preconditioner.

Helpful background. Rayleigh–Ritz supplies the projected problem; sparse eigensolvers supplies stopping criteria.

Residual correction and the Davidson direction

Let $x$ be a normalized Ritz vector with $\theta=x^\dagger Hx$ and residual $r=(H-\theta I)x$. Seeking $(x+t,\theta+\delta\theta)$ and retaining first-order terms gives

$$(H-\theta I)t\approx-r+\delta\theta x.$$

Ordinary Davidson replaces the difficult shifted operator by a convenient approximation $D-\theta I$, frequently using $D=\operatorname{diag}H$:

$$t_i\approx-\frac{r_i}{D_{ii}-\theta}.$$

The direction is orthogonalized against the existing search basis, normalized, and appended. The expanded projected Hamiltonian is diagonalized again. Cache Hamiltonian actions on retained basis vectors to avoid repeating expensive contractions.

The denominator formula is a model, not an exact solution. Small denominators require a documented regularization, such as a level shift or a bounded inverse. Changing that regularization changes convergence behavior. A nearly zero direction after orthogonalization means that the proposed correction contains little new information; adding a numerically dependent vector makes the projected problem ill-conditioned.

Why an exact unprojected inverse can fail

Since $r=(H-\theta I)x$, solving the unprojected correction equation exactly gives $t=-x$ whenever the inverse exists. Orthogonalization then removes the entire direction. This explains why improving an ordinary Davidson inverse indiscriminately need not improve the method.

Jacobi–Davidson instead solves in the tangent space orthogonal to $x$:

$$P(H-\theta I)Pt=-r,\qquad P=I-xx^\dagger,\qquad x^\dagger t=0.$$

The projection removes the redundant normalization direction. The equation is usually solved approximately; early outer iterations rarely justify a highly accurate inner solve. Near a converged isolated eigenpair, tightening the correction tolerance can accelerate local convergence. For generalized eigenproblems, the metric and projectors must be adapted to the overlap matrix rather than copied from the Euclidean formula.

Block targeting and a two-state example

For

$$H=\begin{pmatrix}0&v\\v&\Delta\end{pmatrix}, \qquad 0<|v|\ll\Delta,$$

with real $v$, start with $x=(1,0)^T$. Then $\theta=0$ and $r=(0,v)^T$. The useful diagonal-preconditioned component is $t_2=-v/\Delta$, matching first-order mixing of the ground state. The first denominator is zero, but its residual component is also zero; the implementation must still avoid evaluating an undefined $0/0$. Once the second independent direction is included, the projected two-dimensional problem is exact.

With several nearby states, expand a block of unconverged residuals and orthogonalize the block against both active and locked vectors. This helps preserve degeneracies and avoids repeatedly converging to one member of a cluster. A restart should retain the desired vectors and enough nearby directions to prevent immediate stagnation.

Convergence and limitations

Monitor original residuals, orthogonality, subspace dimension, and the cost of preconditioner applications. Compare methods at fixed residual accuracy, not only by iteration count: one preconditioned iteration can cost far more than one Lanczos step. Test a deliberately weak preconditioner and a diagonal benchmark to expose denominator handling.

Davidson is especially effective when the chosen diagonal reflects the spectrum. It can perform poorly when strong off-diagonal couplings determine low-energy structure. A successful eigenvalue solve also leaves the determinant-space truncation error untouched; that belongs to configuration interaction.

Exercises

Expand the lower eigenvalue of the two-state matrix through order $v^2$ and explain how the correction direction anticipates it.

Solution

$E_-=(\Delta-\sqrt{\Delta^2+4v^2})/2=-v^2/\Delta+O(v^4/\Delta^3)$. The trial vector $(1,-v/\Delta)^T$, when normalized, has Rayleigh quotient $-v^2/\Delta$ to this order. The preconditioner captures the leading virtual admixture.

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 — preconditioned subspace methods.
• Sleijpen, G. L. G., and H. A. van der Vorst, “A Jacobi–Davidson iteration method for linear eigenvalue problems,” SIAM Journal on Matrix Analysis and Applications 17, 401–425 (1996), doi:10.1137/S0895479894270427 — projected correction equations.