Mathematical analysis

Machine epsilon

Convergence (economics)

Iterative refinement

Eigenvalues and eigenvectors

Economics

Linear system

Physics

Generalized minimal residual method

Double-precision floating-point format

Economic growth

Iterative method

Pure mathematics

Mathematics

Condition number

Mathematical optimization

Political science

Computation

Quantum mechanics

Operating system

Rounding

Norm (philosophy)

Computer science

Algorithm

Single-precision floating-point format

Invertible matrix

Factorization

Solver

Round-off error

Applied mathematics

SIAM Journal on Scientific Computing

We propose a general algorithm for solving an $n\times n$ nonsingular linear system $Ax = b$ based on iterative refinement with three precisions. The working precision is combined with possibly different precisions for solving for the correction term and for computing the residuals. Via rounding error analysis of the algorithm we derive sufficient conditions for convergence and bounds for the attainable forward error and normwise and componentwise backward errors. Our results generalize and unify many existing rounding error analyses for iterative refinement. With single precision as the working precision, we show that by using LU factorization in IEEE half precision as the solver and calculating the residuals in double precision it is possible to solve $Ax = b$ to full single precision accuracy for $\infty$-norm condition numbers $\kappa_{\infty}(A) \le 10^4$, with the $O(n^3)$ part of the computations carried out entirely in half precision. We show further that by solving the correction equations by GMRES preconditioned by the LU factors the restriction on the condition number can be weakened to $\kappa_{\infty}(A) \le 10^8$, although in general there is no guarantee that GMRES will converge quickly. Taking for comparison a standard $Ax = b$ solver that uses LU factorization in single precision, these results suggest that on architectures for which half precision is efficiently implemented it will be possible to solve certain linear systems $Ax = b$ up to twice as fast and to greater accuracy. Analogous results are given with double precision as the working precision.

Accelerating the Solution of Linear Systems by Iterative Refinement in Three Precisions