Remnants

gwModelRem family (Islam, Wadekar & Khanna 2026)

Unified remnant models for aligned-spin, precessing and eccentric binaries, plus the point-particle limit. See arXiv:2608.00934.

gwModelRemS: non-precessing quasi-circular

gwModels.remnants.gwModelRemS.gwModelRemS_mf(q, chi1z, chi2z)[source]

Final mass Mf/M of a non-precessing quasi-circular BBH merger.

Ansatz:

Mf/M = 1 - eta * E_rad_hat
E_rad_hat = (1-4 eta) [1 - E_ISCO(chi_hat)]
            + 4 eta E_EM(chi_hat, chi_a)
            + eta (1-4 eta) R_E(eta, chi_hat)
            + 4 eta delta_m chi_a (u0 + u1 chi_hat + u2 eta)

Parameters:

q – Mass ratio m1/m2 >= 1, scalar or array.
chi1z – Aligned spin of the heavier BH, in [-1, 1].
chi2z – Aligned spin of the lighter BH, in [-1, 1].

Returns:

Final mass in units of the total mass.

Return type:

float or array

gwModels.remnants.gwModelRemS.gwModelRemS_chif(q, chi1z, chi2z)[source]

Final spin chi_f,z of a non-precessing quasi-circular BBH merger.

Ansatz:

chi_f,z = S_tilde + eta * ell_orb
S_tilde = (q^2 chi1z + chi2z)/(1+q)^2 = (1 - 2 eta) chi_hat
ell_orb = (1-4 eta) ell_Kerr(chi_hat)
          + 4 eta ell_EM(chi_hat, chi_a)
          + eta (1-4 eta) R_chi(eta, chi_hat)
          + 4 eta delta_m chi_a (v0 + v1 chi_hat + v2 eta)

Parameters:

q – Mass ratio m1/m2 >= 1, scalar or array.
chi1z – Aligned spin of the heavier BH, in [-1, 1].
chi2z – Aligned spin of the lighter BH, in [-1, 1].

Returns:

Final dimensionless spin along z (signed).

Return type:

float or array

gwModels.remnants.gwModelRemS.gwModelRemS_Lpeak(q, chi1z, chi2z)[source]

Peak GW luminosity of a non-precessing quasi-circular BBH merger.

Ansatz:

L_peak = eta^2 * L_hat
L_hat  = L_EM(chi_hat, chi_a) + (1-4 eta) P_L(eta, chi_hat)
         + delta_m chi_a (nu0 + nu1 chi_hat + nu2 eta)

The eta^2 prefactor enforces L_peak -> 0 in the point-particle limit.

Parameters:

q – Mass ratio m1/m2 >= 1, scalar or array.
chi1z – Aligned spin of the heavier BH, in [-1, 1].
chi2z – Aligned spin of the lighter BH, in [-1, 1].

Returns:

Peak luminosity in geometric units (c^5/G).

Return type:

float or array

gwModels.remnants.gwModelRemS.gwModelRemS_omega_peak(q, chi1z, chi2z)[source]

Peak GW frequency M*omega at the peak of abs(h22).

Ansatz:

omega_peak = f_EM + (1-4 eta) P + delta_m chi_a U

Parameters:

q – Mass ratio m1/m2 >= 1, scalar or array.
chi1z – Aligned spin of the heavier BH, in [-1, 1].
chi2z – Aligned spin of the lighter BH, in [-1, 1].

Returns:

M*omega at peak abs(h22) (dimensionless).

Return type:

float or array

gwModels.remnants.gwModelRemS.gwModelRemS_kick(q, chi1z, chi2z)[source]

Recoil velocity of a non-precessing quasi-circular BBH merger.

Ansatz:

v_kick   = sqrt(V_mass^2 + V_spin^2 + 2 V_mass V_spin cos(xi))
V_mass   = A eta^2 delta_m (1 + B eta + C eta^2)
V_spin   = H eta^2 R_v
xi       = (pi/180)(a_deg + b_deg S_tilde_k + c_deg delta_m Delta_tilde)

The point-particle and equal-mass-equal-spin recoils both vanish, so only the mass-asymmetry and spin-asymmetry vectors contribute. The eta^2 scaling enforces v_kick -> 0 as eta -> 0, and delta_m = 0 removes V_mass at q = 1.

Parameters:

q – Mass ratio m1/m2 >= 1, scalar or array.
chi1z – Aligned spin of the heavier BH, in [-1, 1].
chi2z – Aligned spin of the lighter BH, in [-1, 1].

Returns:

Recoil velocity in km/s.

Return type:

float or array

gwModels.remnants.gwModelRemS.gwModelRemS(q, chi1z, chi2z)[source]

All gwModelRemS remnant properties for a non-precessing quasi-circular BBH.

Parameters:

q – Mass ratio m1/m2 >= 1, scalar or array.
chi1z – Aligned spin of the heavier BH, in [-1, 1].
chi2z – Aligned spin of the lighter BH, in [-1, 1].

