Solving Systems of Linear Equations

Abstract

Solving a system of linear equations \( \mathbf {Ax} = \mathbf {b}\) for \(\mathbf {x}\) plays a key role in enough problems and algorithms to deserve a chapter of its own. We assume here that there are as many scalar equations as there are scalar unknowns (so \(\mathbf {A}\) is square), and that the solution is unique (so \(\mathbf {A}\) is invertible). Mathematically, this solution is then given in closed form as \(\mathbf {x} = \mathbf {A}^{- 1}\mathbf {b}\), but this is not how it should be computed numerically. Various methods best avoided for solving a system of linear equations numerically are pointed out, and numerically robust methods are described. The condition number of \(\mathbf {A}\) is a good indicator of the intrinsic difficulty of the problem, independently of the method used for solving it. Direct methods attempt to evaluate \(\mathbf {x}\) from the numerical values of \(\mathbf {A}\) and \(\mathbf {b}\) by a finite number of steps, whereas classical iterative methods aim at converging towards the solution in an infinite number of steps. Krylov subspace iteration has blurred the lines, as it would converge to the exact solution in a finite number of steps if the computations were carried out exactly, just as direct methods.