  The Bunch-Kaufman algorithm and Aasen's algorithm are two of the most widely
used methods for solving symmetric indefinite linear systems, yet they both are
known to suffer from occasional numerical instability due to potentially
exponential element growth or unbounded entries in the matrix factorization. In
this work, we develop a randomized complete pivoting (RCP) algorithm for
solving symmetric indefinite linear systems. RCP is comparable to the
Bunch-Kaufman algorithm and Aasen's algorithm in computational efficiency, yet
enjoys theoretical element growth and bounded entries in the factorization
comparable to that of complete-pivoting, up to a theoretical failure
probability that exponentially decays with an oversampling parameter. Our
finite precision analysis shows that RCP is as numerically stable as Gaussian
elimination with complete pivoting, and RCP has been observed to be numerically
stable in our extensive numerical experiments.
