RelativisticDynamics.jl: Relativistic Spin-Orbital Dynamics in Julia

RelativisticDynamics.jl is an open-source Julia package for relativistic spin-orbital dynamics in the gravitational strong-field of a Kerr spacetime. Existing codes for modelling the dynamics of spinning objects like pulsars in the strong-field regime are generally lacking, since such systems occupy an intermediate regime that is generally overlooked. At the"low"end of this regime there are post-Newtonian descriptions which neglect the influence of the pulsar spin on the underlying spacetime metric ("spin-curvature"coupling). At the"high"end there are the full numerical relativity solutions which are primarily applicable toe two black holes with a mass ratio $\mathcal{O}(1)$, and are computationally intractable for pulsar systems observed over a large number of orbital cycles. RelativisticDynamics.jl aims to bridge this gap by providing a modern, fast code for accurate numerical evolution of spinning relativistic systems via the Mathisson-Papetrou-Dixon formalism. Julia is a modern language that solves the"two language problem", enabling fast dynamic typing and JIT compilation on conjunction with petaflop performance, comparable with numerical languages that are better known in the astrophysics community such as C or Fortran. RelativisticDynamics.jl is written to be fully type flexible, being able to support arbitrary number formats, and fully differentiable via automatic differentiation.


Summary
Relativistic binaries composed of a millisecond pulsar (MSP) orbiting a much more massive (≳ 10 3  ⊙ ), spinning black hole (BH) are exceptional probes for investigating key questions of fundamental physics and astrophysics.Such systems are natural sources of gravitational waves (GWs) in the mHz regime, expected to be detectable by the next generation of space-based GW detectors such as LISA (Thorpe et al., 2019).The associated radio emission from the companion pulsar raises the possibility of an electromagnetic (EM) counterpart, enabling high precision multimessenger measurements to be made.The description of the orbital dynamics of these systems, and the influence on the resultant observed EM and GW signals, is non-trivial.A proper treatment of the spin-orbital dynamics can be derived from the conservation of the energy-momentum tensor   ; = 0 which when expanded into a set of infinite multipole moments leads to a description of the momentum vector   and the spin tensor for affine parameter , 4-velocity   and Riemann curvature tensor    .The system is closed by providing a spin supplementary condition, equivalent to specifying the observerdependent centre of mass.For this work we take the Tulczyjew-Dixon condition (Dixon, 1964;Tulczyjew, 1959)     = 0 (4) Together, equations 2 -4 form the Mathisson-Papetrou-Dixon (MPD) equations (Dixon, 1964;Mathisson, 1937;Papapetrou, 1951), and describe the spin-orbital evolution in a fully consistent way that is applicable to strong field regimes.

Statement of need
RelativisticDynamics.jl is an open-source Julia package for relativistic spin-orbital dynamics in the gravitational strong field for a Kerr spacetime.Existing codes for modelling the dynamics of spinning objects like pulsars in the strong-field regime are generally lacking, since such systems occupy an intermediate regime that is generally overlooked.At the "low" end of this regime there are post-Newtonian or geodesic descriptions (e.g.Damour & Taylor, 1992) which neglect the influence of the pulsar spin on the underlying spacetime metric ("spincurvature" coupling).At the "high" end there is the full Numerical Relativity (NR) solutions (e.g.Andrade et al., 2021) which are primarily applicable to two BHs with a mass ratio (1), and are computationally intractable for these MSP systems which are observed over a large number of orbital cycles.
RelativisticDynamics.jl aims to bridge this gap by providing a modern, fast code for accurate numerical evolution of spinning relativistic systems, via the MPD formalism.Julia is a modern language that solves the "two language problem", enabling fast dynamic typing and JIT compilation in conjunction with petaflop performance, comparable with numerical languages that are better known in the astrophysics community such as C or Fortran.As a modern language, it also provides a dedicated package manager and a large catalogue of composable packages for scientific computing.This enables RelativisticDynamics.jl to easily leverage and interface with other scientific computing packages.The author and collaborators have used the general methods and mathematics described in this package for multiple research projects (e.g.Kimpson et al., 2019Kimpson et al., , 2020aKimpson et al., , 2020b;;Li et al., 2019) with a particular focus on the radio signals from spinning pulsar systems.This package represents an attempt to create a documented, well-tested, open source resource for public use in this area, that can also be used as a computational playground for exploring techniques that could be applicable to more advanced numerical models.The package has been formulated in terms of ODE integration, rather than using e.g.action-angle variables (Witzany, 2022), to allow for extension to general spacetime metrics and straightforward computation of quantities relevant for pulsar observations e.g.spin axis orientation.
In addition to providing a fast, modern package for strong field spin dynamics, RelativisticDynamics.jl has two additional important features from the perspective of modern relativistic astrophysics.Firstly, it is fully type flexible, being able to support arbitrary number formats.By making use of Julia's type-flexibility the model is written in such a way so as to be able to support hardware accelerated, low precision arithmetic and alternative rounding methods such as stochastic rounding.This enables rapid prototyping and exploration of reduced precision numerical techniques in astrophysics, an approach common in other numerical fields such as weather and climate modelling (e.g.Váňa et al., 2017).Secondly, RelativisticDynamcis.jl is written to be fully differentiable via automatic differentiation.This enables the package to be used for differentiable physics applications in astrophysics, for example gravitational waveform modelling and parameter estimation or training neural networks based on the model.Automatic differentiation also provides a potential avenue for extension of the package to general (i.e.non-Kerr) spacetimes, whereby a user can specify the metric and the associated Christoffel symbols and Riemann tensors -which are simply linear combinations of the metric derivatives -are calculated automatically.
Future potential extensions of this code include taking the dynamics beyond second order in the multipole expansion, and the inclusion of alternative spin conditions and spacetime metrics.The inclusion of a diagnostics tool for extracting gravitational waveforms in the time domain via a numerical kludge method would also be a worthwhile addition.Moreover, we have considered only bound dynamical systems -the ability to also explore hyberbolic systems would also be an interesting development.The pulsar is initialised in the orbital plane with zero inclination.In the absence of spin-curvature coupling the particle would remain in the plane ( = 0).Note the -motion is on the scale of km, not gravitational radii.

Figure 1 :
Figure 1:Example orbital trajectories for a ms-pulsar with eccentricity  = 0.1 (left panels),  = 0.8 (right panels), orbiting a massive BH with extremal spin,  = 0.998.The orbital motion is presented in the  −  plane (top panels) and  −  plane (bottom panels).The pulsar is initialised in the orbital plane with zero inclination.In the absence of spin-curvature coupling the particle would remain in the plane ( = 0).Note the -motion is on the scale of km, not gravitational radii.