A Schur–Padé Algorithm for Fractional Powers of a Matrix

Abstract

Abstract. A new algorithm is developed for computing arbitrary real powers A p of a matrix A ∈ C n×n . The algorithm starts with a Schur decomposition, takes k square roots of the triangular factor T, evaluates an [m=m] Padé approximant of ð1 − xÞ p at I − T 1=2 k , and squares the result k times. The parameters k and m are chosen to minimize the cost subject to achieving double precision accuracy in the evaluation of the Padé approximant, making use of a result that bounds the error in the matrix Padé approximant by the error in the scalar Padé approximant with argument the norm of the matrix. The Padé approximant is evaluated from the continued fraction representation in bottom-up fashion, which is shown to be numerically stable. In the squaring phase the diagonal and first superdiagonal are computed from explicit formulae for T p=2 j , yielding increased accuracy. Since the basic algorithm is designed for p ∈ ð−1; 1Þ, a criterion for reducing an arbitrary real p to this range is developed, making use of bounds for the condition number of the A p problem. How best to compute A k for a negative integer k is also investigated. In numerical experiments the new algorithm is found to be superior in accuracy and stability to several alternatives, including the use of an eigendecomposition and approaches based on the formula A p ¼ expðp logðAÞÞ.

Blocked Schur Algorithms for Computing the Matrix Square Root

“…Because of the importance of the (quasi-) triangular square root, which arises in algorithms for computing the matrix logarithm [2], [3], matrix pth roots [5], [8], and arbitrary matrix powers [13], this computational kernel is a strong contender for inclusion in any future extensions of the BLAS.…”

Rhapsody in fractional

“…In [29] are presented the Schur, Newton and Inverse Newton methods. The Schur-Newton and Schur-Padé algorithms are also discussed in [30]. Some of these methods impose additional conditions for matrix A.…”

Reducing the influence of tiny normwise relative errors on performance profiles

“…Although originally introduced as a tool for evaluating and comparing optimization software, performance profiles are now widely used, and the second author has found them very useful in the context of algorithms for matrix functions (see, for example, Al-Mohy and Higham [2009;, Higham [2005;2008;2009], Higham and Lin [2011]). …”