Returns:

(Mf, chif, Lpeak, omega_peak, vkick) where

Mf: final mass in units of the total mass chif: final dimensionless spin along z (signed) Lpeak: peak GW luminosity in geometric units (c^5/G) omega_peak: M*omega at peak abs(h22) (dimensionless) vkick: recoil velocity in km/s

Return type:

tuple

Example

>>> import gwModels
>>> Mf, chif, Lp, wp, vk = gwModels.remnants.gwModelRemS(3.0, 0.5, -0.2)

gwModelRemP: precessing quasi-circular

gwModels.remnants.gwModelRemP.spin_projections(q, a1, a2, theta1, theta2)[source]

Aligned projections and in-plane spin combinations at r = 8M.

Parameters:

q – Mass ratio m1/m2 >= 1.
a1 – Spin magnitudes at r = 8M.
a2 – Spin magnitudes at r = 8M.
theta1 – Spin tilt angles at r = 8M in radians.
theta2 – Spin tilt angles at r = 8M in radians.

Returns:

(chi1z, chi2z, S_perp, Delta_perp).

Return type:

tuple

gwModels.remnants.gwModelRemP.gwModelRemP(q, a1, a2, theta1, theta2, phi1, phi2)[source]

gwModelRemP remnant properties for a precessing quasi-circular BBH merger.

The aligned-spin baseline is evaluated at the z-projections of the spins at r = 8M, then augmented with in-plane-spin corrections.

Parameters:

q – Mass ratio m1/m2 >= 1, scalar or array.
a1 – Spin magnitude of the heavier BH at r = 8M, in [0, 1].
a2 – Spin magnitude of the lighter BH at r = 8M, in [0, 1].
theta1 – Tilt angle of the heavier BH spin at r = 8M in radians, [0, pi].
theta2 – Tilt angle of the lighter BH spin at r = 8M in radians, [0, pi].
phi1 – Azimuthal spin angles at r = 8M in radians. Accepted for interface completeness; the deterministic models do not use them.
phi2 – Azimuthal spin angles at r = 8M in radians. Accepted for interface completeness; the deterministic models do not use them.

Returns:

(Mf, af_mag, theta_f, Lpeak) where

Mf: final mass in units of the total mass af_mag: final spin magnitude, in [0, 1] theta_f: final spin tilt from the orbital angular momentum, radians Lpeak: peak GW luminosity in geometric units (c^5/G)

Return type:

tuple

Note

The recoil of a precessing binary is not modeled deterministically. Use gwModelRemP_flow for the recoil velocity distribution.

Example

>>> import numpy as np, gwModels
>>> Mf, af, thf, Lp = gwModels.remnants.gwModelRemP(
...     2.0, 0.7, 0.3, np.pi/3, np.pi/4, 0.0, 0.0)

gwModels.remnants.gwModelRemP.gwModelRemP_mf(q, a1, a2, theta1, theta2, phi1, phi2)[source]

Final mass Mf/M of a precessing quasi-circular BBH merger.

See gwModelRemP for the parameter description.

Returns:: Final mass in units of the total mass.
Return type:: float or array

gwModels.remnants.gwModelRemP.gwModelRemP_chif(q, a1, a2, theta1, theta2, phi1, phi2)[source]

Final spin magnitude and tilt of a precessing quasi-circular BBH merger.

See gwModelRemP for the parameter description.

Returns:: (af_mag, theta_f) with the tilt in radians.
Return type:: tuple

gwModels.remnants.gwModelRemP.gwModelRemP_Lpeak(q, a1, a2, theta1, theta2, phi1, phi2)[source]

Peak GW luminosity of a precessing quasi-circular BBH merger.

See gwModelRemP for the parameter description.

Returns:: Peak luminosity in geometric units (c^5/G).
Return type:: float or array

gwModelRemSE: eccentric non-precessing

gwModels.remnants.gwModelRemSE.gwModelRemSE_mf(q, chi1z, chi2z, e_ref, l_ref)[source]

Final mass Mf/M of an eccentric non-precessing BBH merger.

Ansatz:

Mf|_ecc = Mf^circ [1 + P_M (1 + alpha_M e sin(l + phi_M))]

Parameters:

q – Mass ratio m1/m2 >= 1, scalar or array.
chi1z – Aligned spin components, in [-1, 1].
chi2z – Aligned spin components, in [-1, 1].
e_ref – Eccentricity at t = -2500M, in [0, 1).
l_ref – Mean anomaly at t = -2500M in radians.

Returns:

Final mass in units of the total mass.

Return type:

float or array

gwModels.remnants.gwModelRemSE.gwModelRemSE_chif(q, chi1z, chi2z, e_ref, l_ref)[source]

Final spin chi_f,z of an eccentric non-precessing BBH merger.

Ansatz:

chif|_ecc = chif^circ [1 + P_chi (1 + alpha_chi e cos(l + phi_chi))]

