PRECESSION 2.1: black-hole binary spin precession on eccentric orbits

TL;DR

Precession 2.1 is an updated Python tool for simulating the dynamics of precessing black hole binaries, now including eccentric orbits and new features for orbital evolution and gravitational-wave analysis.

Contribution

The paper introduces version 2.1 of the precession code, extending its capabilities to handle eccentric orbits and providing new equations for orbital and gravitational-wave parameters.

Findings

Extended the code to handle eccentric orbits using a semi-automatic method.

Added orbit- and precession-averaged evolution equations for eccentricity.

Revised expressions for converting between separation and gravitational-wave frequency.

Abstract

We present version 2.1 of the public code {\sc precession}, a Python module for studying the post-Newtonian dynamics of precessing black hole binaries. In this release, we extend the code to handle eccentric orbits. This extension leverages the existing numerical infrastructure wherever possible, introducing a semi-automatic method to adapt circular-orbit functions to the eccentric case via a Python decorator. Additional new features include orbit- and precession-averaged evolutionary equations for the eccentricity, as well as revised expressions to convert between post-Newtonian separation and gravitational-wave emission frequency.

Equations

$$a⟶a(1-e^2),t⟶t(1-e^2{)}^{3/2},$$

$$M\omega$$

$${\mathcal{A}}_{1\mathrm{PN}}$$

$${\mathcal{A}}_{1.5\mathrm{PN}}$$

$${\mathcal{A}}_{2\mathrm{PN}}$$

$$+\left(323-216\sqrt{1-e^2}\right)q^3+\left(75-60\sqrt{1-e^2}\right)q^4+e^2[36+60\sqrt{1-e^2}+q(38+$$

$$+216\sqrt{1-e^2})-\left(2-312\sqrt{1-e^2}\right)q^2+\left(38+216\sqrt{1-e^2}\right)q^3+12\left(3+5\sqrt{1-e^2}\right)q^4]+$$

$$+e^4\left(147+563q+835q^2+563q^3+147q^4\right)+3(1+q{)}^2\left(1+e^2\right)[2q{\chi }_1{\chi }_2(2\cos {\theta }_1\cos {\theta }_2+$$

$$-\cos \Delta \Phi \sin {\theta }_1\sin {\theta }_2)+\left(3\cos {\theta }_1^2-1\right){\chi }_1^2+\left(3\cos {\theta }_2^2-1\right)q^2{\chi }_2^2]\},$$

Full text

1 Introduction

The detection of gravitational waves (GWs) emitted during the final orbits of black hole (BH) binaries enables precise measurements of the properties of these systems [1]. Some parameters, such as the BH masses and spin magnitudes, remain constant throughout the binary’s evolution. Others, such as spin orientations and orbital eccentricity, evolve significantly over time [2, 3]. As a result, their values at the time of GW detection generally differ from those at BH formation. While the former class of parameters provides some insight into the origin of BH binaries, it is often the latter that carries the most discriminating information about the formation channels of these systems [4, 5]. Connecting measured quantities to their values at formation [5, 6], or, conversely, to evolve the predictions of population synthesis codes [7, 8, 9, 10] to the regime probed by GW detectors, is essential to gain insight into the formation channels of BH binaries. This requires tracking the evolution of the angular momenta vectors as binaries evolve on (potentially eccentric) orbits [11, 12, 13]. A multi-timescale treatment of BH-binary dynamics forms the foundation of the formalism presented in Ref. [11, 12] and implemented in the numerical code precession [14]. While originally designed to handle spinning BH binaries on quasi-circular orbits, we have recently extended precession to the more general case of eccentric binaries. The formalism underlying this extension has already been published in Ref.[13]. In this short contribution to the proceedings of the 2025 GR+Amaldi conference, we summarize the new features of the precession code, which can now seamlessly perform post-Newtonian (PN) evolutions of spinning, eccentric BH binaries. Our modifications fall into two main categories: those implementing the transformations required to extend the circular formalism to eccentric binaries (Sec. 2), and those that have no circular counterpart and are unique to the eccentric case (Sec. 3). This extension to eccentric binaries constitutes version 2.1 of precession. The code is publicly available at github.com/dgerosa/precession, where users can also find detailed documentation and example tutorials.

2 Mapping eccentric binaries to effective circular sources

On the precessional timescale, the dynamics of an eccentric binary can be mapped onto those of an effective circular source [13]. The orbital separation (equal to the orbital radius $r$ in the circular case) is replaced by the semi-latus rectum, while the time coordinate $t$ must be rescaled by an eccentricity-dependent factor. Together, the mapping reads

[TABLE]

where $a$ is the semi-major axis and $e$ is the orbital eccentricity, ranging from [math] (circle) to $1$ (parabola). This transformation must be applied to all functions in the code that depend on $r$ and $t$, with a few exceptions (see Sec. 3). To implement the mapping of the orbital separation, we developed a utility that dynamically wraps functions in precession using a Python decorator. This decorator, named eccentricize, modifies the user-facing interface of selected functions by replacing the argument $r$ with $a$ and $e$. Internally, the functions continue to operate in terms of $r$, which is computed from $a$ and $e$ via Eq. (1), thereby preserving their original logic while extending compatibility to eccentric binaries. The decorator also updates function signatures and docstrings to reflect the new input parameters. Instead, functions that explicitly depend on time have been manually modified, as they required additional adjustments

3 Functions requiring special treatment

