Conditioning, Residuals, and Convergence

A residual measures how well a computed object satisfies the equations that were actually discretized. It does not, by itself, measure distance to the desired continuum solution. When an operator is ill-conditioned or nonnormal, a tiny residual can accompany a large forward error. Credible digits therefore come from residuals combined with condition estimates, resolution and precision sequences, endpoint variations, and an independent formulation.

Required background. Ordinary Differential Equations supplies boundary problems, eigenvalues, norms, and numerical error concepts.

Helpful background. Formulations of General Relativity provides the continuum systems whose discrete conditioning is assessed, while Partial Differential Equations and Green Functions supplies operator domains, boundary conditions, and non-self-adjoint spectra.

Backward error is not forward error

For a linear solve $Ax=b$, the residual $r=b-A\hat x$ is a backward-error diagnostic: $\hat x$ exactly solves a nearby right-hand side. The forward error is controlled schematically by

$$\frac{\lVert x-\hat x\rVert}{\lVert x\rVert} \le \kappa(A)\frac{\lVert r\rVert}{\lVert b\rVert},$$

for an invertible $A$, nonzero $b$, and compatible induced norm, with $\kappa(A)=\lVert A\rVert\lVert A^{-1}\rVert$. This inequality is exact: $x-\hat x=A^{-1}r$ and $\lVert b\rVert\le\lVert A\rVert\lVert x\rVert$ imply it directly. Matrix perturbations and eigenvalue errors require their own bounds. For eigenvalues, left and right eigenvectors replace the ordinary matrix condition number. Nonnormality can make them nearly orthogonal and the eigenvalue highly sensitive even when the discrete residual is near machine precision.

Pseudospectra make this distinction geometric: the $\varepsilon$-pseudospectrum contains points that become eigenvalues under perturbations of size $\varepsilon$, or equivalently where the resolvent is large. Definitions and the first diagnostic examples are given in Trefethen and Embree 2005, §§ 1–2, pp. 1–19.

A two-by-two sensitivity benchmark

Consider

$$A=\begin{pmatrix}1&K\\0&1+\delta\end{pmatrix}.$$

Its eigenvalues are exactly $1$ and $1+\delta$, yet the eigenvectors become nearly parallel when $K/\lvert\delta\rvert$ is large. For the eigenvalue one, a right eigenvector is $v=(1,0)^T$ and a left eigenvector is $w=(1,-K/\delta)^T$ for real nonzero $K,\delta$. Its absolute eigenvalue condition number in the Euclidean norm is

$$\kappa_\lambda=\frac{\lVert w\rVert_2\lVert v\rVert_2}{\lvert w^Tv\rvert} =\sqrt{1+(K/\delta)^2}.$$

Perturb the initially zero lower-left entry by an absolute amount $\varepsilon$. The new eigenvalues are exactly

$$\lambda_\pm=1+\frac{\delta}{2} \pm\sqrt{\frac{\delta^2}{4}+K\varepsilon}.$$

With $K=10^8$ and $\delta=10^{-6}$, the perturbation $\varepsilon=10^{-16}$ moves the eigenvalues by about $10^{-4}$, even though its norm relative to $A$ is about $10^{-24}$. Linear sensitivity predicts shifts $K\varepsilon/\delta$ only when $\lvert K\varepsilon\rvert\ll\delta^2$; outside that range use the square-root formula rather than extrapolating the first-order estimate. Compute this sweep in 53-bit, 100-bit, and 200-bit arithmetic and compare eigenvalue errors, normalized residuals, and left–right overlaps. Exact triangular input can yield exact diagonal eigenvalues despite severe sensitivity to other perturbations, so precision changes alone are not the stress test.

The benchmark teaches a sharp lesson: stable printed digits are not inferred from a residual alone. They are the digits unchanged under precision increase and perturbations justified by the discretization error. Scaling the basis can improve the matrix representation without changing the underlying operator; record that scaling so condition estimates remain interpretable.

From an exact boundary problem to a resonance

For the manufactured problem $-u''=\lambda u$ on $(0,\pi)$ with homogeneous Dirichlet boundaries, the exact eigenvalues $n^2$ provide direct forward errors. At orders $N=20,40,80$, track eigenvalue error, matrix residual, off-grid differential residual, boundary residual, coefficient tail, and matrix condition estimate. Increase arithmetic precision independently. If the residual falls while the eigenvalue error stalls, conditioning or operator inconsistency is limiting.