Parameters:

q – Mass ratio m1/m2 >= 1, scalar or array.
chi1z – Aligned spin components, in [-1, 1].
chi2z – Aligned spin components, in [-1, 1].
e_ref – Eccentricity at t = -2500M, in [0, 1).
l_ref – Mean anomaly at t = -2500M in radians.

Returns:

Final dimensionless spin along z (signed).

Return type:

float or array

gwModels.remnants.gwModelRemSE.gwModelRemSE_kick(q, chi1z, chi2z, e_ref, l_ref)[source]

Recoil velocity of an eccentric non-precessing BBH merger.

Ansatz:

vk|_ecc = vk^circ [1 + P_k (1 + alpha_k e cos(l + phi_k))]

The correction is multiplicative, so it vanishes at q = 1 for non-spinning binaries because the circular recoil itself vanishes there.

Parameters:

q – Mass ratio m1/m2 >= 1, scalar or array.
chi1z – Aligned spin components, in [-1, 1].
chi2z – Aligned spin components, in [-1, 1].
e_ref – Eccentricity at t = -2500M, in [0, 1).
l_ref – Mean anomaly at t = -2500M in radians.

Returns:

Recoil velocity in km/s.

Return type:

float or array

gwModels.remnants.gwModelRemSE.gwModelRemSE_Lpeak(q, chi1z, chi2z, e_ref, l_ref)[source]

Peak GW luminosity of an eccentric non-precessing BBH merger.

Ansatz:

Lp|_ecc = Lp^circ [1 + P_L (1 - alpha_L e sin(l + phi_L))]

Parameters:

q – Mass ratio m1/m2 >= 1, scalar or array.
chi1z – Aligned spin components, in [-1, 1].
chi2z – Aligned spin components, in [-1, 1].
e_ref – Eccentricity at t = -2500M, in [0, 1).
l_ref – Mean anomaly at t = -2500M in radians.

Returns:

Peak luminosity in geometric units (c^5/G).

Return type:

float or array

gwModels.remnants.gwModelRemSE.gwModelRemSE(q, chi1z, chi2z, e_ref, l_ref)[source]

All gwModelRemSE remnant properties for an eccentric non-precessing BBH.

Parameters:

q – Mass ratio m1/m2 >= 1, scalar or array.
chi1z – Aligned spin components, in [-1, 1].
chi2z – Aligned spin components, in [-1, 1].
e_ref – Eccentricity at t = -2500M, in [0, 1).
l_ref – Mean anomaly at t = -2500M in radians.

Returns:

(Mf, chif, vkick, Lpeak) where

Mf: final mass in units of the total mass chif: final dimensionless spin along z (signed) vkick: recoil velocity in km/s Lpeak: peak GW luminosity in geometric units (c^5/G)

Return type:

tuple

Note

Calibrated for q <= 4 and e0 <= 0.25 on non-spinning systems, and provisional: within that domain it is neutral to slightly worse than the circular gwModelRemS baseline on NR data. Do not extrapolate past e_ref ~ 0.3. See the module header for measured numbers.

Example

>>> import gwModels
>>> Mf, chif, vk, Lp = gwModels.remnants.gwModelRemSE(
...     2.0, 0.0, 0.0, 0.1, 0.0)

gwModelRemPE: eccentric precessing

gwModels.remnants.gwModelRemPE.gwModelRemPE(q, a1, a2, theta1, theta2, phi1, phi2, e_ref, l_ref)[source]

gwModelRemPE remnant properties for an eccentric precessing BBH merger.

Parameters:

q – Mass ratio m1/m2 >= 1, scalar or array.
a1 – Spin magnitudes at r = 8M, in [0, 1].
a2 – Spin magnitudes at r = 8M, in [0, 1].
theta1 – Spin tilt angles at r = 8M in radians, in [0, pi].
theta2 – Spin tilt angles at r = 8M in radians, in [0, pi].
phi1 – Azimuthal spin angles at r = 8M in radians. Accepted for interface completeness; the deterministic models do not use them.
phi2 – Azimuthal spin angles at r = 8M in radians. Accepted for interface completeness; the deterministic models do not use them.
e_ref – Eccentricity at t = -2500M, in [0, 1).
l_ref – Mean anomaly at t = -2500M in radians.

Returns:

(Mf, af_mag, theta_f, Lpeak) where

Mf: final mass in units of the total mass af_mag: final spin magnitude, in [0, 1] theta_f: final spin tilt from the orbital angular momentum, radians Lpeak: peak GW luminosity in geometric units (c^5/G)

Return type:

tuple

Note

The eccentric corrections were calibrated on non-spinning aligned-spin systems and are applied here to a precessing baseline, treating the two effects as independent at leading order. theta_f carries no eccentric correction. Use gwModelRemP_flow for the recoil.

Example

