The matrix exponential is a much-studied matrix function having many applications. The Frechet derivative of the matrix exponential describes the first-order sensitivity of $e^A$ to perturbations in $A$ and its norm determines a condition number for $e^A$. Among the numerous methods for computing $e^A$ the scaling and squaring method is the most widely used. We show that the implementation of the method in [N. J. Higham, The scaling and squaring method for the matrix exponential revisited, SIAM J. Matrix Anal. Appl., 26 (2005), pp. 1179-1193] can be extended to compute both $e^A$ and the Frechet derivative at $A$ in the direction $E$, denoted by $L(A,E)$, at a cost about three times that for computing $e^A$ alone. The algorithm is derived from the scaling and squaring method by differentiating the Pade approximants and the squaring recurrence, reusing quantities computed during the evaluation of the Pade approximant, and intertwining the recurrences in the squaring phase. To guide the choice of algorithmic parameters, an extension of the existing backward error analysis for the scaling and squaring method is developed which shows that, modulo rounding errors, the approximations obtained are $e^{A+\Delta A}$ and $L(A+\Delta A,E+\Delta E)$, with the same $\Delta A$ in both cases, and with computable bounds on $\|\Delta A\|$ and $\|\Delta E\|$. The algorithm for $L(A,E)$ is used to develop an algorithm that computes $e^A$ together with an estimate of its condition number. In addition to results specific to the exponential, we develop some results and techniques for arbitrary functions. We show how a matrix iteration for $f(A)$ yields an iteration for the Frechet derivative and show how to efficiently compute the Frechet derivative of a power series. We also show that a matrix polynomial and its Frechet derivative can be evaluated at a cost at most three times that of computing the polynomial itself and give a general framework for evaluating a matrix function and its Frechet derivative via Pade approximation.
[1]
Larry J. Stockmeyer,et al.
On the Number of Nonscalar Multiplications Necessary to Evaluate Polynomials
,
1973,
SIAM J. Comput..
[2]
R. Ward.
Numerical Computation of the Matrix Exponential with Accuracy Estimate
,
1977
.
[3]
A. Laub,et al.
Condition Estimates for Matrix Functions
,
1989
.
[4]
Charles Kenney,et al.
Polar Decomposition and Matrix Sign Function Condition Estimates
,
1991,
SIAM J. Sci. Comput..
[5]
R. Mathias.
Evaluating the Frechet derivative of the matrix exponential
,
1992
.
[6]
Timothy F. Havel,et al.
Derivatives of the Matrix Exponential and Their Computation
,
1995
.
[7]
Roy Mathias,et al.
A Chain Rule for Matrix Functions and Applications
,
1996,
SIAM J. Matrix Anal. Appl..
[8]
A. Laub,et al.
A Schur-Fréchet Algorithm for Computing the Logarithm and Exponential of a Matrix
,
1998,
SIAM J. Matrix Anal. Appl..
[9]
James Demmel,et al.
LAPACK Users' Guide, Third Edition
,
1999,
Software, Environments and Tools.
[10]
ENCYCLOPEDIA OF MATHEMATICS AND ITS APPLICATIONS
,
1999
.
[11]
Nicholas J. Higham,et al.
A Block Algorithm for Matrix 1-Norm Estimation, with an Application to 1-Norm Pseudospectra
,
2000,
SIAM J. Matrix Anal. Appl..
[12]
Cleve B. Moler,et al.
Nineteen Dubious Ways to Compute the Exponential of a Matrix, Twenty-Five Years Later
,
1978,
SIAM Rev..
[13]
Giles Jewitt,et al.
Computing estimates of continuous time macroeconometric models on the basis of discrete data
,
2005,
Comput. Stat. Data Anal..
[14]
Nicholas J. Higham,et al.
Functions of matrices - theory and computation
,
2008
.
[15]
N. Higham.
The Scaling and Squaring Method for the Matrix Exponential Revisited
,
2005,
SIAM J. Matrix Anal. Appl..