Numerical Linear Algebra

Numerical linear algebra is the branch of linear algebra concerned with computing solutions to problems involving vectors and matrices using finite-precision arithmetic. Its principal tasks include solving linear systems, fitting least-squares models, computing eigenvalues and eigenvectors, and constructing matrix factorizations. Unlike purely symbolic treatments, it studies how rounding errors, sensitivity to input data, matrix structure, and computational resources affect the reliability and efficiency of an algorithm. Matrix decompositions provide a central framework for these computations. (netlib.org)

Conditioning and numerical stability

Computers commonly represent numerical data through floating-point arithmetic, so intermediate results generally incur rounding errors. Two distinct ideas govern their interpretation: conditioning describes the sensitivity of a mathematical problem to perturbations in its inputs, whereas stability describes how an algorithm propagates computational errors. A well-designed algorithm cannot eliminate sensitivity intrinsic to the problem. (netlib.org)

For an invertible matrix $A$ , its condition number in a compatible induced norm is

\kappa(A)=\|A\|\,\|A^{-1}\|.

A large condition number indicates that some small relative perturbations can cause much larger relative changes in the solution. Forward error measures the difference between a computed answer and the exact answer; backward error measures the perturbation of the input needed to make the computed answer exact. A backward-stable method produces the exact solution of a nearby problem, but its forward error can still be substantial when the original problem is ill-conditioned. (netlib.org)

For a linear system $Ax=b$ , the residual of an approximation $\hat{x}$ is $r=b-A\hat{x}$ . A small residual does not by itself establish that $\hat{x}$ is close to the true solution: conditioning also matters. Solver diagnostics therefore often combine residual-based error estimates with condition estimation. (netlib.org)

Direct methods and matrix factorizations

Direct methods reduce a problem through a finite sequence of operations, rather than primarily improving an approximation until a stopping criterion is met. Matrix factorization transforms a matrix into factors whose structure makes subsequent computations simpler. (cs.cornell.edu)

Gaussian elimination leads to LU decomposition, commonly written

PA=LU,

where $P$ represents row permutations, $L$ is lower triangular, and $U$ is upper triangular. Forward and backward substitution then solve the system. Partial pivoting selects a largest-magnitude available entry in the current column, helping avoid small pivots and excessive element growth; it is not an unconditional guarantee against every unfavorable case. (netlib.org)

For a real symmetric positive-definite matrix, Cholesky decomposition gives $A=LL^T$ without requiring pivoting. Another fundamental decomposition is QR decomposition, $A=QR$ , where $Q$ has orthonormal columns and $R$ is upper triangular. Multiplication by an orthogonal matrix preserves Euclidean norms, a property important for stable transformations. (netlib.org)

The computational cost of conventional dense LU factorization grows cubically with matrix dimension. Once the factors are available, they can be reused to solve systems with additional right-hand sides. Such factor-and-solve procedures avoid explicitly constructing the inverse matrix merely to solve $Ax=b$ . (johnfoster.pge.utexas.edu)

Least squares and singular values

Linear least squares seeks

\min_x\|Ax-b\|_2.

It is central to ordinary least squares and linear regression. For an overdetermined system with full column rank, the minimizing solution is unique. QR factorization converts the problem into a triangular solve while preserving the relevant norm. (netlib.org)

The normal equations $A^TAx=A^Tb$ are mathematically equivalent under these assumptions, but

\kappa_2(A^TA)=\kappa_2(A)^2.

Consequently, forming and solving the normal equations can amplify numerical difficulties compared with a QR-based approach. (courses.seas.harvard.edu)

Singular value decomposition writes $A=U\Sigma V^*$ , where the diagonal entries of $\Sigma$ are nonnegative singular values and $^{*}$ denotes conjugate transpose. It supports minimum-norm solutions and the treatment of rank-deficient problems. In computation, deciding the effective rank requires a threshold because very small singular values may reflect either genuine structure or perturbations. Tikhonov regularization modifies the fitting problem by adding a penalty, reducing sensitivity at the expense of changing the solution being sought. (netlib.org)

Eigenvalue computations

Computing eigenvalues and eigenvectors is another major task. For general dense matrices, the QR algorithm uses successive transformations to obtain a Schur form, from which eigenvalues can be extracted. Real Schur form is quasi-upper triangular, allowing two-by-two blocks for complex-conjugate eigenvalue pairs. (netlib.org)

Practical QR implementations use orthogonal or unitary transformations and can be backward stable. Nevertheless, eigenvalue sensitivity depends on the matrix: backward stability alone does not ensure small errors in each individual eigenvalue or eigenvector. (netlib.org)

Sparse systems and iterative methods

For large sparse matrices, algorithms exploit the small proportion of nonzero entries. Iterative solvers often rely chiefly on matrix–vector products, avoiding a complete dense factorization. Krylov subspace methods construct approximations from spaces generated by repeated applications of the matrix to an initial residual. (netlib.org)

The conjugate gradient method applies to symmetric positive-definite systems. GMRES addresses general nonsymmetric systems by minimizing the residual over the current Krylov space. Its unrestarted form requires growing storage and orthogonalization work; restarting limits those costs but can affect convergence. Preconditioning transforms the system using an auxiliary operator intended to improve convergence, with effectiveness balanced against construction and application costs. (netlib.org)

Software and implementation

BLAS supplies standardized vector, matrix–vector, and matrix–matrix operations; LAPACK builds higher-level solvers and decompositions on these kernels. Blocked algorithms organize work around matrix–matrix operations to improve data reuse and support parallel computing. Performance therefore depends not only on arithmetic counts, but also on memory access and the efficiency of the underlying computational kernels. (netlib.org)