>>> import numpy as np, gwModels
>>> Mf, af, thf, Lp = gwModels.remnants.gwModelRemPE(
...     2.0, 0.7, 0.3, np.pi/3, np.pi/4, 0.0, 0.0, 0.1, 0.0)

gwModels.remnants.gwModelRemPE.gwModelRemPE_mf(q, a1, a2, theta1, theta2, phi1, phi2, e_ref, l_ref)[source]

Final mass Mf/M of an eccentric precessing BBH merger.

See gwModelRemPE for the parameter description.

Returns:: Final mass in units of the total mass.
Return type:: float or array

gwModels.remnants.gwModelRemPE.gwModelRemPE_chif(q, a1, a2, theta1, theta2, phi1, phi2, e_ref, l_ref)[source]

Final spin magnitude and tilt of an eccentric precessing BBH merger.

See gwModelRemPE for the parameter description.

Returns:: (af_mag, theta_f) with the tilt in radians.
Return type:: tuple

gwModels.remnants.gwModelRemPE.gwModelRemPE_Lpeak(q, a1, a2, theta1, theta2, phi1, phi2, e_ref, l_ref)[source]

Peak GW luminosity of an eccentric precessing BBH merger.

See gwModelRemPE for the parameter description.

Returns:: Peak luminosity in geometric units (c^5/G).
Return type:: float or array

gwModelRemP_flow: probabilistic precessing recoil

class gwModels.remnants.gwModelRemP_flow.gwModelRemP_flow(datadir=None)[source]

Bases: object

Conditional normalizing-flow model for the precessing BBH recoil velocity.

Samples recoil velocities for a given precessing configuration, marginalizing over the spin-orientation information not retained in the context vector.

Context features are computed with the gwModelRemS and gwModelRemP models in this package, so the flow stays consistent with the deterministic remnant models rather than carrying private copies of them.

Parameters:: datadir (str) – Directory containing gwModelRemP_flow.pt. Defaults to the data directory shipped with gwModels.

best_nll

Test negative log likelihood from the checkpoint.

Type:: float or None

Example

>>> import numpy as np, gwModels
>>> flow = gwModels.remnants.gwModelRemP_flow()
>>> med, p5, p95 = flow.predict(2.0, 0.7, 0.3, np.pi/3, np.pi/4, 0.0, 0.0)

compute_context(q, a1, a2, theta1, theta2)[source]

Build the five-dimensional context vector for the flow.

Parameters:

q – Mass ratio m1/m2 >= 1, scalar or array.
a1 – Spin magnitudes at r = 8M, in [0, 1].
a2 – Spin magnitudes at r = 8M, in [0, 1].
theta1 – Spin tilt angles at r = 8M in radians, in [0, pi].
theta2 – Spin tilt angles at r = 8M in radians, in [0, pi].

Returns:

Shape (N, 5), columns Mf, abs(chi_f), eta, S_perp, Delta_perp.

Return type:

array

sample(q, a1, a2, theta1, theta2, phi1, phi2, n_samples=1000)[source]

Draw recoil velocity samples from the conditional distribution.

Parameters:

q – Mass ratio m1/m2 >= 1, scalar or array.
a1 – Spin magnitudes at r = 8M, in [0, 1].
a2 – Spin magnitudes at r = 8M, in [0, 1].
theta1 – Spin tilt angles at r = 8M in radians, in [0, pi].
theta2 – Spin tilt angles at r = 8M in radians, in [0, pi].
phi1 – Spin azimuthal angles at r = 8M in radians. Accepted for interface completeness; the context marginalizes over them.
phi2 – Spin azimuthal angles at r = 8M in radians. Accepted for interface completeness; the context marginalizes over them.
n_samples – Number of samples to draw per binary.

Returns:

Shape (N, n_samples) recoil velocities in km/s, where N is: the number of input systems (1 for scalar inputs).

Return type:

array

predict(q, a1, a2, theta1, theta2, phi1, phi2, n_samples=5000)[source]

Point estimate: median recoil velocity and 90% credible interval.

Parameters:

q – Mass ratio m1/m2 >= 1, scalar or array.
a1 – Spin magnitudes at r = 8M, in [0, 1].
a2 – Spin magnitudes at r = 8M, in [0, 1].
theta1 – Spin tilt angles at r = 8M in radians, in [0, pi].
theta2 – Spin tilt angles at r = 8M in radians, in [0, pi].
phi1 – Spin azimuthal angles at r = 8M in radians. Not used.
phi2 – Spin azimuthal angles at r = 8M in radians. Not used.
n_samples – Number of samples used to compute the statistics.

Returns:

(median, p5, p95) recoil velocities in km/s.

Return type:

tuple

gwModelEMRI: point-particle limit

gwModels.remnants.gwModelEMRI.gwModelEMRI(q, chi, theta_inc, e_sep, warn_unconverged=True, return_converged=False)[source]

Point-particle remnant mass and spin for an extreme-mass-ratio inspiral.

Evaluates the leading-order extreme-mass-ratio limit

