Accuracy and stability of the hyperbolic model time integration scheme revisited

TL;DR

The paper reexamines the hyperbolic-model (HM) time integration scheme for parabolic problems, which augments the PDE with a small second-order time derivative to enable explicit stepping. It derives detailed local-error bounds and shows that the HM local error is fourth-order in the time step $\tau$ for fixed $\varepsilon$, and third-order when $\varepsilon=\mathcal{O}(\tau)$, leading to global convergence of order $2$ or $3$ with respect to the original IVP. A complete stability analysis reveals that for symmetric positive definite $A$, the eigenvalues of the HM amplification matrix $S$ satisfy $|\xi|<1$ if and only if $\varepsilon > \tau^2\lambda_{\max}/4$ (i.e., $\tau<\sqrt{4\varepsilon/\lambda_{\max}}$); however, $\|S^n\|$ can grow significantly when eigenvalues of blocks $S_j$ are nearly equal, potentially compromising convergence despite eigenvalues lying inside the unit disk. Numerical tests corroborate these findings, showing non-monotone behavior and a practical link between amplification growth and the eigenvalue spacing $|\xi_1-\xi_2|$. The work highlights the need for adaptive or more sophisticated parameter choices and diagnostic indicators when employing the HM scheme for stiff parabolic problems.

Abstract

The hyperbolic model (HM) time integration scheme tackles parabolic problems by adding a small artificial second order time derivative term. Described by Samarskii in his 1971 book, the scheme reappeared as the generalized Du Fort-Frankel scheme in a 1976 paper by Gottlieb and Gustafsson. In this note we revisit accuracy and stability properties of the scheme. In particular, we show that the stability condition, formulated by Samarskii based on operator inequalities, coincides with the requirement that the eigenvalues of the amplification matrix (the stability function operator) are smaller than one in absolute value. However, under this condition, the norm of this matrix may exceed one and this, as recently pointed out by Corem and Ditkowski (2012), may corrupt convergence of the scheme. Hence, we also discuss whether this eventual stability lack can be detected and mitigated in practice.

Figures (6)

• Figure 1: Top plots: the first row entries $S_{11}$, $S_{12}$ of the amplification matrix $S$ versus a range of $\mu=\tau\lambda$ for $\tau=3\cdot10^{-5}$, $\varepsilon=2\cdot10^{-4}$ (left plot) and $\varepsilon=1\cdot10^{-4}$ (right plot). Bottom plots: the real parts of the eigenvalues $\xi_{1,2}$ of the amplification matrix versus the same $\mu$ range, for the same $\varepsilon$ and $\tau$ values as at the top plots. The curves for $\xi_1$ and $\xi_2$ coincide if $\xi_{1,2}\in\mathbb{C}$.
• Figure 2: Left plot: the HM error $\tilde{e}(t)$ and its upper estimate \ref{['et_est']} versus time. Right plot: $\varepsilon\tilde{y}"(t)$ versus time. The plots are obtained for $\varepsilon=2\cdot10^{-4}$, $\omega=\lambda_{\max}=1000$.
• Figure 3: Local (left plot) and global (right plot) error values of the HM scheme versus the time step size $\tau$. The plots are obtained for $\varepsilon=2\cdot10^{-4}$, $\lambda_{\max}=1000$, $T=3\cdot10^{-3}$.
• Figure 4: Left plot: norms $\|S^n\|$ for the amplification matrix $S$ versus time. Right plot: $\max_{1\leqslant n\leqslant T/\tau}\|S^n\|$, with $n$ being the time step number, versus the estimate \ref{['diff_est']} for $|\xi_1-\xi_2|^{-1}$. The plots are obtained for $\varepsilon=2\cdot10^{-4}$, $\lambda_{\max}=1000$, $\tau=3\cdot10^{-5}$.
• Figure 5: Left plot: norms $\|S^n\|$ for the amplification matrix $S$ versus time for the 1D heat equation. Right plot: error values of the HM scheme versus the time step size $\tau$ for the 1D heat equation. The plots are obtained for $\varepsilon=2\cdot10^{-4}$, $\lambda_{\max}=1033$, $\tau=3\cdot10^{-5}$, $T=3\cdot10^{-3}$.

Theorems & Definitions (8)

• Proposition 1
• proof
• Proposition 2
• proof
• Proposition 3
• Remark 1
• proof
• Remark 2