Opposite-Side Residuals in Working-Precision Newton-Schulz Refinement of LU-Based Matrix Inverses: Measurements and an Empirical Predictor

An inverse X of a matrix A computed through an LU factorization comes with a guarantee on only one of its two residuals, I - AX or I - XA, according to whether it was obtained by solving AX = I or XA = I. In practice the two are usually of the same size, as Wilkinson observed. A Newton-Schulz step, X ← X + X(I - AX), squares both residuals in exact arithmetic. We measure, in exact rational arithmetic, what one such step does when it is carried out in the working precision, on about 2,000 double-precision matrices with sizes 8 to 1024 and target condition numbers up to 1012. The step does not harm the residual it is built from, and often improves it, and it makes the other one larger. For matrices that are ill-conditioned through their singular values the other residual grows with about twice the conditioning exponent, to within a factor of 10 to 25, at the median, of the largest disparity between the two sides that any approximate inverse can have. Under row or column grading there is no such damage on the corresponding side. One empirical predictor describes both behaviours: the off-side residual is close to (u/2n) ‖(|X||A|)2‖F, or its mirror image, where u is the unit roundoff. With constants fixed on three families of matrices and committed in advance, it predicted the off-side residual on three further families to within a factor of about 1.6 in 90% of matrices, over 14 decades. It is looser on classical test matrices. The damage scales with the precision of the residual: for condition numbers of 106 and above, 11 more bits reduce it by about 211, and with an accurate residual one step gave the entrywise correctly rounded inverse over a wide range of conditioning. A heuristic rounding model with no fitted parameter accounts for the size of the effect under left-to-right summation, and its predictions for three other arithmetics, made in advance, were met: pairwise summation, a fused multiply-add, and the blocked accumulation of an optimised BLAS, where the constant is smaller. The thresholds of every experiment were committed before its code was written, and the predictions that failed are reported. The mechanism is elementary and the qualitative fact is known. We did not find the size of the effect stated in the references we could consult; one that may contain it could not be obtained, and another was read only in part.