Mf/M = 1 - eta (1 - E_sep) chi_f = chi + eta (L_sep - 2 chi (E_sep - 1))

with E_sep and L_sep the conserved quantities of a test particle at the Kerr geodesic separatrix. No terms are fitted to NR data.

Parameters:

q – Mass ratio m1/m2 >= 1, scalar or array. Intended for q >> 1000.
chi – Dimensionless spin of the primary Kerr BH, in [-1, 1].
theta_inc – Orbital inclination in radians, in [0, pi]. 0 is prograde equatorial and pi retrograde equatorial.
e_sep – Eccentricity at the separatrix, in [0, 1).
warn_unconverged – Emit a RuntimeWarning if the solver fails anywhere.
return_converged – If True, also return the boolean convergence mask.

Returns:

(Mf, chif) or (Mf, chif, converged) if return_converged is True,: with Mf in units of the total mass and chif the signed final spin.

Return type:

tuple

Note

Results for systems where the solver has not converged are not physically meaningful. Check the convergence mask when scanning high spin at near-polar inclination; see the module header for rates.

Example

>>> import gwModels
>>> Mf, chif = gwModels.remnants.gwModelEMRI(1e4, 0.9, 0.0, 0.0)

Non-precessing kick (Islam & Wadekar 2025)

gwModels.remnants.IW2025_kick_nonprecessing.gwModel_kick_q200(q, chi1z, chi2z, return_std=False)[source]

Kick velocity for aligned-spin binaries. From Islam & Wadekar (2025), https://arxiv.org/abs/2511.11536

Trained on NR (SXS + RIT, q <= 32) and BHPT data (q <= 200). Valid for 1 <= q <= 1000.

Parameters:

q – Mass ratio m1/m2 >= 1
chi1z – Dimensionless spin of primary along z, in [-1, 1]
chi2z – Dimensionless spin of secondary along z, in [-1, 1]
return_std – If True, also return parameter uncertainty estimate

Returns:

Kick velocity in km/s V_kick_std (optional): Estimated uncertainty in km/s

Return type:

V_kick

GPR kick (Islam & Wadekar 2025)

class gwModels.remnants.IW2025_kick_gpr.gwModel_kick_q200_GPR(model_path)[source]

Bases: object

GPR-based aligned-spin kick velocity model. From Islam & Wadekar (2025), https://arxiv.org/abs/2511.11536

Gaussian Process trained on kick data with features [log2(q), chi_hat, chi_a].

For the analytical kick model, use gwModel_kick_q200 from IW2025_kick_nonprecessing.

Parameters:: model_path (str) – Path to the .pkl model file.

predict(q, chi1z, chi2z)[source]

GPR kick prediction.

Parameters:

q – Mass ratio m1/m2 >= 1
chi1z – Dimensionless spin of primary along z, in [-1, 1]
chi2z – Dimensionless spin of secondary along z, in [-1, 1]

Returns:

Kick velocity in km/s vk_std: GPR uncertainty in km/s

Return type:

vk

info()[source]

Precessing kick flow (Islam & Wadekar 2025)

class gwModels.remnants.IW2025_kick_precessing.gwModel_kick_prec_flow(datadir)[source]

Bases: object

Normalizing-flow model for precessing-spin kick velocity distributions.

Samples kick velocities for given (q, a1, a2), marginalizing over spin orientation angles.

Parameters:: datadir (str) – Directory containing gwModel_kick_prec_flow.pt and gwModel_kick_prec_flow_config.npy

sample(q, a1, a2, num_samples=5000)[source]

Sample kick velocities from the flow model.

Accepts scalar or array (q, a1, a2). Array inputs (length N) are evaluated in a single batched call over the flow rather than one binary at a time. For non-spinning systems (a1==0 and a2==0), falls back to the aligned-spin model gwModel_kick_q200 (applied per element).

Parameters:

q – Mass ratio(s) (m1/m2 >= 1), scalar or array
a1 – Primary spin magnitude(s), scalar or array
a2 – Secondary spin magnitude(s), scalar or array
num_samples – Number of samples to draw per binary

Returns:

kick velocities in km/s. Shape conventions:

scalar (q, a1, a2) -> (num_samples,)
array of length N, num_samples=1 -> (N,)
array of length N, num_samples>1 -> (N, num_samples)

Return type:

samples

Note

Batched and per-binary calls are statistically equivalent but not element-wise identical (the RNG is consumed in a different order).

log_prob(v_kick, q, a1, a2)[source]

Conditional log-density log p(v_kick | q, a1, a2) of the kick magnitude.

Evaluates the flow’s density of the recoil magnitude v = abs(x) by folding the two signed branches of the underlying variable and including the normalization Jacobian, so that exp(log_prob) integrates to 1 over v >= 0. All inputs broadcast to a common shape and are scored in a single batched call.

