[Paper Review] Toward exponentially-convergent simulations of extreme-mass-ratio inspirals: A time-domain solver for the scalar Teukolsky equation with singular source terms

Manas Vishal, Scott E. Field+3 more|arXiv (Cornell University)|Jul 3, 2023

Pulsars and Gravitational Waves ResearchPhysics and Astronomy3 citations

TL;DR

This paper presents a multi-domain discontinuous Galerkin solver for the scalar Teukolsky equation in Kerr spacetime, using spherical harmonic decomposition and a first-order symmetric hyperbolic formulation to achieve global spectral accuracy despite singular Dirac delta source terms. The method enables exponentially convergent simulations of extreme-mass-ratio inspirals with accurate waveforms at future null infinity and correct treatment of the point-particle singularity.

ABSTRACT

Gravitational wave signals from extreme mass ratio inspirals are a key target for space-based gravitational wave detectors. These systems are typically modeled as a distributionally-forced Teukolsky equation, where the smaller black hole is treated as a Dirac delta distribution. Time-domain solvers often use regularization approaches that approximate the Dirac distribution that often introduce small length scales and are a source of systematic error, especially near the smaller black hole. We describe a multi-domain discontinuous Galerkin method for solving the distributionally-forced Teukolsky equation that describes scalar fields evolving on a Kerr spacetime. To handle the Dirac delta, we expand the solution in spherical harmonics and recast the sourced Teukolsky equation as a first-order, one-dimensional symmetric hyperbolic system. This allows us to derive the method's numerical flux to correctly account for the Dirac delta. As a result, our method achieves global spectral accuracy even at the source's location. To connect the near field to future null infinity, we use the hyperboloidal layer method, allowing us to supply outer boundary conditions and providing direct access to the far-field waveform. We document several numerical experiments where we test our method, including convergence tests against exact solutions, energy luminosities for circular orbits, the scheme's superconvergence properties at future null infinity, and the late-time tail behavior of the scalar field. We also compare two systems that arise from different choices of the first-order reduction variables, finding that certain choices are numerically problematic in practice. The methods developed here may be beneficial when computing gravitational self-force effects, where the regularization procedure has been developed for the spherical harmonic modes and high accuracy is needed at the Dirac delta's location.

Motivation & Objective

• To develop a time-domain solver for the scalar Teukolsky equation with singular source terms arising in extreme-mass-ratio inspirals (EMRIs).
• To eliminate regularization-induced errors from approximating the Dirac delta by directly incorporating it into the numerical scheme.
• To achieve global spectral accuracy, even at the location of the point-particle source, using a discontinuous Galerkin method with proper flux treatment.
• To enable direct computation of waveforms at future null infinity via the hyperboloidal layer method, avoiding artificial boundary conditions.

Proposed method

• The scalar field is expanded in spherical harmonics, reducing the 3+1 Teukolsky equation to a 1+1 system in tortoise coordinate and time.
• The second-order wave equation is recast as a first-order symmetric hyperbolic system using variables πₗₘ = -∂ψₗₘ/∂τ and φₗₘ = ∂ψₗₘ/∂ρ.
• A multi-domain discontinuous Galerkin method is applied, with numerical fluxes specifically derived to correctly handle the Dirac delta source term.
• The hyperboloidal layer method is used to map future null infinity to a finite boundary, enabling direct access to the far-field waveform with trivial outer boundary conditions.
• The method ensures strong hyperbolicity and stability, with careful treatment of the wavespeeds at the right boundary (ρ = s), where one wave speed vanishes.
• Two different first-order reduction choices are compared, with one found to be numerically unstable despite theoretical equivalence.

Experimental results

Research questions

• RQ1 Can a discontinuous Galerkin method achieve global spectral convergence for the Teukolsky equation with a Dirac delta source term, even at the singularity location?
• RQ2 Does the proposed symmetric hyperbolic formulation with proper flux treatment eliminate regularization errors from approximating the delta function?
• RQ3 How does the choice of first-order reduction variables affect numerical stability and convergence in the presence of a point-particle source?
• RQ4 Can the hyperboloidal layer method be effectively combined with DG to extract accurate waveforms at future null infinity without artificial boundary conditions?
• RQ5 What is the behavior of the late-time tail of the scalar field, and does the scheme correctly capture this power-law decay?

Key findings

• The proposed DG method achieves global spectral accuracy, even at the location of the Dirac delta source, due to the correct treatment of the singular source in the numerical flux.
• Convergence tests against exact solutions confirm exponential convergence rates, validating the method’s high-order accuracy.
• The scheme correctly captures the late-time power-law tail behavior of the scalar field, consistent with analytical expectations.
• The hyperboloidal layer method successfully enables direct computation of waveforms at future null infinity with trivial boundary conditions.
• One choice of first-order reduction variables leads to numerical instability and loss of accuracy after ~4–5 digits, while the other maintains full convergence.
• Energy luminosities for circular orbits computed with the method show agreement with known results, confirming physical consistency.