Krylov Subspace Methods

The core mechanism of Krylov methods is the iterative construction of the Krylov subspace K_k(A, v) = span{v, Av, A²v, ..., A^(k-1)v}. This subspace is generated by repeatedly applying the matrix A to an initial vector v. The power of the method lies in the fact that for many practical problems, a good approximation to the desired solution (e.g., an eigenvector or the solution to Ax = b) lies within a relatively low-dimensional Krylov subspace, even when A is enormous and sparse.

A defining feature is their matrix-free nature. The algorithm only requires the ability to compute the matrix-vector product A * v, not explicit access to or storage of the full matrix A. This is critical for:

• Extremely large systems where storing A is impossible.
• Implicit operators where A represents a function or a differential operator.
• Sparse matrices where the product can be computed efficiently without forming the dense matrix. This property makes them the algorithm of choice for problems arising from finite element methods, computational fluid dynamics, and quantum chemistry.

To ensure numerical stability, the raw power basis {v, Av, A²v, ...} is transformed into an orthonormal basis for the Krylov subspace. This is achieved via processes like:

• Arnoldi iteration (for general matrices).
• Lanczos algorithm (for Hermitian/symmetric matrices). These processes produce a Hessenberg matrix (Arnoldi) or a tridiagonal matrix (Lanczos) that represents the projection of A onto the Krylov subspace. The eigenvalues of this much smaller projected matrix (Ritz values) approximate the extreme eigenvalues of the original large matrix A.

Krylov methods are foundational for solving large, sparse linear systems Ax = b. Two seminal algorithms are:

• Conjugate Gradient (CG): Used for symmetric positive definite (SPD) matrices. It finds the solution by minimizing the error in the A-norm over the Krylov subspace. It has optimality properties and finite termination in exact arithmetic.
• Generalized Minimal Residual (GMRES): Used for general non-symmetric matrices. It finds the solution that minimizes the residual norm ||b - Ax||₂ over the Krylov subspace. Both methods implicitly build and utilize the Krylov subspace K_k(A, r₀), where r₀ = b - Ax₀ is the initial residual.

The primary cost of Krylov methods grows with the subspace dimension k:

• Storage: O(n*k) for the orthonormal basis vectors.
• Computation: O(k²*n) for full orthogonalization in Arnoldi/GMRES. To manage this cost for problems requiring many iterations, restarting is essential. After m iterations, the algorithm discards the subspace, uses the current best approximation as a new starting vector, and begins building a fresh Krylov subspace. While this can slow convergence, it makes solving very large problems feasible within memory constraints.

The convergence rate of Krylov methods depends heavily on the spectral properties of the matrix A (e.g., eigenvalue distribution). Preconditioning is the critical technique to accelerate convergence by transforming the original system Ax = b into an equivalent one with better spectral properties. A preconditioner M ≈ A⁻¹ is applied, leading to systems like M⁻¹Ax = M⁻¹b. Effective preconditioning (e.g., incomplete LU factorization, multigrid, domain decomposition) is often the difference between a solution that is obtained in seconds versus one that fails to converge.

The Lanczos algorithm is a specialized Krylov subspace method for approximating the extreme eigenvalues and eigenvectors of a large, sparse, symmetric matrix. It constructs an orthonormal basis for the Krylov subspace and reduces the original matrix to a symmetric tridiagonal form, which is much easier to analyze.

• Core Mechanism: Iteratively builds a basis using the recurrence A * q_{k-1} = β_{k-1} * q_{k-2} + α_k * q_{k-1} + β_k * q_k.
• Primary Use: Fundamental for eigenvalue problems in Principal Component Analysis (PCA) and low-rank approximations via Truncated SVD.
• Connection to Compression: Used to efficiently compute the dominant singular values/vectors of large weight matrices for Low-Rank Factorization.

Arnoldi iteration generalizes the Lanczos algorithm to non-symmetric (or non-Hermitian) matrices. It builds an orthonormal basis for the Krylov subspace and reduces the matrix to an upper Hessenberg form.

• Core Mechanism: Employs the modified Gram-Schmidt process to orthogonalize each new vector against all previous basis vectors.
• Primary Use: Forms the foundation for the GMRES (Generalized Minimal Residual) method for solving non-symmetric linear systems Ax = b.
• Application Context: While less common in symmetric factorization for compression, it is crucial for solving the large linear systems that can arise during model training or fine-tuning phases.

The Rayleigh quotient is a scalar function R(A, x) = (x* A x) / (x* x) for a given vector x and a Hermitian matrix A. It provides an estimate for an eigenvalue of A.

• Key Property: For an eigenvector, the Rayleigh quotient equals the corresponding eigenvalue. For an approximation, it provides an optimal eigenvalue estimate given the vector.
• Role in Krylov Methods: In the Lanczos algorithm, the Rayleigh quotients of the basis vectors with respect to the tridiagonal matrix yield Ritz values, which are approximations to the eigenvalues of A. This is central to monitoring convergence.
• Optimization Link: The minimization/maximization of the Rayleigh quotient is directly related to finding extreme eigenvectors, connecting to variational principles in PCA.

The Conjugate Gradient (CG) method is a Krylov subspace algorithm for solving systems of linear equations Ax = b where A is symmetric and positive-definite. It is optimal in that it finds the exact solution in at most n iterations (in exact arithmetic).

• Krylov Subspace: The k-th iterate x_k lies in the Krylov subspace K_k(A, r_0), where r_0 is the initial residual.
• Optimization Perspective: CG can be viewed as an iterative method for minimizing the quadratic form (1/2)xᵀAx - bᵀx.
• System Solution in ML: While direct compression uses eigenvalue methods, efficient solvers like CG are vital for the linear subproblems in Alternating Least Squares (ALS) for tensor factorization or certain training procedures.

Power iteration is the simplest Krylov subspace method, used to find the dominant eigenvector (associated with the eigenvalue of largest magnitude) of a matrix. It repeatedly applies the matrix to a vector and normalizes: x_{k+1} = (A * x_k) / ||A * x_k||.

• Krylov Subspace: The sequence {x_k} spans the Krylov subspace K_k(A, x_0).
• Limitations: Only finds one eigenpair and convergence can be slow if the eigenvalue gap is small.
• Foundation Role: It is the conceptual precursor to more sophisticated methods like Lanczos and Arnoldi, which leverage the entire subspace for better approximations. Used in simple Truncated SVD implementations.

GMRES is a Krylov subspace method for solving non-symmetric linear systems Ax = b. It finds the vector in the Krylov subspace that minimizes the norm of the residual ||b - Ax||.

• Subspace Construction: Uses Arnoldi iteration to build an orthonormal basis for the Krylov subspace K_m(A, r_0).
• Minimization: Solves a least squares problem in the smaller Hessenberg matrix generated by Arnoldi to find the optimal approximation.
• Engineering Relevance: In machine learning pipelines, large, sparse, non-symmetric systems can arise from constraints or specific model architectures. GMRES provides a robust, iterative solver for these scenarios, ensuring numerical stability where direct methods fail.