This is the primitive for inverting the model (inferring progenitor q, a1, a2 from an observed kick via Bayes). It represents the smooth flow density and does not apply the non-spinning delta-function fallback used by sample(); intended for spinning systems.

Parameters:

v_kick – Observed kick magnitude(s) [km/s], scalar or array
q – Progenitor parameters, scalar or array
a1 – Progenitor parameters, scalar or array
a2 – Progenitor parameters, scalar or array

Returns:

log-density, scalar if all inputs are scalar else a 1-D array.

info()[source]

HLZ kick

gwModels.remnants.HLZ_kick.spin_difference_vector(small_q, a1, a2, theta1, theta2, delta_phi)[source]

Calculate the antisymmetric spin combination Delta = (S2/m2 - S1/m1)/(1+q).

Parameters:

small_q – mass ratio m2/m1 <= 1 (internal convention)
a1 – dimensionless spin magnitudes
a2 – dimensionless spin magnitudes
theta1 – angles between L and spin vectors (radians)
theta2 – angles between L and spin vectors (radians)
delta_phi – azimuthal angle difference between spins (radians)

Returns:

Delta_parallel, Delta_perp

gwModels.remnants.HLZ_kick.total_spin_vector(small_q, a1, a2, theta1, theta2, delta_phi)[source]

Calculate the symmetric spin combination chi = (S1 + S2)/M^2.

Parameters:

small_q – mass ratio m2/m1 <= 1 (internal convention)
a1 – dimensionless spin magnitudes
a2 – dimensionless spin magnitudes
theta1 – angles between L and spin vectors (radians)
theta2 – angles between L and spin vectors (radians)
delta_phi – azimuthal angle difference between spins (radians)

Returns:

chi_tilde_parallel, chi_tilde_perp

gwModels.remnants.HLZ_kick.calculate_kick_components(small_q, a1, a2, theta1, theta2, delta_phi, Theta=None)[source]

Calculate the three kick velocity components.

Parameters:

small_q – mass ratio m2/m1 <= 1 (internal convention)
a1 – dimensionless spin magnitudes
a2 – dimensionless spin magnitudes
theta1 – angles between L and spin vectors (radians)
theta2 – angles between L and spin vectors (radians)
delta_phi – azimuthal angle difference between spins (radians)
Theta – angle between Delta x L and fiducial infall direction (radians). If None, a random value between 0 and 2pi is used.

Returns:

Vm, Vs_perp, Vs_parallel, Theta_used

gwModels.remnants.HLZ_kick.bbh_final_kick_precessing_CLZM2007(q, a1, a2, theta1, theta2, delta_phi, Theta=None, debug=False)[source]

Calculate total kick velocity magnitude.

Parameters:

q – mass ratio q = m1/m2 >= 1
a1 – dimensionless spin magnitudes
a2 – dimensionless spin magnitudes
theta1 – angles between L and spin vectors (radians)
theta2 – angles between L and spin vectors (radians)
delta_phi – azimuthal angle difference between spins (radians)
Theta – angle between Delta x L and fiducial infall direction (radians). If None, a random value between 0 and 2pi is used.
debug – if True, return all components

Returns:

total kick velocity in km/s If debug=True: V_kick, Vm, Vs_perp, Vs_parallel, Theta_used

Return type:

V_kick (float)

gwModels.remnants.HLZ_kick.bbh_final_kick_nonprecessing_HLZ2014(q, chi1z, chi2z)[source]

RIT aligned-spin recoil (kick) for binaries with spins along +/-z. Coefficients from arXiv:1406.7295.

Parameters:

q – mass ratio q = m1/m2 >= 1
chi1z – dimensionless spin of primary along z in [-1, 1]
chi2z – dimensionless spin of secondary along z in [-1, 1]

Returns:

kick velocity in km/s

Return type:

V_kick

HBR final mass and spin

gwModels.remnants.HBR_mass_spin.energy_at_isco(a)[source]: Dimensionless specific energy at the ISCO: E_ISCO(a).

gwModels.remnants.HBR_mass_spin.angular_momentum_at_isco(a)[source]: Dimensionless specific angular momentum at ISCO: L_ISCO(a).

gwModels.remnants.HBR_mass_spin.angle_between_spins(theta1, theta2, delta_phi)[source]: Angle alpha between the two spin vectors using spherical law of cosines.

gwModels.remnants.HBR_mass_spin.angle_correction(theta, eps)[source]: Angle remapping from Eq. (18): tan(theta’/2) = (1+eps) tan(theta/2).

gwModels.remnants.HBR_mass_spin.bbh_final_mass_precessing_BMR2012(q, a1, a2, theta1, theta2, verbose=False)[source]

Final remnant mass using Barausse-Morozova-Rezzolla (2012) fit.

Parameters:

q – Mass ratio q = m1/m2 >= 1
a1 – Dimensionless spin magnitudes (0 <= a <= 1)
a2 – Dimensionless spin magnitudes (0 <= a <= 1)
theta1 – Angles (radians) between orbital momentum and spins
theta2 – Angles (radians) between orbital momentum and spins
verbose – Print intermediate calculations (default: False)