For a non-self-adjoint resonance, no exact spectrum is generally available. Let $M$ denote physical mass and $GM$ its geometrized length. Use the spin-2, $\ell=2$, $n=0$ Schwarzschild gravitational target, common to the Regge–Wheeler and Zerilli parities,

$$GM\,\omega=0.37367168-0.08896232i$$

only as an external benchmark after solving with the declared $e^{-i\omega t}$ convention, harmonic indices, master variable, boundary factors, and spectral sheet. Leaver states the radiative boundary-value problem and continued-fraction spectrum in Leaver 1985, § 2, Eqs. 2–14 and Table 1; the displayed higher-precision value is given with its mode identification in Dias, Godazgar, and Santos 2022, Figure 4 and surrounding text. Repeat $N=40,60,80$, 53-bit, 100-bit, and 200-bit arithmetic, map-scale and endpoint sweeps, and compare with an independent continued-fraction or time-domain result. Jansen’s examples illustrate resolution/precision cross-comparison and filtering in Jansen 2017, § 3, Eqs. 3.1–3.2 and Figures 2–3; they do not supply residual plots or a residual formulation.

Tests of convergence

Begin with at least three resolutions. Estimate the observed order only in a regime where the error decreases regularly; do not fit through a roundoff plateau. Next raise precision at fixed resolution. Then vary domain maps and finite endpoints. Finally change the numerical method. Each axis tests a different possible failure.

Report solution error when an exact solution exists, differential and boundary residuals in physical norms, and stable digits rather than a long raw number. Singular values or a resolvent estimate can expose nonnormal sensitivity. LAPACK’s distinction among balancing, condition estimation, and generalized eigenvalue algorithms is documented in Anderson et al. 1999, Chapter 4 §§ 4.1–4.3.

Agreement between successive resolutions can be false when both lie below the same unresolved boundary layer. Agreement between double- and high-precision results can also be false when discretization dominates. The evidence hierarchy works because neither comparison is treated as sufficient alone.

Reproducible record and stopping rules

Pin the numerical library, eigensolver, balancing and scaling options, basis, domain map, boundary enforcement, norms, arithmetic precision, and matching algorithm across resolutions. A clean rerun emits matrices or their builders, solutions, residuals, singular-value or eigenvector-condition estimates, coefficient tails, stable-digit tables, and independent comparisons in machine-readable form.

The validity regime is the declared operator domain, spectral branch, and resolution range before saturation. Failure modes include spurious modes, endpoint contamination, ill-conditioned matrices, tolerance saturation, branch swapping, and a residual computed in a norm that suppresses the region of largest error. Stop adding digits when any required diagnostic stops improving. A second plot generated from the same matrix and eigenvector is not an independent check.

Check your understanding

An eigenpair has normalized residual $10^{-14}$ and condition estimate $10^{10}$. Is it reasonable to quote fourteen accurate digits of the eigenvalue?

Solution

No. The residual is a backward-error measure; sensitivity can amplify it by roughly the condition estimate, leaving far fewer forward-accurate digits. Verify through higher precision, resolution and endpoint changes, and an independent method. Quote only digits stable under those tests, together with the condition information.

References

• Anderson, E., et al. (1999), LAPACK Users’ Guide, 3rd edition, SIAM. DOI.
• Dias, O. J. C., Godazgar, M., and Santos, J. E. (2022), “Eigenvalue Repulsions and Quasinormal Mode Spectra of Kerr–Newman: An Extended Study,” Journal of High Energy Physics 2022, 076. DOI. Open preprint.
• Jansen, A. (2017), “Overdamped Modes in Schwarzschild–de Sitter and a Mathematica Package for the Numerical Computation of Quasinormal Modes,” European Physical Journal Plus 132, 546. arXiv.
• Leaver, E. W. (1985), “An Analytic Representation for the Quasi-Normal Modes of Kerr Black Holes,” Proceedings of the Royal Society of London A 402, 285–298. DOI. Open PDF.
• Trefethen, L. N., and Embree, M. (2005), Spectra and Pseudospectra: The Behavior of Nonnormal Matrices and Operators, Princeton University Press. DOI.