Gram-Schmidt Orthogonalization

The process of turning a set of vectors into an orthonormal basis, and the numerical care needed to do it stably.

Building an orthonormal basis

Gram-Schmidt takes a set of linearly independent vectors and produces an orthonormal set spanning the same space. Each vector has its projections onto the already-processed directions subtracted off, leaving a component orthogonal to all of them, which is then normalized. This process is the conceptual heart of QR factorization and of the Krylov orthogonalization inside Arnoldi and GMRES.

Classical vs modified

The classical Gram-Schmidt algorithm computes all projections against the original vector at once. It is mathematically correct but numerically unstable: rounding errors cause the computed vectors to lose orthogonality, sometimes severely. The modified Gram-Schmidt algorithm subtracts each projection sequentially, updating the working vector before computing the next projection. This small reordering greatly improves stability at the same operation count, and it is the version used in practice.

import numpy as np

def modified_gram_schmidt(A):
    A = A.astype(float).copy()
    n = A.shape[1]
    Q = np.zeros_like(A); R = np.zeros((n, n))
    for j in range(n):
        v = A[:, j]
        for i in range(j):
            R[i, j] = Q[:, i] @ v
            v = v - R[i, j] * Q[:, i]
        R[j, j] = np.linalg.norm(v)
        Q[:, j] = v / R[j, j]
    return Q, R

Reorthogonalization

When vectors are nearly dependent, even modified Gram-Schmidt loses orthogonality. The remedy is reorthogonalization: applying the projection step a second time recovers full orthogonality to machine precision. The twice-is-enough rule states that two passes suffice in almost all cases, which is why iterative eigensolvers use selective or full reorthogonalization to combat the loss of orthogonality that plagues Lanczos.

Alternatives

For computing a QR factorization directly, Householder reflections are more stable than Gram-Schmidt and are the default in dense linear algebra libraries. Gram-Schmidt remains preferred when vectors arrive one at a time (as in Krylov methods) and only a few need to be orthogonalized against a growing set, where Householder's all-at-once structure is inconvenient.