Returns:

Remnant mass as fraction of initial total mass M_f/M

Return type:

Mfin

Reference: Barausse, Morozova & Rezzolla (2012), ApJ 758, 63, arXiv:1206.3803

gwModels.remnants.HBR_mass_spin.bbh_final_spin_precessing_HBR2016(q, a1, a2, theta1, theta2, delta_phi, model='HBR16_34corr', verbose=False)[source]

Final spin magnitude using Hofmann, Barausse & Rezzolla (2016) fit.

Parameters:

q – Mass ratio q = m1/m2 >= 1
a1 – Dimensionless spin magnitudes (0 <= a <= 1)
a2 – Dimensionless spin magnitudes (0 <= a <= 1)
theta1 – Angles (radians) between orbital momentum and spins
theta2 – Angles (radians) between orbital momentum and spins
delta_phi – Angle between spin projections on orbital plane
model – Fit model selection (default: “HBR16_34corr”)
verbose – Print intermediate calculations

Returns:

Final spin magnitude, abs(a_final) <= 1

Return type:

chi_final

Reference: Hofmann, Barausse & Rezzolla (2016), ApJL 825, L19, arXiv:1605.01938

UIB2016 final mass and spin

gwModels.remnants.UIB2016_mass_spin.bbh_final_mass_non_precessing_UIB2016(q, chi1z, chi2z, version='v2')[source]

Final mass fraction Mf/M using the UIB2016 fit, parameterized by mass ratio.

Parameters:

q – Mass ratio q = m1/m2 >= 1
chi1z – Dimensionless spins along z in [-1, 1]
chi2z – Dimensionless spins along z in [-1, 1]
version – “v1” or “v2” (default “v2”)

Returns:

Final mass as a fraction of total mass

Return type:

Mf_over_M

gwModels.remnants.UIB2016_mass_spin.bbh_final_spin_non_precessing_UIB2016(q, chi1z, chi2z, version='v2')[source]

Final spin using the UIB2016 fit, parameterized by mass ratio.

Parameters:

q – Mass ratio q = m1/m2 >= 1
chi1z – Dimensionless spins along z in [-1, 1]
chi2z – Dimensionless spins along z in [-1, 1]
version – “v1” or “v2” (default “v2”)

Returns:

Final spin magnitude

Return type:

chif

Kerr geodesic quantities

gwModels.remnants.Kerr.clip_spin(chi, limit=0.9999)[source]

Clip a dimensionless spin away from the extremal values +/-1.

At spin magnitude 1 the Kerr ISCO degenerates: r_ISCO -> 1 for prograde orbits, and the denominator sqrt(r^1.5 - 3 sqrt(r) + 2 chi) vanishes, so kerr_isco_angular_momentum divides by zero. Effective spin combinations such as chi_hat reach exactly +/-1 when both component spins are extremal, so calibrated fits clip before evaluating the ISCO.

This is a fit convention, not physics, and is therefore applied explicitly by the models rather than being folded into the functions below.

Parameters:

chi – Dimensionless spin, scalar or array.
limit – Clipping magnitude (default SPIN_CLIP = 0.9999).

Returns:

Spin clipped to [-limit, limit].

Return type:

float or array

gwModels.remnants.Kerr.kerr_isco_radius(a)[source]

Boyer-Lindquist r_ISCO(a) for equatorial orbits of a Kerr black hole.

r_ISCO = 3 + Z2 - sign(a) sqrt((3 - Z1)(3 + Z1 + 2 Z2)) Z1 = 1 + (1 - a^2)^(1/3) [(1 + a)^(1/3) + (1 - a)^(1/3)] Z2 = sqrt(3 a^2 + Z1^2)

Limits: r_ISCO(0) = 6, r_ISCO(1) = 1, r_ISCO(-1) = 9.

Parameters:: a – Dimensionless spin parameter, abs(a) <= 1. Positive for prograde, negative for retrograde orbits.
Returns:: ISCO radius in units of GM/c^2.
Return type:: float or array

gwModels.remnants.Kerr.kerr_isco_energy(a)[source]

Specific energy of a test particle at the Kerr ISCO.

E_ISCO = sqrt(1 - 2/(3 r_ISCO))

Limits: E_ISCO(0) = sqrt(8/9) = 0.942809, E_ISCO(1) = 1/sqrt(3) = 0.577350.

Parameters:: a – Dimensionless spin parameter, abs(a) <= 1.
Returns:: Specific energy at the ISCO.
Return type:: float or array

gwModels.remnants.Kerr.kerr_isco_angular_momentum(a)[source]

Specific axial angular momentum magnitude at the Kerr ISCO.

L_ISCO = (r^2 - 2 d a sqrt(r) + a^2)
/ (r^(3/4) sqrt(r^(3/2) - 3 sqrt(r) + 2 d a))

