Conjugate gradient for ill-posed problems: regularization by preconditioning, preconditioning by regularization

This paper investigates using the conjugate gradient iterative solver for ill-posed problems. We show that preconditioner and Tikhonov-regularization work in conjunction. In particular when they employ the same symmetric positive semi-definite operator, a powerful Ritz analysis allows one to estimate at negligible computational cost the solution for any Tikhonov’s weight. This enhanced linear solver is applied to the boundary data completion problem and as the inner solver for the optical flow estimator.

💡 Research Summary

The paper addresses the challenging class of ill‑posed linear systems that arise frequently in mechanics, inverse problems, and imaging, where the operator A is symmetric positive semi‑definite (SPSD) and its spectrum contains many small eigenvalues that amplify noise. Classical remedies—truncation, spectral shift (Tikhonov regularization), and filtering—are usually applied separately, and iterative Krylov solvers are often avoided because poor conditioning leads to unstable iterates. The authors propose a unified framework that simultaneously uses the same symmetric positive (semi‑)definite matrix M as a preconditioner for the Conjugate Gradient (CG) method and as the regularization operator in a Tikhonov term λM.

The core algorithm is a preconditioned CG (PCG) where the inner product is defined by M (i.e., ⟨u,v⟩_M = uᵀMv). At iteration i the algorithm builds the Krylov subspace K_i(M⁻¹A, M⁻¹r₀) and maintains two orthogonal bases: Z_i, which is M‑orthogonal, and W_i, which is A‑orthogonal. By normalising Z_i one obtains a basis ˆZ_i that is equivalent to the Arnoldi basis for the matrix M⁻¹A. The projected matrix T_i = ˆZ_iᵀAˆZ_i is tridiagonal; its eigen‑decomposition yields Ritz values Θ_i and Ritz vectors V_i = ˆZ_iΞ_i, where Ξ_i contains the eigenvectors of T_i. Crucially, V_i are M‑orthonormal and satisfy V_iᵀAV_i = Θ_i, so they approximate the generalized eigenpairs (A,M).

When the regularized system (A+λM)x = b is solved for a particular λ using the PCG, the stored Ritz data allow a post‑processing step that reconstructs an approximation of the solution for any other λ at negligible cost. Because V_iᵀ(A+λM)V_i = Θ_i + λI, the solution can be expressed as

 x(λ) ≈ x₀ + V_i (Θ_i + λI)⁻¹ V_iᵀ r₀,

where r₀ = b – Ax₀. Thus the whole λ‑parameter space can be explored without re‑solving the linear system, which is a major computational advantage for parameter‑selection techniques such as L‑curve or Morozov’s discrepancy principle.

The authors also show how the M‑norms naturally provide error and solution‑size measures that are more appropriate for ill‑posed problems than the Euclidean norm. The residual norm ‖r_i‖_{M⁻¹} and the solution norm ‖x_i – x₀‖_M give a physically meaningful L‑curve: as iterations progress, the residual decreases while the solution norm grows, tracing a monotone curve toward the upper‑left corner. Moreover, the Ritz values give a direct way to construct a Picard plot in the reduced basis, allowing a comparison of eigenvalue decay versus right‑hand‑side coefficient decay.

Three stopping criteria are proposed: (1) a relative reduction of the M⁻¹‑residual, (2) a balance between residual size and growth of the M‑norm of the iterate (‖r_i‖_{M⁻¹} < ε‖T_i‖_F‖x_i – x₀‖_M), and (3) stagnation detection via the quantity γ_i²/δ_i. The second criterion is especially valuable because it prevents over‑iteration that would otherwise amplify noise in the ill‑conditioned subspace.

The methodology is validated on two problems. The first is a linear boundary‑data completion problem based on the Steklov‑Poincaré formulation. Because the problem is linear, the same M can be used throughout, and the Ritz‑based post‑processing yields accurate reconstructions with fewer CG iterations than standard MINRES or CGLS. The second application is the estimation of optical flow between two images, a nonlinear inverse problem solved by an outer Gauss‑Newton loop. Inside each outer iteration, the linearised normal equations are solved with the proposed PCG‑Ritz scheme. By re‑using the same preconditioner/regularizer M and by post‑processing the Ritz data for a whole family of λ values, the authors avoid costly line‑searches and achieve a noticeable speed‑up while maintaining or improving accuracy.

In summary, the paper makes the following contributions:

1. Demonstrates that a single SPD operator can serve both as a preconditioner and as a Tikhonov regularizer, aligning physical intuition with numerical stabilization.
2. Shows that the PCG algorithm naturally yields Ritz pairs that separate the influence of A and M, enabling inexpensive λ‑independent post‑processing.
3. Provides M‑based error and solution norms that give meaningful L‑curve and Picard diagnostics during iteration.
4. Introduces practical, inexpensive stopping criteria tailored to ill‑posed problems.
5. Extends the approach to nonlinear problems by embedding the linear solver within an outer iteration, achieving acceleration and λ‑parameter flexibility without extra solves.

Overall, the work offers a theoretically sound and practically efficient framework for solving ill‑posed linear and mildly nonlinear problems, turning the traditionally problematic combination of preconditioning and regularization into a synergistic tool.