DiscoTEX 1.0: Discontinuous collocation and implicit-turned-explicit (IMTEX) integration symplectic, symmetric numerical algorithms with higher order jumps for differential equations I: numerical black hole perturbation theory applications

TL;DR

A generic numerical algorithm which constructs discontinuous spatial and temporal discretisations by operating on discontinuous Lagrange and Hermite interpolation formulae, respectively is described and it is demonstrated that numerical weak-form solutions can be recovered to high-order accuracy by solving a first-order reduced system of ODEs.

Abstract

Dirac $δ-$ distributionally sourced differential equations emerge in many dynamical physical systems from machine learning, finance, neuroscience, and seismology to black hole perturbation theory. These systems lack exact analytical solutions and are thus best tackled numerically. We describe a generic numerical algorithm which constructs discontinuous spatial and temporal discretisations by operating on discontinuous Lagrange and Hermite interpolation formulae, respectively. By solving the distributionally sourced wave equation, possessing analytical solutions, we demonstrate that numerical weak-form solutions can be recovered to high-order accuracy by solving a first-order reduced system of ODEs. The method-of-lines framework is applied to the \texttt{DiscoTEX} algorithm i.e. through \underline{dis}continuous \underline{co}llocation with implicit\underline{-turned-explicit} integration methods which are symmetric and conserve symplectic structure. Furthermore, the main application of the algorithm is proved by calculating the amplitude at any desired location within the numerical grid, including at the position (and at its right and left limit) where the wave- (or wave-like) equation is discontinuous via interpolation using \texttt{DiscoTEX}. This is demonstrated, firstly by solving the wave- (or wave-like) equation and comparing the numerical weak-form solution to the exact solution. We further demonstrate how to reconstruct the gravitational metric perturbations from weak-form numerical solutions of a non-rotating black hole, which do not have known exact analytical solutions, and compare them against state-of-the-art frequency domain results. We conclude by motivating how \texttt{DiscoTEX}, and related numerical algorithms, both open a promising new alternative waveform generation route for modelling highly asymmetric binaries and complement current frequency domain methods.

Figures (18)

• Figure 1: Top: Plot showing the $l_{\infty}$ of the relative error from the exact and numerical solution with the number of Chebyshev nodes $N$ versus the number of jumps $J$. It is readily clear that for optimal implementations of the algorithm, we need to study both these numerical algorithm optimisation control factors. Bottom: From the plot on the top it is clear the optimal number of nodes is $N=32$ Chebyshev nodes with about $J=12$ jumps.
• Figure 2: Plots of the exact spatial discontinuous function as given by equation \ref{['ch3_simplefunction_exact']} against the numerical solution is obtained with traditional Lagrangian smooth interpolation methods and followed with the correction introduced by the discontinuous algorithm. We demonstrate the accuracy of the algorithm by considering both the first- and fourth-order spatial derivatives.
• Figure 3: Fractional energy error from the evolution of equation (25) as given in Ref.da2023hyperboloidal with no RHS terms i.e $\mathcal{S} = S_{lm}(\tau,\sigma) = 0$ and potential on the LHS is $V_{l}= 0$. We observe both HF EX RK and IMTEX Hermite integrators errors are bounded, with lowest order IMTEX integrators significantly outperforming even the highest HF EX RK schemes. Furthermore, it is observed at this juncture, that there is no point in using a Hermite IMTEX scheme of order higher than 4, this is likely due to the build-up of round-off error and the fact most of the residual error comes from the spatial discretisation which becomes more significant as the order of Hermite IMTEX increases.
• Figure 4: Numerical error associated with the evaluation of the integral in equation \ref{['ch3_legendreQP_time']} and with the numerical scheme highlighted in equation \ref{['ch33_hermiteintegration_hermite6rule']} corrected by incorporating the discontinuous behaviour through the discontinuous time-integration rule given in \ref{['sec2_disco_time_h6']}. For an order-6 integration scheme, sixth-order convergence is observed as given by the DH6 line. For a smooth integrator as given by line SH6 inaccurate results are recorded. In a companion paper i.e. discotexII we further compare results up to order-12th.
• Figure 5: Numerical weak-form solution to $\Psi(t,x)$ obtained via DiscoTEX with moving radiation boundary conditions. Left: Numerical field $\Psi(t,x)$ for a point-particle in time-dependent linear motion $r_{p}(t) = v t$, where $v$ is the particle's velocity. Specifically, here $r_{p} \approx 0.0634$ at the coordinate-time $t=0.1902$ and $v=1/4$. Right: Waveform for the point-particle computed on the numerical domain $x \in [-4,4]$ and $t \in [-1.521, 6.602]$.