Logo
ScienceStack
Table of Contents
Fetching ...
Sign In

Stable Hermite transforms via the Golub-Welsch algorithm

Marcus Webb, Georg Maierhofer

Abstract

We introduce an efficient stable algorithm for transforms associated with expansions in Hermite functions interpolated at Hermite polynomial roots. The Hermite transform matrix can be factorised into a diagonal component and an orthogonal matrix, leading to a form which allows both the forward and inverse Hermite transforms to be computed stably. Our novel algorithm computes this factorisation based on the eigendecomposition of the Jacobi matrix associated with Hermite functions. Through numerical experiments, we demonstrate the stability and efficiency gains of this novel method over prior work. Numerical experiments show that the new approach matches or improves on the accuracy of existing stabilized methods, is substantially faster in practice, and enables reliable use of large Hermite expansions in downstream PDE computations. We also provide an open-source implementation, together with reference implementations of previous methods, to facilitate adoption by the community.

Stable Hermite transforms via the Golub-Welsch algorithm

Abstract

We introduce an efficient stable algorithm for transforms associated with expansions in Hermite functions interpolated at Hermite polynomial roots. The Hermite transform matrix can be factorised into a diagonal component and an orthogonal matrix, leading to a form which allows both the forward and inverse Hermite transforms to be computed stably. Our novel algorithm computes this factorisation based on the eigendecomposition of the Jacobi matrix associated with Hermite functions. Through numerical experiments, we demonstrate the stability and efficiency gains of this novel method over prior work. Numerical experiments show that the new approach matches or improves on the accuracy of existing stabilized methods, is substantially faster in practice, and enables reliable use of large Hermite expansions in downstream PDE computations. We also provide an open-source implementation, together with reference implementations of previous methods, to facilitate adoption by the community.

Paper Structure

This paper contains 14 sections, 2 theorems, 23 equations, 6 figures, 2 algorithms.

Table of Contents

  1. Key words:
  2. Mathematics Subject Classification (MSC2020):
  3. Introduction
  4. Hermite transforms
  5. Contributions and paper structure
  6. Algorithms for the Hermite transform
  7. Bunck's stabilized recurrence
  8. Orthogonal Decomposition of the transform matrix
  9. New approach via the Golub--Welsch algorithm
  10. Numerical examples
  11. Performance metrics of the stable Hermite transforms
  12. Downstream task: the Gross--Pitaevskii equation
  13. Conclusions
  14. MATLAB Code

Key Result

Lemma 2.1

The Hermite transform matrix $T \in \mathbb{R}^{N\times N}$, with the nodes $\mathbf{x}$ taken to be the Gauss--Hermite quadrature nodes, satisfies where $Q \in \mathbb{R}^{N\times N}$ is orthogonal and $D \in \mathbb{R}^{N\times N}$ is diagonal, with $\blacktriangleleft$$\blacktriangleleft$

Figures (6)

  • Figure 1: Depiction of the matrix $Q$ and the diagonal of $D$ for $N = 100$.
  • Figure 2: Error in Algorithm \ref{['alg:GolubWelschHermiteTransform']} using Clenshaw's algorithm vs the asymptotic expression described above (cf. townsend2016fast). We have excluded errors that are larger than $10^{-1}$. We speculate that the slowly growing error for large $N$ (in the blue lines) is due to the error in the computation of the Airy function and its derivative.
  • Figure 3: Assembly time for the matrix $T$.
  • Figure 4: Error in the transform matrices ($2$-norm). We have excluded errors that are larger than $10^{-1}$.
  • Figure 5: Performance of the Hermite spatial discretisation on the Gross--Pitaevskii equation \ref{['eqn:Gross-Pitaevskii_equation']}.
  • ...and 1 more figures

Theorems & Definitions (4)

  • Lemma 2.1
  • proof
  • Lemma 2.2
  • proof