with r = r_ISCO(a) and d = sign(a). Limits: L_ISCO(0) = 2 sqrt(3) = 3.464102, L_ISCO(1) = 2/sqrt(3) = 1.154701, L_ISCO(-1) = 22/(3 sqrt(3)) = 4.232809.

Note that the denominator vanishes at extremal spin, where this is a removable singularity: evaluating at exactly +/-1 divides by zero. Pass the spin through clip_spin first if extremal values are possible. An equivalent form regular at the endpoints is (2/(3 sqrt(3)))(1 + 2 sqrt(3 r - 2)); the two agree to about 1e-12 and this one is retained because the calibrated remnant fits were built on it.

Parameters:: a – Dimensionless spin parameter, abs(a) <= 1.
Returns:: Specific angular momentum magnitude at the ISCO.
Return type:: float or array

gwModels.remnants.Kerr.kerr_ell(a)[source]

Effective angular momentum combination at the Kerr ISCO.

ell_Kerr = L_ISCO - 2 a (E_ISCO - 1)

This is the combination that enters remnant-spin fits anchored to the point-particle limit, following Hofmann, Barausse & Rezzolla (2016).

Parameters:: a – Dimensionless spin parameter, abs(a) <= 1.
Returns:: Effective angular momentum at the ISCO.
Return type:: float or array

gwModels.remnants.Kerr.separatrix_energy(e_s, chi)[source]

Specific energy at the separatrix, with a phenomenological Kerr extension.

E_sep = sqrt(8/(9 - e_s^2)) + [E_ISCO(chi) - E_ISCO(0)] 9/(9 - e_s^2)

The first term is the exact Schwarzschild separatrix (p_sep = 6 + 2 e_s); the second interpolates toward the Kerr ISCO. Reduces to E_ISCO(chi) in the circular limit e_s -> 0.

Parameters:

e_s – Eccentricity at the separatrix, in [0, 1).
chi – Dimensionless spin parameter.

Returns:

Specific energy at the separatrix.

Return type:

float or array

gwModels.remnants.Kerr.separatrix_angular_momentum(e_s, chi)[source]

Specific angular momentum at the separatrix, with a Kerr extension.

L_sep = (6 + 2 e_s)/sqrt(3 + 2 e_s - e_s^2)

[L_ISCO(chi) - L_ISCO(0)] sqrt(3)/sqrt((3 - e_s)(1 + e_s))

Reduces to L_ISCO(chi) in the circular limit e_s -> 0.

Parameters:

e_s – Eccentricity at the separatrix, in [0, 1).
chi – Dimensionless spin parameter.

Returns:

Specific angular momentum at the separatrix.

Return type:

float or array

gwModels.remnants.Kerr.separatrix_ell(e_s, chi)[source]

Effective angular momentum at the separatrix.

ell_sep = L_sep - 2 chi (E_sep - 1)

Reduces to kerr_ell(chi) in the circular limit e_s -> 0.

Parameters:

e_s – Eccentricity at the separatrix, in [0, 1).
chi – Dimensionless spin parameter.

Returns:

Effective angular momentum at the separatrix.

Return type:

float or array

gwModels.remnants.Kerr.separatrix_EL(a, e, x, max_iter=20, tol=1e-12, return_converged=False)[source]

Energy and axial angular momentum at the generic Kerr separatrix.

Parameters:

a – Dimensionless Kerr spin, in [-1, 1].
e – Orbital eccentricity, in [0, 1). Values below 1e-8 are treated as circular and use the R = R’ = R’’ = 0 branch.
x – Cosine of the orbital inclination, in [-1, 1]. x = 1 is prograde equatorial and x = -1 retrograde equatorial.
max_iter – Maximum Newton iterations.
tol – Convergence tolerance on the residuals of the defining conditions.
return_converged – If True, also return a boolean convergence mask.

Returns:

(E, Lz) or (E, Lz, converged) if return_converged is True.

Return type:

tuple

Note

Check the convergence mask before using results, particularly at high spin combined with near-polar inclination. See the module header.

Utilities

gwModels.remnants.remnant_utils.validate_q(q)[source]: Validate that q = m1/m2 >= 1.

gwModels.remnants.remnant_utils.validate_spin_magnitudes(a1, a2)[source]: Validate that spin magnitudes are in [0, 1].

gwModels.remnants.remnant_utils.validate_spin_z(chi1z, chi2z)[source]: Validate that spin z-components are in [-1, 1].

gwModels.remnants.remnant_utils.symmetric_mass_ratio(q)[source]

Calculate symmetric mass ratio eta(q) = q/(1+q)^2.

eta = 0.25 for equal masses (q = 1), eta -> 0 for q -> infinity. The formula is symmetric under q -> 1/q, so works for any q > 0.

Parameters:: q – Mass ratio q = m1/m2 >= 1.
Returns:: Symmetric mass ratio.
Return type:: float or array