Accurate symmetric indefinite linear equation solvers

Abstract

The Bunch-Kaufman factorization is widely accepted as the algorithm of choice for the direct solution of symmetric indefinite linear equations; it is the algorithm employed in both UNPACK and LAPACK. It has also been adapted to sparse symmetric indefinite linear systems. While the Bunch-Kaufman factorization is normwise backward stable, its factors can have unusual scaling, with entries bounded by terms depending both on ∥A∥ and on K(A). This scaling, combined with the block nature of the algorithm, may degrade the accuracy of computed solutions unnecessarily. Overlooking the lack of a triangular factor bound leads to a further complication in LAPACK such that the LAPACK Bunch-Kaufman factorization can be unstable. We present two alternative algorithms, close cousins of the Bunch-Kaufman factorization, for solving dense symmetric indefinite systems. Both share the positive attributes of the Bunch-Kaufman algorithm but provide better accuracy by bounding the triangular factors. The price of higher accuracy can be kept low by choosing between our two algorithms. One is appropriate as the replacement for the blocked LAPACK Bunch-Kaufman factorization; the other would replace the LINPACK-like unblocked factorization in LAPACK. Solving sparse symmetric indefinite systems is more problematic. We conclude that the Bunch-Kaufman algorithm cannot be rescued effectively in the sparse case. Imposing the constraint of bounding the triangular factors leads naturally to one particular version of the Duff-Reid algorithm, which we show gives better accuracy than Liu's sparse variant of the Bunch-Kaufman algorithm. We extend the work of Duff and Reid in two respects that often provide higher efficiency: a more effective procedure for finding pivot blocks and a stable extension to pivot blocks of size larger than two.