Computing the matrix exponential and the Cholesky factor of a related finite horizon Gramian

Abstract In this article, an efficient numerical method for computing both the matrix exponential and a finite horizon controllability Gramian in Cholesky-factored form is proposed. The method is applicable to general dense matrices of moderate size and produces a Cholesky factor of the Gramian without computing the full product. It is a generalization of the scaling-and-squaring approach for approximating the matrix exponential and exploits a similar doubling formula for the Gramian to compute both at once. The required computational effort is thereby kept modest. Most importantly, a rigorous backward error analysis is provided, which guarantees that the approximation is accurate to the round-off error level in double precision. This accuracy is illustrated in practice on a large number of standard test examples. The method has been implemented in the Julia package FiniteHorizonGramians.jl, which is available online under the MIT license. Code for reproducing the experimental results is included in this package, as well as code for determining the optimal method parameters. The analysis can thus easily be adapted to a different finite-precision arithmetic.