pespace: A new tool of GPU-accelerated and auto-differentiable response generation and likelihood evaluation for space-borne gravitational wave detectors

Space-borne gravitational wave detectors will expand the scope of gravitational wave astronomy to the milli-Hertz band in the near future. The development of data analysis software infrastructure at the current stage is crucial. In this paper, we introduce \texttt{pespace} which can be used for the full Bayesian parameter estimation of massive black hole binaries with detectors including LISA, Taiji, and Tianqin. The core computations are implemented using the high-performance parallel programming framework \texttt{taichi-lang} which enables automatic differentiation and hardware acceleration across different architectures. We also reimplement the waveform models \texttt{PhenomXAS} and \texttt{PhenomXHM} in the separate package \texttt{tiwave} to integrate waveform generation within the \texttt{taichi-lang} scope, making the entire computation accelerated and differentiable. To demonstrate the functionality of the tool, we use a typical signal from a massive black hole binary to perform the full Bayesian parameter estimation with the complete likelihood function for three scenarios: including a single detector using the waveform with only the dominant mode; a single detector using the waveform including higher modes; and a detector network with higher modes included. The results demonstrate that higher modes are essential in breaking degeneracies, and coincident observations by the detector network can significantly improve the measurement of source properties. Additionally, automatic differentiation provides an accurate way to obtain the Fisher matrix without manual fine-tuning of the finite difference step size. Using a subset of extrinsic parameters, we show that the approximated posteriors obtained by the Fisher matrix agree well with those derived from Bayesian parameter estimation.

💡 Research Summary

The paper introduces pespace, a novel software framework designed for full Bayesian parameter estimation of massive black‑hole binary (MBHB) signals observed by space‑based gravitational‑wave detectors such as LISA, Taiji, and Tianqin. The core of pespace is built on the high‑performance parallel programming language taichi‑lang, which provides just‑in‑time compilation, automatic differentiation (AD), and hardware acceleration across CPUs, GPUs, and other architectures. By leveraging taichi‑lang, the authors re‑implement the state‑of‑the‑art frequency‑domain waveform models PhenomXAS (dominant mode) and PhenomXHM (higher‑order modes) in a separate package called tiwave. This integration allows waveform generation, detector response calculation, and likelihood evaluation to be performed entirely within the taichi‑lang environment, ensuring that all steps are both GPU‑accelerated and differentiable.

The paper first reviews the distinctive features of space‑based detectors: long arm lengths, orbital motion, and frequency‑dependent transfer functions, which render the long‑wavelength approximation insufficient. The response model adopts the stationary‑phase approximation to construct frequency‑domain responses that incorporate these effects.

For Bayesian inference, the authors employ the full Whittle likelihood under the assumption of stationary Gaussian noise, using a noise‑free likelihood for computational convenience. They conduct three benchmark analyses on a representative MBHB signal: (1) a single detector with only the dominant mode, (2) a single detector including higher‑order modes, and (3) a network of LISA, Taiji, and Tianqin with higher‑order modes. The results demonstrate that (i) higher‑order modes dramatically break degeneracies among intrinsic parameters (mass ratio, spins) and extrinsic parameters (sky location, inclination), and (ii) a multi‑detector network substantially increases the signal‑to‑noise ratio, tightening constraints on distance and sky localization.

A key contribution is the use of automatic differentiation to compute the Fisher information matrix directly, avoiding the delicate step‑size tuning required in finite‑difference approaches. By restricting the Fisher analysis to a subset of extrinsic parameters, the authors show that the multivariate Gaussian approximation derived from the Fisher matrix reproduces the posterior distributions obtained from stochastic sampling with high fidelity, confirming the validity of the Fisher approximation in the high‑SNR regime.

The authors acknowledge several limitations: the current implementation focuses solely on isolated MBHB signals and does not yet support global fitting of overlapping sources or other source classes (e.g., EMRIs, stochastic backgrounds). The noise model assumes ideal stationary Gaussian noise, neglecting realistic complications such as non‑stationarity, data gaps, and glitches. The frequency‑domain response relies on the stationary‑phase approximation, limiting applicability to rapidly chirping signals and precessing binaries. Moreover, while AD is available for extrinsic parameters, derivatives with respect to intrinsic parameters still require numerical differencing because waveform generation is not yet AD‑compatible.

Despite these constraints, pespace offers a unique combination of (i) GPU‑accelerated, differentiable waveform and response generation, (ii) support for higher‑order modes, (iii) multi‑detector network capability, and (iv) seamless Fisher matrix computation. The code and data are publicly released on GitHub and Zenodo, facilitating reproducibility and future development. The authors outline a roadmap that includes extending the framework to global fitting, incorporating realistic noise models, adding precessing waveform support, and enabling full AD for intrinsic parameters. In summary, pespace represents a significant step toward efficient, high‑precision data analysis pipelines required for the upcoming era of space‑based gravitational‑wave astronomy.