Computing a matrix function for exponential integrators

An efficient numerical method is developed for evaluating phi(A), where A is a symmetric matrix and phi is the function defined by phi(x)=(e(x)-1)/x=1+x/2+x(2)/6+.... This matrix function is useful in the so-called exponential integrators for differential equations. In particular, it is related to the exact solution of the ODE system dy/dt=Ay+b, where A and b are t-independent. Our method avoids the eigenvalue decomposition of the matrix A and it requires about 10n(3)/3 operations for a general symmetric n x n matrix. When the matrix is tridiagonal, the required number of operations is only O(n(2)) and it can be further reduced to O(n) if only a column of the matrix function is needed. These efficient schemes for tridiagonal matrices are particularly useful when the Lanczos method is used to calculate the product: of this matrix function (for a large symmetric matrix) with a given vector. (C) 2003 Elsevier B.V. All rights reserved.