A New Analysis of Iterative Refinement and Its Application to Accurate Solution of Ill-Conditioned Sparse Linear Systems

Abstract

Abstract. Iterative refinement is a long-standing technique for improving the accuracy of a computed solution to a nonsingular linear system Ax = b obtained via LU factorization. It makes use of residuals computed in extra precision, typically at twice the working precision, and existing results guarantee convergence if the matrix A has condition number safely less than the reciprocal of the unit roundoff, u. We identify a mechanism that allows iterative refinement to produce solutions with normwise relative error of order u to systems with condition numbers of order u −1 or larger, provided that the update equation is solved with a relative error sufficiently less than 1. A new rounding error analysis is given, and its implications are analyzed. Building on the analysis, we develop a GMRES (generalized minimal residual)-based iterative refinement method (GMRES-IR) that makes use of the computed LU factors as preconditioners. GMRES-IR exploits the fact that even if A is extremely ill conditioned the LU factors contain enough information that preconditioning can greatly reduce the condition number of A. Our rounding error analysis and numerical experiments show that GMRES-IR can succeed where standard refinement fails, and that it can provide accurate solutions to systems with condition numbers of order u −1 and greater. Indeed, in our experiments with such matrices-both random and from the University of Florida Sparse Matrix Collection-GMRES-IR yields a normwise relative error of order u in at most three steps in every case.

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

“…In doing so, we \bullet provide rigorous componentwise forward error bounds and both componentwise and normwise backward error bounds; \bullet make minimal assumptions about the solvers used in steps 1 and 4, so that the analysis is applicable to all the situations mentioned above, as well as to others that can be envisaged; \bullet treat the precisions u f , u s , u, and u r as independent parameters. Our results generalize and unify most existing analyses, including the recent forward error analysis of [9]. We make one omission: we do not try to prove a``one step of iterative refinement in fixed precision implies componentwise backward stability"" result [17], [34], which is of lesser practical importance.…”

“…All that matters is that the factors yield an \widetil A with condition number much smaller than that of A. The convergence condition \phi i \ll 1 from the forward error analysis, where \phi i is defined in (3.9), therefore holds if As mentioned above, and explained in detail in [9], \mu i is much less than 1 in the early iterations, so this condition is in practice dominated by the second term, for which we need f (n)u(\gamma f n ) 2 \kappa \infty (A) 2 \ll 1, and hence certainly \kappa \infty (A) < u - 1/2 u - 1 f ; so (8.4) can hold for \kappa \infty (A) greater than u - 1 f . Then the limiting accuracy is, from (3.10), 4pu r cond(A, x) + u.…”

“…Carson and Higham [9] give a new forward error analysis of iterative refinement that identifies a mechanism that allows accurate solutions to be obtained to systems for which A has condition number as large as the order of the reciprocal of the unit roundoff. Their analysis requires the update equation in line 4 of Algorithm 1.1 to be solved with relative error less than Table 1.1 Summary of existing rounding error analyses for iterative refinement in floating point arithmetic indicating (a) whether the analyses apply to LU factorization only or to an arbitrary solver, (b) whether the backward or forward error analyses are componentwise (``comp"") or normwise (``norm""), and (c) the assumptions on the precisions u f , us, u, ur in Algorithm 1.1 (u f = u and us = u f unless otherwise stated).…”

“…comp comp u \geq ur Tisseur [37] 2001 arb. norm norm u \geq ur Langou et al [24] 2006 LU norm norm u f \geq u = ur Carson and Higham [9] 2017 arb. comp --u \geq ur This work 2017 arb.…”

“…The fourth column shows a bound on \kappa \infty (A) that must hold for the analysis to guarantee convergence (under suitable conditions described in the text) with limiting backward or forward errors of the orders shown in the final three columns. [9] introduce a new form of iterative refinement that corresponds to Algorithm 1.1 with u f = u and u r = u 2 and a special way of solving the correction equation on line 4. The algorithm is intended to handle the situation where A is extremely ill conditioned, possibly even singular to working precision, so that \kappa \infty (A) could exceed u - 1 .…”