Numerical linear algebra solves linear systems and related matrix problems efficiently while managing floating-point errors. Gaussian elimination is still widely used, usually as LU decomposition with pivoting.

Common factorizations

  • LU: Splits a matrix into lower and upper triangular factors, with row swaps for pivoting. For a dense square matrix, factorization costs O(n³), then each new right-hand side can be solved in O(n²).
  • QR: Splits a matrix into orthogonal and upper triangular factors. A numerically stable choice for least-squares problems, though not free of rounding errors.
  • SVD: Splits a matrix into singular vectors and singular values. Useful for rank estimation, PCA, and low-rank approximation, including image compression.

Direct or iterative?

  • Direct methods, such as LU and QR, use a finite sequence of operations. Floating-point results are still approximate.
  • Iterative methods, such as Conjugate Gradient and GMRES, refine an initial guess until a stopping criterion is met. Often useful for large sparse systems, but convergence and memory needs depend on the method and problem. Conjugate Gradient requires a symmetric positive-definite matrix.

In practice

Use established libraries rather than implementing solvers from scratch. BLAS provides basic vector and matrix operations; LAPACK builds on them to provide factorizations and solvers used by tools such as NumPy and SciPy.

Faster matrix multiplication algorithms also exist: Strassen’s method takes about O(n²·⁸¹) operations. A better asymptotic complexity does not always mean better practical performance.

Reference