Not all functions in precession can be extended to the eccentric case using the eccentricize decorator. Some require more careful treatment, either because they do not have a counterpart in the circular case (e.g., the definition and evolution of eccentricity) or because they involve non-trivial modifications. The inspiral of a BH binary using orbit-averaged equations is implemented in inspiral_orbav. For circular binaries, precession employs the 2PN spin-precession equations together with the 3.5PN orbital-velocity evolution equation (including spin terms up to 2PN; see Sec.V of Ref.[12]). For eccentric binaries, we have implemented the 2PN spin-precession equations and the 3PN evolutionary equations for the semi-major axis and eccentricity (with spin terms up to 2PN) [see Eqs. (C1a)–(C2l) of Ref.[15]]. Three additional angular variables appear in this case: $\Psi_{1}$, $\Psi_{2}$, and $\delta\lambda$[16, 15, 17]. The first two describe the angles between the projections of the spin vectors onto the orbital plane and the periastron line, while the third identifies the position of the periastron, as illustrated in Fig.1. Therefore, in addition to tracking the spin orientations, it is also necessary to evolve $\delta\lambda$. In this version of the code, we adopt the evolution equation given in Eq.(17b) of Ref.[17]. The precession-averaged evolution of eccentric binaries (inspiral_precav) has been thoroughly described in Ref. [13], to which we refer for further details. In particular, we have derived a new regularized version of the Peters’ equation for $\differential a/\differential e$ that is well-suited to numerical integrations in the $e\to 0$ limit. Figure 1 shows an example of an eccentric BH binary evolution. We present the orbit-averaged evolution of the tilt angles $\theta_{1}$, $\theta_{2}$, and $\delta\lambda$ for a binary with initial eccentricity $e_{0}=0.5$, initialized at a large separation ($a=100M$) and evolved down to $a=10M$. For comparison, we also show the evolution of the spin angles for the same binary initialized with $e_{0}=0$ (the circular case). We note that, when setting $e=0$, the eccentric and circular evolutions do not coincide, because the derivative of $e$ is not zero. This is a known feature of orbit-averaged equations [18, 15, 16], which encodes the influence of spins on eccentricity. Another aspect that required careful handling is the conversion between the GW emission frequency $f_{\rm GW}$ and the orbital separation, namely pnseparation_to_gwfrequency and gwfrequency_to_pnseparation. For eccentric binaries, GW emission is more complex than in the circular case [19], as it occurs across multiple harmonics of the orbital frequency. Moreover, the distribution of power among these harmonics depends sensitively on the eccentricity [19, 20, 21]. In the version of the code presented here, we proceed as follows: given the semi-major axis $a$ and eccentricity $e$, we first compute the orbital frequency $\omega$ and then map $\omega$ to $f_{\rm GW}$, either considering the second harmonic, which dominates in the circular limit, or the harmonic that carries the maximum power. The latter is identified using the fitting formulae of Refs. [20, 21], which we have reimplemented. If $\omega$ or $f_{\rm GW}$ at a specific harmonic is known, and the eccentricity $e$ is given, the corresponding semi-major axis can be directly computed using Eqs. (4a) and (B2a) of Ref.[15]. To perform the inverse operation—recovering the frequency from the semi-major axis—we analytically invert Eq. (B2a) of Ref. [15], obtaining

[TABLE]

where the coefficients are defined as

[TABLE]

where $M=m_{1}+m_{2}$ is the total mass of the system, $q=m_{2}/m_{1}$ is the mass ratio ($m_{2}\leq m_{1}$), $\chi_{1}$, $\chi_{2}$ are the spin magnitudes, and $\Delta\Phi$ is the difference between the projection of spin vectors onto the orbital plane. In the previous, circular-only version of precession, the conversions between $r$ and $f_{\rm GW}$ were based on Eqs. (4.5) and (4.13) of Ref. [22]. However, these did not include the spin-induced quadrupole term, a contribution that enters at 2PN order [23, 24]. In the current release of the code, we have corrected the expressions for the circular case, now including this missing contribution. While the impact of this term in the frequency-separation conversion equation is expected to be small, it is nonetheless important for consistency and completeness in spin-precession modeling. precession also includes a legacy implementation of analytic fits for the post-merger BH mass, spin, and recoil [14]. These fits were developed only for quasi-circular binaries and, therefore, have not been ported to the eccentric version of the code.

precession is distributed via git (github.com/dgerosa/precession); installation is as easy as pip install precession. Currently, the default version of the code remains that for quasi-circular binaries, as described in Ref. [12]; this is accessible via import precession. The functionalities for studying eccentric binaries described here and in Ref. [13] are implemented in a dedicated submodule, which we recommend importing as import precession.eccentricity as precession.

Acknowledgments

We thank Matteo Boschini for discussions. G.F., D.G., and N.L. are supported by ERC Starting Grant No. 945155–GWmining, Cariplo Foundation Grant No. 2021-0555, MUR PRIN Grant No. 2022-Z9X4XS, Italian-French University (UIF/UFI) Grant No. 2025-C3-386, MUR Grant “Progetto Dipartimenti di Eccellenza 2023-2027” (BiCoQ), and the ICSC National Research Centre funded by NextGenerationEU. D.G. is supported by MSCA Fellowship No. 101064542–StochRewind, MSCA Fellowship No. 101149270–ProtoBH and MUR Young Researchers Grant No. SOE2024-0000125. Computational work was performed at CINECA with allocations through INFN and Bicocca, and at NVIDIA with allocations through the Academic Grant program.