Preconditioning
A preconditioner transforms a linear system into an equivalent one with clustered eigenvalues, so Krylov solvers converge in far fewer iterations.
Reshaping the problem
Krylov methods converge slowly when a matrix is ill-conditioned or has a spread-out spectrum. Preconditioning replaces A x = b with an equivalent system M^{-1} A x = M^{-1} b, where M approximates A but is cheap to invert. If M^{-1} A has clustered eigenvalues, the solver converges in a fraction of the iterations.
The central trade-off
A good preconditioner must satisfy two competing demands: it should approximate A well enough to cluster the spectrum, yet applying M^{-1} each iteration must be cheap. The extreme M = A converges in one step but is as hard as the original problem; M = I does nothing. The art is finding the sweet spot.
Common preconditioners
- Jacobi (diagonal): divide by the diagonal, trivial and cheap, modest effect.
- Incomplete LU/Cholesky (ILU, IC): an approximate factorization that drops small fill, widely effective.
- Multigrid: near optimal for elliptic problems, driving iteration counts nearly independent of grid size.
- Domain decomposition: solve subdomains and couple them, well suited to parallel machines.
Left, right, and split
Preconditioning can be applied on the left, on the right, or split symmetrically. Symmetric positive-definite systems need a symmetric preconditioner to keep CG applicable; split or right preconditioning preserves the residual norm GMRES minimizes. The choice interacts with the Krylov method used.
from scipy.sparse.linalg import spilu, LinearOperator, gmres
import numpy as np
# ILU preconditioner for a sparse matrix A
# ilu = spilu(A.tocsc())
# M = LinearOperator(A.shape, ilu.solve)
# x, info = gmres(A, b, M=M)
Effective preconditioning, often multigrid or incomplete factorizations, is what makes the large sparse solves in breeder Hyperion simulations converge in a practical number of iterations.