A New Scaling and Squaring Algorithm for the Matrix Exponential

摘要: The scaling and squaring method for the matrix exponential is based on approximation $e^A \approx (r_m(2^{-s}A))^{2^s}$, where $r_m(x)$ $[m/m]$ Pade approximant to $e^x$ integers $m$ $s$ are be chosen. Several authors have identified a weakness of existing algorithms termed overscaling, in which value much larger than necessary chosen, causing loss accuracy floating point arithmetic. Building algorithm Higham [SIAM J. Matrix Anal. Appl., 26 (2005), pp. 1179-1193], used by MATLAB function expm, we derive new that alleviates overscaling problem. Two key ideas employed. first, specific triangular matrices, compute diagonal elements phase as exponentials instead from powers $r_m$. second idea base backward error analysis underlies members sequence $\{\|A^k\|^{1/k}\}$ $\|A\|$, since nonnormal matrices it possible $\|A^k\|^{1/k}$ smaller indeed this likely when occurs algorithms. terms estimated without computing $A$ using 1-norm estimator conjunction with bound form $\|A^k\|^{1/k} \le \max\bigl( \|A^p\|^{1/p}, \|A^q\|^{1/q} \bigr)$ holds certain fixed $p$ $q$ less $k$. improvements truncation bounds balanced potential large $\|A\|$ cause inaccurate evaluation $r_m$ We employ rigorous along some heuristics ensure rounding errors kept under control. Our numerical experiments show generally provides at least good no higher cost, while or usually yields significant accuracy, both.