The solution to a differential Lyapunov equation can be expressed in closed form as a matrix-valued integral, the so-called finite-horizon Gramian. Such Gramians also have applications in many other areas, such as Gauss-Markov regression. The Gramian is positive semi-definite, and often it is more useful to have a Cholesky factorization of it rather than the Gramian itself. I will present a new efficient numerical method for computing such Cholesky factors of finite-horizon Gramians without first computing the full Gramian. In contrast to other methods applicable to this task, this method is a generalization of the scaling-and-squaring approach for approximating the matrix exponential. It exploits a similar doubling formula for the Gramian, and thereby keeps the required computational effort modest. Most importantly, we have performed a rigorous backward error analysis that guarantees that the approximation is accurate to the round-off error level in double precision if the method parameters are chosen appropriately. I will sketch the outline of this proof. I will also show the results of supporting numerical experiments on a large number of standard test cases, which illustrate that this accuracy is indeed achieved in practice.

[link to pdf] [back to Numdiff-17]