Random Projection Technique For Solving Symmetric Indefinite Linear Systems

It is well known that the Bunch-Kaufinan algorithm,bounded Bunch-Kaufinan al-gorithm and Aasen's algorithm are three of the most widely used methods for solving symmetric indefinite linear systems Ax ?b,yet they all are known to suffer from occa-sional numerical instabilities.For example,the Bunch-Kaufman algorithm may produce potentially exponential element growth or unbounded entries in the unit lower triangular matrix L.The flops of the bounded Bunch-Kaufinan algorithm may achieve the worst case,that is,its flops may be the same as that of the Bunch-Parlett algorithm.Moreover,there also exists exponential element growth factors for the bounded Bunch-Kaufman al-gorithm.Researchers have found that there exist pathological matrices for the Bunch-Kaufman algorithm and bounded Bunch-Kaufman algorithm respectively.There are no known pathological matrices for Aasen's algorithm,but the element growth factor for Aasen's algorithm may grow exponentially by analysis and cause numerical instability problems.In this work,we develop a randomized complete pivoting(RCP)algorithm for solv-ing symmetric indefinite linear systems,and an efficient pivoting strategy,simplified Bunch-Kaufinan partial pivoting(SBKP)strategy,for selecting symmetric pivots.The RCP algorithm is based on this strategy and obtains a block LDLT factorization.In this pa-per,we also show a detailed theoretical analysis of the efficiency and stability of the RCP algorithm,and compare it with that of three other non-randomized algorithms(Bunch-Kaufman algorithm,bounded Bunch-Kaufinan algorithm and Aasen's algorithm)to illus-trate its effectiveness and reliability.In terms of theoretical analysis,we first perform probabilistic analysis on the column pivot reliability of the random projection of the RCP algorithm.Then the upper bounds on the entries of the L matrix and the element growth factor are analyzed base on the previous result.The RCP algorithm enjoys theoretical element growth and bounded entries in the factorization comparable to that which results from complete pivoting,up to a theoretical failure probability that exponentially decays with an oversampling parameter.Moreover,in our numerical experiments,the RCP algorithm is better than the Bunch-Kaufinan algo-rithm,and comparable to the bounded Bunch-Kaufman algorithm and Aasen's algorithm in terms of backward error,element growth factor,and the L entry upper bound.We analyze the operation counts of the RCP algorithm to illustrate its efficiency.The RCP algorithm is slightly slower than the Bunch-Kaufinan algorithm,bounded Bunch-Kaufman algorithm and Aasen's algorithm from theoretical analysis and numerical results,but this difference can decrease to a negligible amount as matrix size increases.Further-more,we compare two blocked versions of the computations for implementing the RCP algorithm,and find that picking one pivot every time is more reliable than picking multiple pivots in each loop.In terms of stability,we perform finite precision analysis on the RCP algorithm,which shows that RCP is as numerically stable as Gaussian elimination with complete pivoting.In addition,RCP has been observed to be numerically stable in our extensive numerical experiments.According to the results of the finite precision analysis,updating random projections could potentially result in numerical instability.On the other hand,the input matrix A could itself be nearly rank deficient,leading to potentially low-quality column pivots.Therefore we proposed solving these two problems by performing a modi-fied algorithm.However,we did not experimentally see any need for correcting numerical instability arising from updating random projections.Additionally,the other work in this paper is that we show a detailed analysis on the element growth factor of Aasen's algorithm and obtain its upper bound of 2n-1,which is much smaller than the 4n-2 upper bound proposed by Higham.We also analyze and prove that this upper bound is not tight when the matrix size is greater than or equal to 6.We construct matrix examples of n = 4,5,6,respectively,to validate our analysis results.