Improved Backward Error Bounds for LU and Cholesky Factorizations

Assuming standard floating-point arithmetic (in base β, precision p) and barring underflow and overflow, classical rounding error analysis of the LU or Cholesky factorization of an n× n matrix A provides backward error bounds of the form |Δ A| ≤ γₙ | L| | U| or |Δ A| ≤ γₙ₊₁ | RT| | R|. Here, L, U, and R denote the computed factors, and γₙ is the usual fraction nu/(1-nu) = nu + O(u²) with u the unit roundoff. Similarly, when solving an n× n triangular system $Tx = b$ by substitution, the computed solution x satisfies (T+Δ T) x = b with |Δ T| ≤ γₙ |T|. All these error bounds contain quadratic terms in u and limit n to satisfy either $nu<1$ or $(n+1)u < 1$. We show in this paper that the constants γₙ and γₙ₊₁ can be replaced by $nu$ and $(n+1)u$, respectively, and that the restrictions on n can be removed. To get these new bounds the main ingredient is a general framework for bounding expressions of the form |ρ-s|, where s is the exact sum of a floating-point number and $n-1$ real numbers and where ρ is a real number approximating the computed sum s. By instantiating this framework with suitable values of ρ, we obtain improved versions of the well-known Lemma 8.4 from [N. J. Higham, Accuracy and Stability of Numerical Algorithms, 2nd ed., SIAM, Philadelphia, 2002] (used for analyzing triangular system solving and LU factorization) and of its Cholesky variant. All our results hold for rounding to nearest with any tie-breaking strategy and whatever the order of summation.