The Lanczos Algorithm

A Krylov method that reduces a large symmetric matrix to a small tridiagonal one, yielding extreme eigenvalues cheaply.

Tridiagonalizing on the fly

The Lanczos algorithm builds an orthonormal basis for the Krylov subspace of a symmetric matrix A. Because A is symmetric, the projection of A onto this basis is tridiagonal, and it is generated by a three-term recurrence: each new basis vector depends only on the previous two. This makes each step cheap in both work and memory, needing one matrix-vector product per iteration.

From tridiagonal to eigenvalues

After m steps the algorithm has a small m-by-m symmetric tridiagonal matrix T whose eigenvalues (the Ritz values) approximate the eigenvalues of A. The extreme eigenvalues, largest and smallest, converge first and fastest, often to high accuracy in far fewer than n iterations. This is exactly what is needed for stability analysis, where the most unstable mode corresponds to an extreme eigenvalue.

import numpy as np

def lanczos(A, v0, m):
    n = len(v0)
    V = np.zeros((n, m))
    alpha = np.zeros(m); beta = np.zeros(m)
    V[:, 0] = v0 / np.linalg.norm(v0)
    w = A @ V[:, 0]
    alpha[0] = V[:, 0] @ w
    w = w - alpha[0] * V[:, 0]
    for j in range(1, m):
        beta[j] = np.linalg.norm(w)
        V[:, j] = w / beta[j]
        w = A @ V[:, j]
        alpha[j] = V[:, j] @ w
        w = w - alpha[j] * V[:, j] - beta[j] * V[:, j-1]
    return alpha, beta[1:], V

The loss of orthogonality

In finite-precision arithmetic the basis vectors lose orthogonality as Ritz values converge, producing spurious duplicate eigenvalues. Practical implementations use selective or full reorthogonalization, or postprocessing to detect and remove ghosts. The implicitly restarted Lanczos method (as in ARPACK) restarts the recurrence to keep the subspace small and orthogonality manageable.

Relation to conjugate gradient

Lanczos and conjugate gradient are two views of the same underlying process: CG solves linear systems, Lanczos finds eigenvalues, both driven by the symmetric three-term recurrence. For nonsymmetric matrices the analogous tool is Arnoldi iteration.