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 one has solved a problem with the same physical connectivity as the original problem. Theorems of Oettli, Prager and Skeel 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. We incorporate these results into the stopping criterion of the iterative refinement step of a direct sparse matrix solver and verify by numerical experiments 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. We also discuss a condition estimator corresponding to this new backward error which 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.