Functions | |
def | L_vec_from_masses_spins_and_f (m1, m2, f, PNorder=2, S1_vec=np.zeros(3), S2_vec=np.zeros(3), no_mass_spin_input_check=False) |
Compute the PN orbital angular momentum from the binary's masses and spins and m = 2 GW frequency using the expression in Eq. More... | |
def | L_from_masses_spins_and_f (m1, m2, f, PNorder=2, S1_vec=np.zeros(3), S2_vec=np.zeros(3), no_mass_spin_input_check=False) |
Compute the magnitude of the PN orbital angular momentum from the binary's masses and spins and m = 2 GW frequency using the expression in Eq. More... | |
def | L_from_masses_a_and_ecc (m1, m2, a, e) |
Compute the magnitude of the Newtonian orbital angular momentum from the binary's masses, semimajor axis, and eccentricity. More... | |
def | E_ratio_term (m) |
Compute m*E_ratio(m) [cf. More... | |
def | eq_coeffs (u, kappaxiq, q, S1, S2, xi, S0sq) |
Compute coefficients for the cubic equation in bar{S}^2 to be solved [Eq. More... | |
def | solve_cubic (coeffs, u, kappaxiq, q, S1, S2, imag_tol, use_mpmath=False, root_tol=1.e-6, polyroots_maxsteps=50, polyroots_extraprec=50) |
Solve the cubic equation, check for complex roots, and order the roots. More... | |
def | dkappaxiq_du (u, kappaxiq, q, S1, S2, xi, S0sq, imag_tol, root_tol, lin_tol, root_tol_raise_err, use_mpmath=False, polyroots_maxsteps=50, polyroots_extraprec=50) |
Compute (<S^2>prec - S0^2)/(1 - q) [= dkappa{xi q}/du] using Eq. More... | |
def | kappaxiq_evol (q, S1, S2, S1L, xi, S0sq, u0, uf, du, method, imag_tol, root_tol, lin_tol, ode_atol, ode_rtol, root_tol_raise_err, use_solve_ivp=False, solve_ivp_lsoda_min_step=0., use_mpmath=False, polyroots_maxsteps=50, polyroots_extraprec=50, odefun_degree=None) |
Compute kappaxiq when the binary has u = uf starting from u = u0. More... | |
def | prec_avg_tilt_comp_vec_inputs (m1, m2, chi1_vec, chi2_vec, fref, L0_vec=None, Lf=None, du=1.e-3, ode_atol_base=1.e-8, ode_atol_inplane_spin_scaling=True, ode_atol_floor=1.e-13, ode_rtol=-1, method="lsoda", use_solve_ivp=False, solve_ivp_lsoda_min_step=0., use_fallback=True, du_fallback=1.e-5, ode_atol_fallback=1.e-13, ode_rtol_fallback=-1, method_fallback="lsoda", use_mpmath_fallback=True, mpmath_dps=30, ode_tol_mpmath=1.e-15, odefun_degree=None, polyroots_maxsteps=50, polyroots_extraprec=50, LPNorder=2, LPNspins=True, imag_tol=1.e-6, root_tol=1.e-6, lin_tol=-1, lin_tol_fallback=-1, root_tol_raise_err=False, failure_mode='None', version='v1') |
Compute tilt angles at infinite separation or bounds and average value of tilt angles at finite separation from masses, spins, tilt angles, and in-plane spin angle at some frequency, taken to be small enough that the precession averaged spin evolution from Gerosa et al. More... | |
def | prec_avg_tilt_comp (m1, m2, chi1, chi2, tilt1, tilt2, phi12, fref, Lf=None, **kwargs) |
Compute tilt angles at infinite separation or bounds and average value of tilt angles at finite separation from masses, spins, tilt angles, and in-plane spin angle at some frequency, taken to be small enough that the precession averaged spin evolution from Gerosa et al. More... | |
Variables | |
bool | solve_ivp_avail = False |
bool | mpmath_avail = False |
list | LPN_orders = [0, 1, 1.5, 2, 2.5] |
string | LPN_order_string = ', '.join(map(str, LPN_orders)) |
def lalsimulation.tilts_at_infinity.calc_tilts_prec_avg_regularized.L_vec_from_masses_spins_and_f | ( | m1, | |
m2, | |||
f, | |||
PNorder = 2 , |
|||
S1_vec = np.zeros(3) , |
|||
S2_vec = np.zeros(3) , |
|||
no_mass_spin_input_check = False |
|||
) |
Compute the PN orbital angular momentum from the binary's masses and spins and m = 2 GW frequency using the expression in Eq.
(4) of https://dcc.ligo.org/public/0122/T1500554/023/dLdS.pdf.
We orbit-average the contributions from the spins as in Eq. (17) of that reference, but the spins are set to zero by default.
Inputs:
m1, m2: Detector frame masses in kg f: m = 2 GW frequency in Hz
PNorder: Order of PN corrections to include–choices are 0, 1, 1.5, 2, or 2.5, where the half-integer orders are the spin-orbit contributions, default: 2 S1_vec, S2_vec: 3d numpy arrays giving the components of the dimensionful spins (with total mass = 1) in a coordinate system where the Newtonian orbital angular momentum is in the z-direction, default: [0., 0., 0.] no_mass_spin_input_check: Turn off mass and spin input checking for when this is called by other functions that have already performed this check, default: False
Output: 3d array giving the PN orbital angular momentum in total mass = 1 units in a coordinate system where the Newtonian orbital angular momentum is in the z-direction
Definition at line 72 of file calc_tilts_prec_avg_regularized.py.
def lalsimulation.tilts_at_infinity.calc_tilts_prec_avg_regularized.L_from_masses_spins_and_f | ( | m1, | |
m2, | |||
f, | |||
PNorder = 2 , |
|||
S1_vec = np.zeros(3) , |
|||
S2_vec = np.zeros(3) , |
|||
no_mass_spin_input_check = False |
|||
) |
Compute the magnitude of the PN orbital angular momentum from the binary's masses and spins and m = 2 GW frequency using the expression in Eq.
(4) of https://dcc.ligo.org/public/0122/T1500554/023/dLdS.pdf.
We orbit-average the contributions from the spins as in Eq. (17) of that reference, but the spins are set to zero by default.
This wraps L_vec_from_masses_spins_and_f()
Inputs:
m1, m2: Detector frame masses in kg f: m = 2 GW frequency in Hz
PNorder: Order of PN corrections to include–choices are 0, 1, 1.5, 2, or 2.5, where the half-integer orders are the spin-orbit contributions, default: 2 S1_vec, S2_vec: 3d numpy arrays giving the components of the dimensionful spins (with total mass = 1) in a coordinate system where the Newtonian orbital angular momentum is in the z-direction, default: [0., 0., 0.] no_mass_spin_input_check: Turn off mass and spin input checking for when this is called by other functions that have already performed this check, default: False
Output: Magnitude of the PN orbital angular momentum in total mass = 1 units
Definition at line 170 of file calc_tilts_prec_avg_regularized.py.
def lalsimulation.tilts_at_infinity.calc_tilts_prec_avg_regularized.L_from_masses_a_and_ecc | ( | m1, | |
m2, | |||
a, | |||
e | |||
) |
Compute the magnitude of the Newtonian orbital angular momentum from the binary's masses, semimajor axis, and eccentricity.
Inputs:
m1, m2: Source frame masses in kg a: Semimajor axis in m e: eccentricity
Output: Magnitude of orbital angular momentum in total mass = 1 units
Definition at line 188 of file calc_tilts_prec_avg_regularized.py.
def lalsimulation.tilts_at_infinity.calc_tilts_prec_avg_regularized.E_ratio_term | ( | m | ) |
Compute m*E_ratio(m) [cf.
Eq. (B13) in the paper https://dcc.ligo.org/P2100029, arXiv:2107.11902, which gives m^2*E_ratio(m)]
Definition at line 214 of file calc_tilts_prec_avg_regularized.py.
def lalsimulation.tilts_at_infinity.calc_tilts_prec_avg_regularized.eq_coeffs | ( | u, | |
kappaxiq, | |||
q, | |||
S1, | |||
S2, | |||
xi, | |||
S0sq | |||
) |
Compute coefficients for the cubic equation in bar{S}^2 to be solved [Eq.
(20) in the paper https://dcc.ligo.org/P2100029, arXiv:2107.11902]
Inputs: u: 1/(2L), where L is the orbital angular momentum (i.e., the integration variable) kappaxiq: The rescaled variable defined in Eq. (14) of the paper q: Mass ratio of the binary, with the < 1 convention S1, S2: Magnitudes of the dimensionful spins of the binary (with total mass = 1) xi: Conserved effective spin, defined in Eq. (2) of the paper S0sq: Square of magnitude of the initial total spin (at u = u0, and with total mass = 1)
Output: The coefficients of the equation (vareps, Btransf, Ctransf, Dtransf)
Definition at line 235 of file calc_tilts_prec_avg_regularized.py.
def lalsimulation.tilts_at_infinity.calc_tilts_prec_avg_regularized.solve_cubic | ( | coeffs, | |
u, | |||
kappaxiq, | |||
q, | |||
S1, | |||
S2, | |||
imag_tol, | |||
use_mpmath = False , |
|||
root_tol = 1.e-6 , |
|||
polyroots_maxsteps = 50 , |
|||
polyroots_extraprec = 50 |
|||
) |
Solve the cubic equation, check for complex roots, and order the roots.
The inputs besides the coefficients are just needed in order to be printed if one of the checks is triggered
Inputs:
coeffs: Array of coefficients of cubic [vareps, Btransf, Ctransf, Dtransf] u: 1/(2L), where L is the orbital angular momentum (i.e., the integration variable) kappaxiq: The rescaled variable defined in Eq. (14) of the paper https://dcc.ligo.org/P2100029, arXiv:2107.11902 q: Mass ratio of the binary, with the < 1 convention S1, S2: Magnitudes of the dimensionful spins of the binary (with total mass = 1) imag_tol: Tolerance for the imaginary part of roots of the cubic
use_mpmath: If True, use mpmath.polyroots to solve the cubic instead of numpy's Polynomial roots method (slower but more accurate), default: False root_tol: Tolerance for the error estimate from mpmath.polyroots, if use_mpmath == True, default: 1.e-6 polyroots_maxsteps: The maximum number of steps allowed in mpmath.polyroots if use_mpmath == True, default: 50 polyroots_extraprec: The number of extra digits to use in mpmath.polyroots if use_mpmath == True, default: 50
Output: The roots of the equation ordered from smallest to largest (Ssq3, Ssqm, Ssqp)
Definition at line 302 of file calc_tilts_prec_avg_regularized.py.
def lalsimulation.tilts_at_infinity.calc_tilts_prec_avg_regularized.dkappaxiq_du | ( | u, | |
kappaxiq, | |||
q, | |||
S1, | |||
S2, | |||
xi, | |||
S0sq, | |||
imag_tol, | |||
root_tol, | |||
lin_tol, | |||
root_tol_raise_err, | |||
use_mpmath = False , |
|||
polyroots_maxsteps = 50 , |
|||
polyroots_extraprec = 50 |
|||
) |
Compute (<S^2>prec - S0^2)/(1 - q) [= dkappa{xi q}/du] using Eq.
(42) in Chatziioannou, Klein, Yunes, and Cornish PRD 95, 104004 (2017) [rewritten as bar{S}_+^2 + [E(m)/K(m) - 1](bar{S}_+^2 - bar{S}_3^2), where bar{S}_*^2 := (S_*^2 - S_0^2)/(1 - q)]. Here kappa_{xi q} := (kappa - S_0^2/(2 L) - q xi/(1 + q))/(1 - q) - xi/4, and S_0 is the magnitude of the initial total spin (i.e., at u = u0).
Inputs:
u: 1/(2L), where L is the orbital angular momentum (i.e., the integration variable) kappaxiq: The rescaled variable defined above and in Eq. (14) of the paper https://dcc.ligo.org/P2100029, arXiv:2107.11902 q: Mass ratio of the binary, with the < 1 convention S1, S2: Magnitudes of the dimensionful spins of the binary (with total mass = 1) xi: Conserved effective spin, defined in Eq. (2) of the paper https://dcc.ligo.org/P2100029, arXiv:2107.11902 S0sq: Square of magnitude of the initial total spin (at u = u0, and with total mass = 1) imag_tol: Tolerance for the imaginary part of roots of the cubic root_tol: Tolerance for the error in the roots of the cubic (is an absolute tolerance for small values of the roots, as generally occurs during the beginning of the evolution, and is a fractional tolerance when the roots are large, which generally occurs towards the end of the evolution; also is reduced automatically when the roots are close together) lin_tol: Tolerance for errors on the r.h.s. of dkappa/du due to linearization in m (is a fractional tolerance when m is large, as generally occurs during the beginning of the evolution, and is an absolute tolerance when m is small, as generally occurs towards the end of the evolution) root_tol_raise_err: Whether to raise an error (True) or just print a warning (False) when the error in the roots of the cubic is larger than root_tol
use_mpmath: If True, use mpmath functions instead of numpy functions to solve the cubic and evaluate the elliptic integrals (slower but more accurate), default: False polyroots_maxsteps: The maximum number of steps allowed in mpmath.polyroots if use_mpmath == True, default: 50 polyroots_extraprec: The number of extra digits to use in mpmath.polyroots if use_mpmath == True, default: 50
Output: (<S^2>prec - S0^2)/(1 - q) [= dkappa{xi q}/du]
Definition at line 395 of file calc_tilts_prec_avg_regularized.py.
def lalsimulation.tilts_at_infinity.calc_tilts_prec_avg_regularized.kappaxiq_evol | ( | q, | |
S1, | |||
S2, | |||
S1L, | |||
xi, | |||
S0sq, | |||
u0, | |||
uf, | |||
du, | |||
method, | |||
imag_tol, | |||
root_tol, | |||
lin_tol, | |||
ode_atol, | |||
ode_rtol, | |||
root_tol_raise_err, | |||
use_solve_ivp = False , |
|||
solve_ivp_lsoda_min_step = 0. , |
|||
use_mpmath = False , |
|||
polyroots_maxsteps = 50 , |
|||
polyroots_extraprec = 50 , |
|||
odefun_degree = None |
|||
) |
Compute kappaxiq when the binary has u = uf starting from u = u0.
Inputs:
q: Mass ratio of the binary, with the < 1 convention S1, S2: Magnitudes of the dimensionful spins of the binary (with total mass = 1) S1L: Component of the dimensionful spin of hole 1 along the orbital angular momentum (with total mass = 1) xi: Conserved effective spin, defined in Eq. (2) of the paper https://dcc.ligo.org/P2100029, arXiv:2107.11902 S0sq: Square of initial dimensionful total spin (in total mass = 1 units) u0: Initial value of u uf: Final value of u du: Step in u for the evolution method: method to use for scipy integration ("vode", "lsoda", "dopri5", or "dop853" for the ode interface and "RK45", "RK23", "DOP853", or "LSODA" for the solve_ivp interface); not used if use_mpmath == True imag_tol: Tolerance for the imaginary part of roots of the cubic root_tol: Tolerance for the error in the roots of the cubic (is an absolute tolerance for small values of the roots, as generally occurs during the beginning of the evolution, and is a fractional tolerance when the roots are large, which generally occurs towards the end of the evolution; also is reduced automatically when the roots are close together) lin_tol: Tolerance for errors on the r.h.s. of dkappa/du due to linearization in m (is a fractional tolerance when m is large, as generally occurs during the beginning of the evolution, and is an absolute tolerance when m is small, as generally occurs towards the end of the evolution) ode_atol: Absolute tolerance for the default scipy ode solver, also the single tolerance setting for the mpmath ode solver if it is used ode_rtol: Relative tolerance for the default scipy ode solver root_tol_raise_err: Whether to raise an error (True) or just print a warning (False) when the error in the roots of the cubic is larger than root_tol
use_solve_ivp: If True, use scipy's solve_ivp function for the ODE integration, which gives access to a different set of integrators and allows tighter tolerances for LSODA than the default ode function, but can hang for certain parameters, default: False solve_ivp_lsoda_min_step: Minimum step to use in the solve_ivp LSODA integrator, default: 0.
use_mpmath: If True, use mpmath functions and ode integration instead of numpy/scipy (slower but more accurate), default: False polyroots_maxsteps: The maximum number of steps allowed in mpmath.polyroots if use_mpmath == True, default: 50 polyroots_extraprec: The number of extra digits to use in mpmath.polyroots if use_mpmath == True, default: 50 odefun_degree: Degree of polynomial used by mpmath.odefun if use_mpmath == True; calculated automatically by mpmath.odefun if None, default: None
Output: kappaxiq at Lf
Definition at line 706 of file calc_tilts_prec_avg_regularized.py.
def lalsimulation.tilts_at_infinity.calc_tilts_prec_avg_regularized.prec_avg_tilt_comp_vec_inputs | ( | m1, | |
m2, | |||
chi1_vec, | |||
chi2_vec, | |||
fref, | |||
L0_vec = None , |
|||
Lf = None , |
|||
du = 1.e-3 , |
|||
ode_atol_base = 1.e-8 , |
|||
ode_atol_inplane_spin_scaling = True , |
|||
ode_atol_floor = 1.e-13 , |
|||
ode_rtol = -1 , |
|||
method = "lsoda" , |
|||
use_solve_ivp = False , |
|||
solve_ivp_lsoda_min_step = 0. , |
|||
use_fallback = True , |
|||
du_fallback = 1.e-5 , |
|||
ode_atol_fallback = 1.e-13 , |
|||
ode_rtol_fallback = -1 , |
|||
method_fallback = "lsoda" , |
|||
use_mpmath_fallback = True , |
|||
mpmath_dps = 30 , |
|||
ode_tol_mpmath = 1.e-15 , |
|||
odefun_degree = None , |
|||
polyroots_maxsteps = 50 , |
|||
polyroots_extraprec = 50 , |
|||
LPNorder = 2 , |
|||
LPNspins = True , |
|||
imag_tol = 1.e-6 , |
|||
root_tol = 1.e-6 , |
|||
lin_tol = -1 , |
|||
lin_tol_fallback = -1 , |
|||
root_tol_raise_err = False , |
|||
failure_mode = 'None' , |
|||
version = 'v1' |
|||
) |
Compute tilt angles at infinite separation or bounds and average value of tilt angles at finite separation from masses, spins, tilt angles, and in-plane spin angle at some frequency, taken to be small enough that the precession averaged spin evolution from Gerosa et al.
PRD 92, 064016 (2015) [implemented following Chatziioannou, Klein, Yunes, and Cornish PRD 95, 104004 (2017)] is accurate, using the regularized expressions from https://dcc.ligo.org/P2100029, arXiv:2107.11902. This version takes vectors of spins as input (as opposed to spin angles) and allows one optionally to specify the initial orbital angular momentum instead of computing it from fref (and the binary's parameters).
Inputs:
m1, m2: Detector frame masses of the binary, in kg chi1_vec, chi2_vec: 3d numpy arrays with the dimensionless spins of the binary; these must be in coordinates where the Newtonian orbital angular momentum is in the z-direction if L0_vec is None (the default), otherwise all that is necessary is that chi1_vec, chi2_vec, and L0_vec are all in the same coordinate system fref: Reference frequency, in Hz; ignored when L0_vec is not None
L0_vec: 3d numpy array with the initial orbital angular momentum in units where the binary's total mass = 1 and in the same coordinates as the spins; this is computed from fref and the binary parameters if None is passed and in that case given in coordinates where the Newtonian orbital angular momentum is in the z-direction, default: None Lf: Final magnitude of orbital angular momentum in total mass = 1 units for output at finite separation, None gives the output at infinity, default: None
du: Step in u for the evolution (external to the ODE integrator), default: 1.e-3 ode_atol_base: Base value of absolute tolerance for ode solver, default: 1.e-8 ode_atol_inplane_spin_scaling: Scale ode_atol_base by min(chi1*np.sin(tilt1),chi2*np.sin(tilt2)), with a floor of ode_atol_floor, to give good accuracy for small in-plane spin cases, default: True ode_atol_floor: Floor value for absolute tolerance for ode solver used when ode_atol_inplane_spin_scaling == True; should be several orders of magnitude less than ode_atol_base, default: 1.e-13 ode_rtol: Relative tolerance for ode solver, -1 sets this to 0.1*ode_atol, default: -1 method: method to use for integration (either from scipy's ode interface, where "vode", "lsoda", "dopri5", and "dop853" have been tested, though "lsoda" is the recommended one; "RK45", "RK23", "DOP853", or "LSODA" for scipy's solve_ivp interface, if use_solve_ivp == True ("lsoda" and "dop853" are equivalent to "LSODA" and "DOP853" here, for compatibility); or "odefun" to use mpmath.odefun with its default Taylor method, which is all that is available as of v1.2.0; this uses the mpmath fallback settings and disables the fallback evolution–it is also quite slow), default: "lsoda" use_solve_ivp: If True, use scipy's solve_ivp interface and its integrators instead of the default ode interface; note that the solve_ivp integrators apparently hang for certain difficult cases, default: False solve_ivp_lsoda_min_step: Minimum step to use in the solve_ivp LSODA integrator, default: 0. use_fallback: Whether to use the fallback integration if the primary integration fails (True) or not (False), not available if method == "odefun", default: True du_fallback: Step in u for the evolution if the first evolution fails, default: 1.e-5 ode_atol_fallback: Absolute tolerance for ode solver if the first evolution fails, default: 1.e-13 ode_rtol_fallback: Relative tolerance for ode solver if the first evolution fails, -1 sets this to 0.1*ode_atol_fallback, default: -1 method_fallback: method to use for integration if the first evolution fails (any of the scipy integrators available for method), default: "lsoda" use_mpmath_fallback: Use the slow but accurate mpmath integration if the first fallback integration fails (True) or not (False), default: True mpmath_dps: The number of decimal places to use for the mpmath computations, if the mpmath integration is triggered or method == "odefun", default: 30 ode_tol_mpmath: Tolerance for the mpmath ode solver, only used if the mpmath fallback evolution is triggered or method == "odefun", default: 1.e-15 odefun_degree: Degree of polynomial used by mpmath.odefun if the mpmath fallback evolution is triggered or method == "odefun"; calculated automatically by mpmath.odefun if None, default: None polyroots_maxsteps: The maximum number of steps allowed in mpmath.polyroots if the mpmath fallback evolution is triggered or method == "odefun", default: 50 polyroots_extraprec: The number of extra digits to use in mpmath.polyroots if the mpmath fallback evolution is triggered or method == "odefun", default: 50 LPNorder: PN order to use in computing the intial orbital angular momentum from fref–choices are 0, 1, 1.5, 2, or 2.5, where the half-integer orders are the spin contributions; ignored when L0_vec is not None, default: 2 LPNspins: Whether to include the orbit-averaged spin contributions in the computation of the initial orbital angular momentum; ignored when L0_vec is not None, default: True imag_tol: Tolerance for the imaginary part of roots of the cubic, default: 1.e-6 root_tol: Tolerance for the error in the roots of the cubic (is an absolute tolerance for small values of the roots, as generally occurs during the beginning of the evolution, and is a fractional tolerance when the roots are large, which generally occurs towards the end of the evolution; also is reduced automatically when the roots are close together), default: 1.e-6 lin_tol: Tolerance for errors on the r.h.s. of dkappa/du due to linearization in m (is a fractional tolerance when m is large, as generally occurs during the beginning of the evolution, and is an absolute tolerance when m is small, as generally occurs towards the end of the evolution); -1 sets this to be the same as ode_atol, except for the mpmath evolution, where it is always set to the ode_tol_mpmath value, default: -1 lin_tol_fallback: Tolerance for errors on the r.h.s. of dkappa/du due to linearization in m if the first evolution fails, -1 sets this to be the same as ode_atol_fallback, default: -1 root_tol_raise_err: Whether to raise an error (True) or just print a warning (False) when the error in the roots of the cubic is larger than root_tol, default: False failure_mode: How the code behaves when the evolution fails (after the fallback evolution, if use_fallback is True, the default). 'Error' means that the code raises a RuntimeError, while 'NAN' and 'None' mean that the code raises a RuntimeWarning and returns np.nan or None for the output, respectively, default: 'None' version: Version of the calculation to use–currently only the initial version, v1, is available, default: "v1"
Output: dictionary with entries 'tilt1_inf', 'tilt2_inf' for evolution to infinity and entries 'tilt1_sep_min', 'tilt1_sep_max', 'tilt1_sep_avg', 'tilt2_sep_min', 'tilt2_sep_max', 'tilt2_sep_avg' for evolution to a finite separation (i.e., a finite orbital angular momentum)
Definition at line 870 of file calc_tilts_prec_avg_regularized.py.
def lalsimulation.tilts_at_infinity.calc_tilts_prec_avg_regularized.prec_avg_tilt_comp | ( | m1, | |
m2, | |||
chi1, | |||
chi2, | |||
tilt1, | |||
tilt2, | |||
phi12, | |||
fref, | |||
Lf = None , |
|||
** | kwargs | ||
) |
Compute tilt angles at infinite separation or bounds and average value of tilt angles at finite separation from masses, spins, tilt angles, and in-plane spin angle at some frequency, taken to be small enough that the precession averaged spin evolution from Gerosa et al.
PRD 92, 064016 (2015) [implemented following Chatziioannou, Klein, Yunes, and Cornish PRD 95, 104004 (2017)] is accurate, using the regularized expressions. This takes in spin angles and a reference frequency to initialize the evolution.
Inputs:
m1, m2: Detector frame masses of the binary, in kg chi1, chi2: Dimensionless spin magnitudes of the binary tilt1, tilt2: Tilt angles of the binary's spins (w.r.t. the Newtonian orbital angular momentum) at fref phi12: Angle between the in-plane components of the spins at fref fref: Reference frequency, in Hz
Lf: Final magnitude of orbital angular momentum in total mass = 1 units for output at finite separation None gives the output at infinity, default: None
du: Step in u for the evolution (external to the ODE integrator), default: 1.e-3 ode_atol_base: Base value of absolute tolerance for ode solver, default: 1.e-8 ode_atol_inplane_spin_scaling: Scale ode_atol_base by min(chi1*np.sin(tilt1),chi2*np.sin(tilt2)), with a floor of ode_atol_floor, to give good accuracy for small in-plane spin cases, default: True ode_atol_floor: Floor value for absolute tolerance for ode solver used when ode_atol_inplane_spin_scaling == True; should be several orders of magnitude less than ode_atol_base, default: 1.e-13 ode_rtol: Relative tolerance for ode solver, -1 sets this to 0.1*ode_atol, default: -1 method: method to use for integration (either from scipy's ode interface, where "vode", "lsoda", "dopri5", and "dop853" have been tested, though "lsoda" is the recommended one; "RK45", "RK23", "DOP853", or "LSODA" for scipy's solve_ivp interface, if use_solve_ivp == True ("lsoda" and "dop853" are equivalent to "LSODA" and "DOP853" here, for compatibility); or "odefun" to use mpmath.odefun with its default Taylor method, which is all that is available as of v1.2.0; this uses the mpmath fallback settings and disables the fallback evolution–it is also quite slow), default: "lsoda" use_solve_ivp: If True, use scipy's solve_ivp interface and its integrators instead of the default ode interface; note that the solve_ivp integrators apparently hang for certain difficult cases, default: False solve_ivp_lsoda_min_step: Minimum step to use in the solve_ivp LSODA integrator, default: 0. use_fallback: Whether to use the fallback integration if the primary integration fails (True) or not (False), not available if method == "odefun", default: True du_fallback: Step in u for the evolution if the first evolution fails, default: 1.e-5 ode_atol_fallback: Absolute tolerance for ode solver if the first evolution fails, default: 1.e-13 ode_rtol_fallback: Relative tolerance for ode solver if the first evolution fails, -1 sets this to 0.1*ode_atol_fallback, default: -1 method_fallback: method to use for integration if the first evolution fails (any of the scipy integrators available for method), default: "lsoda" use_mpmath_fallback: Use the slow but accurate mpmath integration if the first fallback integration fails (True) or not (False), default: True mpmath_dps: The number of decimal places to use for the mpmath computations, if the mpmath integration is triggered or method == "odefun", default: 30 ode_tol_mpmath: Tolerance for the mpmath ode solver, only used if the mpmath fallback evolution is triggered or method == "odefun", default: 1.e-15 odefun_degree: Degree of polynomial used by mpmath.odefun if the mpmath fallback evolution is triggered or method == "odefun"; calculated automatically by mpmath.odefun if None, default: None polyroots_maxsteps: The maximum number of steps allowed in mpmath.polyroots if the mpmath fallback evolution is triggered or method == "odefun", default: 50 polyroots_extraprec: The number of extra digits to use in mpmath.polyroots if the mpmath fallback evolution is triggered or method == "odefun", default: 50 LPNorder: PN order to use in computing the intial orbital angular momentum from fref–choices are 0, 1, 1.5, 2, or 2.5, where the half-integer orders are the spin contributions, default: 2 LPNspins: Whether to include the orbit-averaged spin contributions in the computation of the initial orbital angular momentum, default: True imag_tol: Tolerance for the imaginary part of roots of the cubic, default: 1.e-6 root_tol: Tolerance for the error in the roots of the cubic (is an absolute tolerance for small values of the roots, as generally occurs during the beginning of the evolution, and is a fractional tolerance when the roots are large, which generally occurs towards the end of the evolution; also is reduced automatically when the roots are close together), default: 1.e-6 lin_tol: Tolerance for errors on the r.h.s. of dkappa/du due to linearization in m (is a fractional tolerance when m is large, as generally occurs during the beginning of the evolution, and is an absolute tolerance when m is small, as generally occurs towards the end of the evolution); -1 sets this to be the same as ode_atol, except for the mpmath evolution, where it is always set to the ode_tol_mpmath value, default: -1 lin_tol_fallback: Tolerance for errors on the r.h.s. of dkappa/du due to linearization in m if the first evolution fails, -1 sets this to be the same as ode_atol_fallback, default: -1 root_tol_raise_err: Whether to raise an error (True) or just print a warning (False) when the error in the roots of the cubic is larger than root_tol, default: False failure_mode: How the code behaves when the evolution fails (after the fallback evolution, if use_fallback is True, the default). 'Error' means that the code raises a RuntimeError, while 'NAN' and 'None' mean that the code raises a RuntimeWarning and returns np.nan or None for the output, respectively, default: 'None' version: Version of the calculation to use–currently only the initial version, v1, is available, default: "v1"
Output: dictionary with entries 'tilt1_inf', 'tilt2_inf' for evolution to infinity and entries 'tilt1_sep_min', 'tilt1_sep_max', 'tilt1_sep_avg', 'tilt2_sep_min', 'tilt2_sep_max', 'tilt2_sep_avg' for evolution to a finite separation (i.e., a finite orbital angular momentum)
Definition at line 1402 of file calc_tilts_prec_avg_regularized.py.
bool lalsimulation.tilts_at_infinity.calc_tilts_prec_avg_regularized.solve_ivp_avail = False |
Definition at line 26 of file calc_tilts_prec_avg_regularized.py.
bool lalsimulation.tilts_at_infinity.calc_tilts_prec_avg_regularized.mpmath_avail = False |
Definition at line 33 of file calc_tilts_prec_avg_regularized.py.
list lalsimulation.tilts_at_infinity.calc_tilts_prec_avg_regularized.LPN_orders = [0, 1, 1.5, 2, 2.5] |
Definition at line 39 of file calc_tilts_prec_avg_regularized.py.
string lalsimulation.tilts_at_infinity.calc_tilts_prec_avg_regularized.LPN_order_string = ', '.join(map(str, LPN_orders)) |
Definition at line 41 of file calc_tilts_prec_avg_regularized.py.