Logo
ScienceStack
Table of Contents
Fetching ...
Sign In

Computing the matrix exponential with the double exponential formula

Fuminori Tatsuoka, Tomohiro Sogabe, Tomoya Kemmochi, Shao-Liang Zhang

TL;DR

To deal with non-Hermitian matrices, the double exponential (DE) formula specialized to Fourier integrals is considered and another integral representation including an oscillatory term is used and the truncation error is analyzed.

Abstract

This paper considers the computation of the matrix exponential $\mathrm{e}^A$ with numerical quadrature. Although several quadrature-based algorithms have been proposed, they focus on (near) Hermitian matrices. In order to deal with non-Hermitian matrices, we use another integral representation including an oscillatory term and consider applying the double exponential (DE) formula specialized to Fourier integrals. The DE formula transforms the given integral into another integral whose interval is infinite, and therefore it is necessary to truncate the infinite interval. In this paper, to utilize the DE formula, we analyze the truncation error and propose two algorithms. The first one approximates $\mathrm{e}^A$ with the fixed mesh size which is a parameter in the DE formula affecting the accuracy. Second one computes $\mathrm{e}^A$ based on the first one with automatic selection of the mesh size depending on the given error tolerance.

Computing the matrix exponential with the double exponential formula

TL;DR

To deal with non-Hermitian matrices, the double exponential (DE) formula specialized to Fourier integrals is considered and another integral representation including an oscillatory term is used and the truncation error is analyzed.

Abstract

This paper considers the computation of the matrix exponential with numerical quadrature. Although several quadrature-based algorithms have been proposed, they focus on (near) Hermitian matrices. In order to deal with non-Hermitian matrices, we use another integral representation including an oscillatory term and consider applying the double exponential (DE) formula specialized to Fourier integrals. The DE formula transforms the given integral into another integral whose interval is infinite, and therefore it is necessary to truncate the infinite interval. In this paper, to utilize the DE formula, we analyze the truncation error and propose two algorithms. The first one approximates with the fixed mesh size which is a parameter in the DE formula affecting the accuracy. Second one computes based on the first one with automatic selection of the mesh size depending on the given error tolerance.
Paper Structure (12 sections, 1 theorem, 25 equations, 6 figures, 2 algorithms)

This paper contains 12 sections, 1 theorem, 25 equations, 6 figures, 2 algorithms.

Table of Contents

  1. Introduction
  2. The DE formula for Fourier integrals
  3. Computing the matrix exponential with the DE formula
  4. Truncation error
  5. Algorithms
  6. Numerical examples
  7. Computing the scalar exponential with the DE formula
  8. Accuracy dependence on the shift parameter
  9. Accuracy of automatic quadrature
  10. A comparison of the DE formula with the Cauchy integral for non-Hermitian matrices
  11. Conclusion
  12. Derivative of the change of variables

Key Result

Proposition 3.1

Suppose that all the eigenvalues of $A\in\mathbb{C}^{n\times n}$ lie in the left half plane. Let $u(t) = hx_h(t)/\pi - 1$, where $x_h(t)$ is a double exponential transformation for Fourier integrals satisfying $x_h'(t) \le 2\pi/h$. For given $l, r \in \mathbb{Z}$ with $l < r$ and $x_h(lh) \le 1/\sqr where $\|\cdot\|$ is any subordinate norm.

Figures (6)

  • Figure 1: The absolute value of the summands in the upper bound on the truncation error. The left figure illustrates $|x'_h(kh)|$ in \ref{['eq:ub_trunc_err']}, and the right one illustrates $|ku(kh)|$.
  • Figure 2: The error of the double exponential formula for the scalar exponential $\mathrm{e}^z$. The mesh size is set to $h=0.2, 0.1, 0.05$. The left column shows the error for the area $\{z\in\mathbb{C}\colon -30 \le \mathrm{Re}(z) \le 10,\ -20 \le \mathrm{Im}(z) \le 20\}$ and the right one shows the error for $\{z\in\mathbb{C}\colon -5000 \le \mathrm{Re}(z) \le 0,\ -2500 \le \mathrm{Im}(z) \le 2500\}$.
  • Figure 3: Accuracy dependence of the double exponential formula on the shift parameter $\sigma$.
  • Figure 4: The accuracy of Algorithm \ref{['alg:2']}. The left figure shows the accuracy dependence on the input tolerance $\epsilon$ and the error of the computational results. The right figure shows the Error of the case when $A = A_1$ and $\epsilon=10^{-10}$. The horizontal axis is the inverse of the selected mesh size, and the vertical axis is the error.
  • Figure 5: Eigenvalues of the test matrix $M^{-1}K$ in Subsection \ref{['sect:ex_comparison']}.
  • ...and 1 more figures

Theorems & Definitions (2)

  • Proposition 3.1
  • proof