Solving sparse linear systems with sparse backward error

When solving sparse linear systems, it is desirable to produce the solution of a nearby sparse problem with the same sparsity structure. This kind of backward stability helps guarantee, for example, that a problem with the same physical connectivity as the original has been solved. Theorems of Oettli, Prager [Number Math., 6 (1964), pp. 405-409] and Skeel [Math. Comput., 35 (1980), pp. 817-832] show that one step of iterative refinement, even with single precision accumulation of residuals, guarantees such a small backward error if the final matrix is not too ill-conditioned and the solution components do not vary too much in magnitude. These results are incorporated into the stopping criterion of the iterative refinement step of a direct sparse matrix solver, and numerical experiments verify that the algorithm frequently stops after one step of iterative refinement with a componentwise relative backward error at the level of the machine precision. Furthermore, calculating this stopping criterion is very inexpensive. A condition estimator corresponding to this new backward error is discussed that provides an error estimate for the computed solution. This error estimate is generally tighter than estimates provided by standard condition estimators. We also consider the effects of using a drop tolerance during the LU decomposition.