PN expressions

𝓔(pnsystem)
binding_energy(pnsystem)

Compute the gravitational binding energy of a compact binary.

Note that this may not be as useful as its derivative, 𝓔′, which is used as part of the right-hand side for orbital evolutions.

The nonspinning orbital binding energy is known through 4pN. The expressions through 3.5pN here come from Eq. (233) of Blanchet (2014).

The 4pN term from Eq. (5.2d) of Jaranowski and Schäfer is known exactly, now that the $ν$-linear piece is given as Eq. (32) of Bini and Damour (2013a). The remaining terms are not known exactly, but Bini and Damour (2013b) have derived some terms, though there is incomplete information, which are noted as the constants in this code.

The spin-orbit terms in the energy are now complete to 4.0pN (the last term is zero). These terms come from Eq. (4.6) of Bohé et al. (2012).

The spin-squared terms (by which we mean both spin-spin and spin-orbit squared terms) in the energy are known to 3pN order, and given in Eq. (3.33) of Bohé et al. (2015).

The tidal-coupling terms come in to the binding energy at relative 5pN order, and are known to 6pN order. These terms come from Eq. (2.11) of Vines et al. (2011). Note their unusual convention for mass ratios, where $χ₁ = m₁/m$ in their notation; in particular, $χ$ is not a spin parameter. Also note that $λ̂ = λ₂ v^{10}/(m₁+m₂)^5$, and we need to add the coupling terms again with $1 ↔ 2$. Finally, note the normalization difference, where a different overall factor is used, leading to a sign difference.

𝓔′(pnsystem)
binding_energy_deriv(pnsystem)

Compute the derivative with respect to $v$ of the binding energy of a compact binary.

This is computed automatically (via FastDifferentiation) from 𝓔; see that function for details of the PN formulas.

𝓕(pnsystem)
gw_energy_flux(pnsystem)

Compute the gravitational-wave energy flux to infinity

The nonspinning flux terms are complete to 4.5pN order, and are given in Eq. (6.11) of Blanchet et al. (2023).

The spin-orbit terms in the flux are now known to 4.0pN. These terms come from Eq. (4.9) of Marsat et al. (2013)

The spin-squared terms (by which we mean both spin-spin and spin-orbit squared terms) in the flux are known to 3pN order, and given most conveniently in Eq. (4.14) of Bohé et al. (2015).

Beyond 4.5pN, terms are only known in the extreme-mass-ratio limit. These terms are given in Appendix A of Fujita (2012). He computed them up to 22pN. That seems like overkill, so we'll just go up to 6pN.

For systems with matter, the tidal-heating terms come in at relative 5pN order, and are known partially at 6pN order. These terms come from Eq. (3.6) of Vines et al. (2011). Note their unusual convention for mass ratios, where $χ₁ = m₁/m$ in their notation; in particular, $χ$ is not a spin parameter. Also note that $λ̂ = λ₂ v^{10}/(m₁+m₂)^5$, and we need to add the coupling terms again with $1 ↔ 2$. Finally, note the normalization difference, where a different overall factor is used, leading to a sign difference.

tidal_heating(pnsystem)

Compute the rate of energy and angular-momentum absorption into each black hole in a binary.

The returned quantity is a tuple (Ṡ₁, Ṁ₁, Ṡ₂, Ṁ₂), representing the corresponding rates of change of spin (magnitude) and mass on black hole 1 and black hole 2. These apply to general — possibly precessing — non-eccentric binaries. This collection of terms comes from Alvi (2001). It probably wouldn't be too hard to extend Alvi's analysis to eccentric systems.

Note that the validity of the result depends not only on the PN parameter $v$, but also on the angles of the spins relative to the separation vector $n̂$: the smaller the angle, the lower the $v$ at which the approximations should be expected to break down.

See also

• Tagoshi and Sasaki (1994)
• Poisson and Sasaki (1995)
• Tagoshi et al. (1997)
• Porto (2007)
• Chatziioannou et al. (2012)

See the documentation section on "Horizons" for details about the computation of horizon quantities used in this function.

r(pnsystem)
separation(pnsystem)

Compute the separation between the two black holes. This is essentially the inverse of γₚₙ, with some factors of G and M thrown in.

Note that there should be a factor of 1/c^2 in this expression; we reserve it to use explicitly in PN expansions. That is, for every factor of 1/r, we explicitly include a factor of 1/c^2 in the expansion.

r′(pnsystem)
separation_deriv(pnsystem)

Compute the derivative of the separation between the two black holes with respect to v.

r⁻¹(r, pnsystem)
separation_inverse(r, pnsystem)

Return v such that r = r(v) when pnsystem is evaluated at v.

Note that the value of v in the input pnsystem is ignored; you may use any value. It may also be convenient to know that you can set the value of v in pnsystem to the returned value using PostNewtonian.vindex as in

pnsystem.state[PostNewtonian.vindex] = r⁻¹(r, pnsystem)

See also γₚₙ⁻¹.

γₚₙ(pnsystem)
inverse_separation(pnsystem)

Compute the post-Newtonian parameter

\[\gamma_{\mathrm{PN}} \equiv \frac{G\, M}{r\, c^2},\]

where $r$ is the magnitude of the orbital separation. This quantity has PN order 1, and is given by Eq. (4.3) of Bohé et al. (2013) and Eq. (3.32) of Bohé et al. (2015).

