Computing the Matrix Exponential with an Optimized Taylor Polynomial Approximation

Abstract

A new way to compute the Taylor polynomial of a matrix exponential is presented which reduces the number of matrix multiplications in comparison with the de-facto standard Paterson-Stockmeyer method for polynomial evaluation. Combined with the scaling and squaring procedure, this reduction is sufficient to make the Taylor method superior in performance to Padé approximants over a range of values of the matrix norms. An efficient adjustment to make the method robust against overscaling is also introduced. Numerical experiments show the superior performance of our method to have a similar accuracy in comparison with state-of-the-art implementations, and thus, it is especially recommended to be used in conjunction with Lie-group and exponential integrators where preservation of geometric properties is at issue.

Computing the matrix sine and cosine simultaneously with a reduced number of products

“…The search of fast algorithms for evaluating matrix polynomials has received considerable interest in the recent literature [2,3,17,19,22,24,28,29]. We next briefly summarize how to approximate the matrix sine and cosine functions by means of certain polynomials involving a reduced number of matrix products.…”

“…We next briefly summarize how to approximate the matrix sine and cosine functions by means of certain polynomials involving a reduced number of matrix products. This reduction essentially follows the same approach used in [9] to minimise the number of commutators appearing in different Liegroup integrators and was successfully adapted to the Taylor expansion of the exponential matrix in [2] and especially in [3].…”

“…The same test matrices have been generated as in Remark 5 of [3] and all matrices are adjusted to have 1-norms over (10 −4 , 10 from the Matrix Function Toolbox [11]. The reference solutions of the matrix cosine and sine have been calculated with Mathematica with 100 digits of precision.…”

“…There are different techniques to compute efficiently the exponential of a matrix [2,3,5,12,21,25,26,27]. However, using any of these general algorithms to approximate the unitary matrix e −itA in (1) involves products of complex matrices making them computationally expensive.…”

Numerical Picard Iteration Methods for Simulation of Non-Lipschitz Stochastic Differential Equations

“…• Symmetries: We consider symmetries in the operators and decompose to symmetrical operators, such that we could apply fast solver methods for each symmetric operator-part, see [7][8][9]. • Deterministic-Stochastic splitting: We decompose into deterministic and stochastic parts, while we apply fast deterministic and fast stochastic solver, see [10,11].…”