Note that there is a 3PN gauge term of $-22ν\ln(r/r₀')/3$ that is simply ignored here, as it should cancel out of any physical quantity.

γₚₙ′(pnsystem)
inverse_separation_deriv(pnsystem)

Compute the derivative of γₚₙ with respect to v.

γₚₙ⁻¹(γ, pnsystem)
inverse_separation_inverse(γ, pnsystem)

Return v such that γₚₙ(pnsystem) = γ when pnsystem is evaluated at v.

Note that the value of v in the input pnsystem is ignored; you may use any value. It may also be convenient to know that you can set the value of v in pnsystem to the returned value using PostNewtonian.vindex as in

pnsystem.state[PostNewtonian.vindex] = γₚₙ⁻¹(γ, pnsystem)

See also r⁻¹.

This module contains a few expressions from Kidder (1995).

This is mostly here for testing, because these expressions are not directly used in this package: they are somewhat outdated and describe quantities that are not actually used in this formulation. However, they were used in the SpEC code as an initial guess for eccentricity reduction, so we want to make sure that results from this package are consistent with those from SpEC.

r(pnsystem)

Eq. (4.13).

ṙ(pnsystem)

Eq. (4.12), computed as ṙ = (dE/dt) / (dE/dr), re-expanded and truncated.

aₗ(pnsystem)

Eq. (4.4) of Bohé et al. (2013). This term contributes to 𝛡.

Ω⃗ᵪ(Mⱼ, Mₖ, χ⃗ⱼ, χ⃗ₖ, R)

Compute the angular velocity of precession of spin vector χ⃗ⱼ.

In the approximation that the spin magnitude is constant, the time derivative of $χ⃗ⱼ$ is due to its rotation alone, and is given by $Ω⃗ᵪ × χ⃗ⱼ$.

Note that this function is called by Ω⃗ᵪ₁ and Ω⃗ᵪ₂ with the appropriate parameters; you probably want to use those instead of this one.

The spin-spin term is given by Eq. (2.4) of Kidder (1995); the spin-orbit terms by Eq. (4.5) of Bohé et al. (2013); and the quadrupole-monopole term by Eq. (2.7) Racine (2008).

Ω⃗ᵪ₁(pnsystem)

Compute the angular velocity of precession of χ⃗₁

In the approximation that the spin magnitude is constant, the time derivative of $χ⃗₁$ is due to its rotation alone, and is given by $Ω⃗ᵪ₁ × χ⃗₁$.

Note that this function simply calls Ω⃗ᵪ with the appropriate parameters.

Ω⃗ᵪ₂(pnsystem)

Compute the angular velocity of precession of χ⃗₂

In the approximation that the spin magnitude is constant, the time derivative of $χ⃗₂$ is due to its rotation alone, and is given by $Ω⃗ᵪ₂ × χ⃗₂$.

Note that this function simply calls Ω⃗ᵪ with the appropriate parameters.

Ω⃗ₚ(pnsystem)

Compute the angular velocity of orbital precession.

This is the angular velocity of the orbital angular velocity direction unit vector $ℓ̂$; the time derivative of that unit vector is $Ω⃗ₚ × ℓ̂$.

At the moment, this is computed solely by expressions from Bohé et al. (2013). See 𝛡 for details.

See also R.

𝛡(pnsystem)

Compute the angular velocity of orbital precession according to Bohé et al.

As Bohé et al. (2013) explain above their Eq. (4.1), the orbital precession is given by the time derivative of the orbital axis: 𝓵̇ = 𝛡 × 𝓵, where the angular velocity is along the separation vector 𝓷, so that 𝛡 = ϖ 𝓷. And in turn, they define aₗ ≔ r ω ϖ, where r is the separation and ω is the orbital angular frequency. Then, they define the PN parameter γₚₙ≔M/r and we have Mω = v³ so that ϖ = γₚₙ aₗ / v³. The parameters γₚₙ and aₗ are given by Eqs. (4.3) and (4.4), and given here by the functions γₚₙ and aₗ.

The spin-squared terms (by which we mean both spin-spin and spin-orbit squared terms) in the energy are known to 3pN order, and given in Eq. (3.32) of Bohé et al. (2015).

See also R.

h!(h, pnsystem; ℓₘᵢₙ=0, ℓₘₐₓ=typemax(Int))
mode_weights!(h, pnsystem; ℓₘᵢₙ=0, ℓₘₐₓ=typemax(Int))

Compute mode weights of gravitational waves emitted by pn system, modifying h in place.

These modes are computed in the "co-orbital" frame, in which the larger object lies on the positive $x$ axis, the smaller lies on the negative $x$ axis, and the instantaneous angular velocity is in the positive $z$ direction.

The modes are stored in h in order of increasing $ℓ$ and increasing $m$, with $m$ iterating fastest, all the way up to the highest available mode, $(8,8)$.

Because gravitational waves have spin weight -2, the $(ℓ,m)=(0,0)$, $(1,-1)$, $(1,0)$, and $(1,1)$ modes are always 0. By default, we assume that these modes are nonetheless included in h. If that is not the case, set ℓₘᵢₙ to the smallest $ℓ$ value that should be present in the output data — ℓₘᵢₙ=2 being the most reasonable alternative.

These results come most directly from Eqs. (A5) of Boyle et al. (2014), with the exception of errors in the 2PN spin-spin terms, in which cases we must multiply by $\nu/2$ and make the substitutions $S_1 \mapsto S_0^+$ and $S_2 \mapsto S_0^-$. In turn, those expressions are synthesized from the following: Non-spinning terms are taken from Blanchet (2014), except for the highest-pN term in the (2,±2) mode, which are taken from Blanchet et al. (2023), and the $m=0$ modes, which are taken from Favata (2008). The 1PN spin-orbit term is from Eq. (3.22d) of Kidder (1995). The 1.5PN spin-orbit term is from Eq. (3.22f) of Kidder (1995) and Eq. (F15b) of Will and Wiseman (1996). The 2PN spin-orbit term is from Eq. (4.13) of Buonanno, Faye, Hinderer (2013), while the 2PN spin-spin term is from Eq. (4.15) of that reference.