Trajectory Package

The trajectory package houses modules that turn initial orbital parameters into phase-space trajectories (orbital phases and orbital parameters). These are determined by numerical integrators that interpolate over the relevant quantities.

class few.trajectory.base.TrajectoryBase(*args, **kwargs)

Bases: Citable, ABC

Base class used for trajectory modules.

This class provides a flexible interface to various trajectory implementations. Specific arguments to each trajectory can be found with each associated trajectory module discussed below.

abstractmethod classmethod get_inspiral(*args, **kwargs)

Inspiral Generator

@classmethod that requires a child class to have a get_inspiral method.

Raises:

NotImplementedError – The child class does not have this method.

__call__(*args, in_coordinate_time=True, dt=10.0, T=1.0, new_t=None, spline_kwargs=None, DENSE_STEPPING=False, buffer_length=1000, upsample=False, err=1e-11, fix_t=False, integrate_backwards=False, max_step_size=None, **kwargs)

Call function for trajectory interface.

This is the function for calling the creation of the trajectory. Inputs define the output time spacing.

Parameters:

  • *args – Input of variable number of arguments specific to the inspiral model (see the trajectory class’ get_inspiral method). Important Note: M must be the first parameter of any model that uses this base class.

  • in_coordinate_time (bool) – If True, the trajectory will be outputted in coordinate time. If False, the trajectory will be outputted in units of M. Default is True.

  • dt (float) – Time step for output waveform in seconds. Also sets initial step for integrator. Default is 10.0.

  • T (float) – Total observation time in years. Sets the maximum time for the integrator to run. Default is 1.0.

  • new_t (ndarray | None) – If given, this represents the final time array at which the trajectory is analyzed. This is performed by using a cubic spline on the integrator output. Default is None.

  • spline_kwargs (dict | None) – If using upsampling, spline_kwargs provides the kwargs desired for scipy.interpolate.CubicSpline. Default is {}.

  • DENSE_STEPPING (bool) – If True, the trajectory used in the integrator will be densely stepped at steps of dt. If False, the integrator will determine its stepping. Default is False.

  • buffer_length (int) – Sets the allocation of memory for trajectory parameters. This should be the maximum length expected for a trajectory. If it is reached, output arrays will be extended, but this is more expensive than allocating a larger array initially. Trajectories with default settings will be ~100 points. Default is 1000.

  • upsample (bool) – If True, upsample, with a cubic spline, the trajectories from 0 to T in steps of dt. Default is False.

  • err (float) – Tolerance for integrator. Default is 1e-10. Decreasing this parameter will give more steps over the trajectory, but if it is too small, memory issues will occur as the trajectory length will blow up. We recommend not adjusting this parameter.

  • fix_T – If upsampling, this will affect excess points in the t array. If True, it will shave any excess on the trajectories arrays where the time is greater than the overall time of the trajectory requested.

  • integrate_backwards (bool) – If True, the integrator will run backwards in time. Default is False.

  • max_step_size (float | None) – If given, this will set the maximum step size for the integrator (matching the convention according to in_coordinate_time). Default is None (no maximum step size).

  • **kwargs – kwargs passed to trajectory module. Default is {}.

  • fix_t (bool)

Returns:

Tuple of (t, p, e, Phi_phi, Phi_r, flux_norm).

Raises:

ValueError – If input parameters are not allowed in this model.

Return type:

tuple[ndarray]

classmethod all_citations()

Return all the citations from the registry as printable BibTeX string

Return type:

str

classmethod citation()

Return the module references as a printable BibTeX string.

Return type:

str

classmethod module_references()

Method implemented by each class to define its list of references

Return type:

Iterable[REFERENCE | str]

Generic Inspiral Generator

class few.trajectory.inspiral.EMRIInspiral(*args, func, integrate_constants_of_motion=False, enforce_schwarz_sep=False, convert_to_pex=True, rootfind_separatrix=True, **kwargs)

Bases: TrajectoryBase

EMRI trajectory module.

This module implements generic trajectories by integrating with a DOP853 Runge-Kutta integrator (see few.trajectory.integrate). Both adaptive (default) and fixed timesteps are supported. For an adaptive integration, a continuous solution is generated that can be evaluated at any point in time.

The trajectory operates on generic ODE functions that are defined in few.trajectory.ode. Flux-based trajectories (which interpolate data grids) can be found in few.trajectory.ode.flux. A generic Post-Newtonian trajectory is also provided in few.trajectory.ode.pn5. Users can subclass few.trajectory.ode.base.ODEBase to define their own ODE functions. See the documentation for examples on how to do this.

Parameters:

  • func (str | Type[ODEBase]) – Function name for the ode to use in the integration. This must be given as a keyword argument, even though it is required. To get inbuilt stock options for this argument, check few.trajectory.ode._STOCK_TRAJECTORY_OPTIONS.

  • integrate_constants_of_motion (bool) – If True, the trajectory will integrate the constants of motion (E, L, Q). Default is False.

  • enforce_schwarz_sep (bool) – Enforce the separatrix of Schwarzschild spacetime. This can mitigate issues at higher spin and/or higher eccentricity where (e.g.) PN approximations become increasingly inaccurate. Default is False.

  • convert_to_pex (bool) – Convert the output from ELQ to pex coordinates (only used if integrating constants of motion). Default is True.

  • rootfind_separatrix (bool) – Finish trajectory by performing a numerical root-finding operation to place the final point of the trajectory. If False, performs Euler integration to the final point. Default is True.

  • *args – Any arguments for parent class few.trajectory.base.TrajectoryBase

  • **kwargs – Any keyword arguments for the integrator and ODE classes few.trajectory.integrate.Integrate and few.trajectory.ode.ODEBase.

Raises:

  • ValueError – func kwarg not given or not available.

  • ValueError – File necessary for ODE not found.

inspiral_generator

Integrator class for the trajectory.

Type:

class

specific_kwarg_keys

Specific keywords that need to transferred to the inspiral function that can be adjusted with each call.

Type:

dict

classmethod module_references()

Return citations related to this module

Return type:

list[REFERENCE]

property npoints

Number of points in the trajectory.

property tolerance: float

Absolute tolerance of the integrator.

property dense_stepping: bool

If True, trajectory is using fixed stepping.

property integrate_backwards: bool

If True, integrate backwards.

bounds_p0(e0, x0, a, integrate_backwards=False)

Get bounds on p for given e, x, a.

bounds_e0(p0, x0, a, integrate_backwards=False)

Get bounds on e for given p, x, a.

bounds_x0(p0, e0, a)

Get bounds on x for given p, e, a.

bounds_a(p0, e0, x0)

Get bounds on a for given p, e, x.

get_inspiral(m1, m2, a, y1, y2, y3, *args, Phi_phi0=0.0, Phi_theta0=0.0, Phi_r0=0.0, **kwargs)

Generate the inspiral.

This is the function for calling the creation of the trajectory. Inputs define the output time spacing.

Parameters:

  • m1 (float) – Mass of massive black hole in solar masses.

  • m2 (float) – Mass of compact object in solar masses.

  • a (float) – Dimensionless spin of massive black hole.

  • p0 – Initial semi-latus rectum in terms units of M (p/M).

  • e0 – Initial eccentricity (dimensionless).

  • x0 – Initial \(\cos{\iota}\). Note: This value is different from \(x_I\)

  • waveforms. (used in the relativistic)

  • *args – Added for flexibility.

  • Phi_phi0 (float) – Initial phase for \(\Phi_\phi\). Default is 0.0.

  • Phi_theta0 (float) – Initial phase for \(\Phi_\Theta\). Default is 0.0.

  • Phi_r0 (float) – Initial phase for \(\Phi_r\). Default is 0.0.

  • **kwargs – kwargs passed from parent.

  • y1 (float)

  • y2 (float)

  • y3 (float)

Returns:

Tuple of (t, p, e, x, Phi_phi, Phi_theta, Phi_r).

Return type:

tuple[ndarray]

__call__(*args, in_coordinate_time=True, dt=10.0, T=1.0, new_t=None, spline_kwargs=None, DENSE_STEPPING=False, buffer_length=1000, upsample=False, err=1e-11, fix_t=False, integrate_backwards=False, max_step_size=None, **kwargs)

Call function for trajectory interface.

This is the function for calling the creation of the trajectory. Inputs define the output time spacing.

Parameters:

  • *args – Input of variable number of arguments specific to the inspiral model (see the trajectory class’ get_inspiral method). Important Note: M must be the first parameter of any model that uses this base class.

  • in_coordinate_time (bool) – If True, the trajectory will be outputted in coordinate time. If False, the trajectory will be outputted in units of M. Default is True.

  • dt (float) – Time step for output waveform in seconds. Also sets initial step for integrator. Default is 10.0.

  • T (float) – Total observation time in years. Sets the maximum time for the integrator to run. Default is 1.0.

  • new_t (ndarray | None) – If given, this represents the final time array at which the trajectory is analyzed. This is performed by using a cubic spline on the integrator output. Default is None.

  • spline_kwargs (dict | None) – If using upsampling, spline_kwargs provides the kwargs desired for scipy.interpolate.CubicSpline. Default is {}.

  • DENSE_STEPPING (bool) – If True, the trajectory used in the integrator will be densely stepped at steps of dt. If False, the integrator will determine its stepping. Default is False.

  • buffer_length (int) – Sets the allocation of memory for trajectory parameters. This should be the maximum length expected for a trajectory. If it is reached, output arrays will be extended, but this is more expensive than allocating a larger array initially. Trajectories with default settings will be ~100 points. Default is 1000.

  • upsample (bool) – If True, upsample, with a cubic spline, the trajectories from 0 to T in steps of dt. Default is False.

  • err (float) – Tolerance for integrator. Default is 1e-10. Decreasing this parameter will give more steps over the trajectory, but if it is too small, memory issues will occur as the trajectory length will blow up. We recommend not adjusting this parameter.

  • fix_T – If upsampling, this will affect excess points in the t array. If True, it will shave any excess on the trajectories arrays where the time is greater than the overall time of the trajectory requested.

  • integrate_backwards (bool) – If True, the integrator will run backwards in time. Default is False.

  • max_step_size (float | None) – If given, this will set the maximum step size for the integrator (matching the convention according to in_coordinate_time). Default is None (no maximum step size).

  • **kwargs – kwargs passed to trajectory module. Default is {}.

  • fix_t (bool)

Returns:

Tuple of (t, p, e, Phi_phi, Phi_r, flux_norm).

Raises:

ValueError – If input parameters are not allowed in this model.

Return type:

tuple[ndarray]

classmethod all_citations()

Return all the citations from the registry as printable BibTeX string

Return type:

str

classmethod citation()

Return the module references as a printable BibTeX string.

Return type:

str

get_rhs_ode(m1, m2, a, y1, y2, y3, *args, Phi_phi0=0.0, Phi_theta0=0.0, Phi_r0=0.0, **kwargs)

Compute the right hand side of the ordinary differential equation.

This is a convenience function for interfacing with the call method of the ODE class.

Parameters:

  • m1 (float) – Mass of massive black hole in solar masses.

  • m2 (float) – Mass of compact object in solar masses.

  • a (float) – Dimensionless spin of massive black hole.

  • p0 – Initial semi-latus rectum in terms units of M (p/M).

  • e0 – Initial eccentricity (dimensionless).

  • x0 – Initial \(\cos{\iota}\). Note: This value is different from \(x_I\)

  • waveforms. (used in the relativistic)

  • *args – Added for flexibility.

  • Phi_phi0 (float) – Initial phase for \(\Phi_\phi\). Default is 0.0.

  • Phi_theta0 (float) – Initial phase for \(\Phi_\Theta\). Default is 0.0.

  • Phi_r0 (float) – Initial phase for \(\Phi_r\). Default is 0.0.

  • **kwargs – kwargs passed from parent.

  • y1 (float)

  • y2 (float)

  • y3 (float)

Returns:

Tuple of (t, p, e, x, Phi_phi, Phi_theta, Phi_r).

Return type:

tuple[ndarray]

few.trajectory.inspiral.get_0PA_frequencies(m1, m2, a, p, e, x, use_gpu=False)

Get frequencies for 0PA phase evolution in Hertz

Parameters:

  • m1 (float) – Mass of the massive black hole in solar masses.

  • m2 (float) – Mass of the secondary body in solar masses.

  • a (float | ndarray) – Dimensionless spin of massive black hole. If other parameters are arrays and the spin is scalar, it will be cast to a 1D array.

  • p (float | ndarray) – Values of separation, \(p\).

  • e (float | ndarray) – Values of eccentricity, \(e\).

  • x (float | ndarray) – Values of cosine of the inclination, \(x=\cos{I}\). Please note this is different from \(Y=\cos{\iota}\).

  • use_gpu (bool)

Returns:

Tuple of (OmegaPhi, OmegaTheta, OmegaR). These are 1D arrays or scalar values depending on inputs.

Return type:

tuple[float | ndarray]

ODE Integrators

class few.trajectory.integrate.Integrate(func, integrate_constants_of_motion=False, downsample=None, buffer_length=10000, enforce_schwarz_sep=False, max_iter=None, **kwargs)

Bases: ABC

Custom integrator class.

Flexible integration options. # TODO: discuss options

Parameters:

  • func (Type[ODEBase]) – ODE base function to integrate.

  • integrate_constants_of_motion (bool) – If True, integrate in ELQ rather than orbital parameters.

  • buffer_length (int) – Initial buffer length for trajectory. Will adjust itself if needed.

  • enforce_schwarz_sep (bool) – If True, enforce the Schwarzschild separatrix as the separatrix. \(p_s = 6 + 2e\).

  • max_iter (Optional[int]) – Maximum number of iterations for the integrator.

  • **kwargs – Keyword arguments for the ODE function, func.

property nparams: int

Dimensionality of ODE

property num_add_args: int

Number of additional arguments to the ODE

property convert_Y

Convert to Y if needed.

property background

Spacetime background (spin or no spin)

property equatorial

Orbit limited to equatorial plane.

property circular

Circular orbit.

property integrate_constants_of_motion

Integrating in ELQ.

property separatrix_buffer_dist

Buffer distance from separatrix.

distance_to_outer_boundary(y)

Distance to outer boundary of interpolation grid.

take_step(t, h, y)

Take a step of the integrator.

Parameters:

  • t (float) – Time.

  • h (float) – Step size.

  • y (ndarray) – Current position of integrator.

Returns:

(new time, new step size, new position)

Return type:

Tuple[float, float, ndarray]

log_failed_step(lastError)

Add to and check failed step count.

Raises:

ValueError – Too many failed tries.

Parameters:

lastError (Exception)

tune_initial_step_size(t0, y0, h_max)

Tune the initial step size of the integrator if using adaptive stepping.

Adapted from scipy.integrate._ivp.common.select_initial_step

Parameters:

  • t0 (float) – Initial time.

  • y0 (ndarray) – Initial state.

  • h_max (float) – Largest step size permitted.

Returns:

Initial step size.

Return type:

float

integrate(t0, y0)

Integrate from (t0, y0).

Parameters:

  • t0 (float) – Initial time.

  • y0 (ndarray) – Initial position.

Returns:

(Time array, Position array).

Return type:

Tuple[ndarray, ndarray]

action_function(t, y)

Act on the integrator.

Parameters:

  • t (float) – Time.

  • y (ndarray) – Current position of integrator.

Return type:

None

initialize_integrator(err=1e-11, DENSE_STEPPING=False, integrate_backwards=False, max_step_size=None, **kwargs)

Setup the integrator.

Parameters:

  • err (float) – Absolute tolerance of the integrator.

  • DENSE_STEPPING (bool) – If True, use fixed stepping.

  • integrate_backwards (bool) – If True, integrate backwards.

  • **kwargs – For future interoperability.

  • max_step_size (float | None)

property tolerance: float

Absolute tolerance of the integrator.

property dense_stepping: bool

If True, trajectory is using fixed stepping.

property npoints: int

Number of points in the trajectory.

property trajectory: ndarray

Trajectory array.

property integrator_t_cache: ndarray

Spline coefficient time array.

property integrator_spline_coeff: ndarray

Spline coefficients from dopr backend.

save_point(t, y, spline_output=None)

Save point in trajectory array.

Parameters:

  • t (float) – Time.

  • y (ndarray) – Current position of integrator.

  • spline_output (ndarray | None) – Spline coefficients from the backend.

Return type:

None

eval_integrator_spline(t_new)

Evaluate integration at new time array.

Parameters:

t_new (ndarray) – New time array.

Returns:

New trajectory.

Return type:

ndarray

eval_integrator_derivative_spline(t_new, order=1)

Evaluate integration derivatives at new time array.

Parameters:

  • t_new (ndarray) – New time array.

  • order (int) – Order of derivative to evaluate. Defaults to 1.

Returns:

New trajectory derivatives.

run_inspiral(m1, m2, a, y0, additional_args, T=1.0, dt=10.0, **kwargs)

Run inspiral integration.

Parameters:

  • m1 (float) – Larger mass in solar masses.

  • m2 (float) – Small mass in solar masses.

  • a (float) – Dimensionless spin of central black hole.

  • y0 (ndarray) – Initial set of coordinates for integration.

  • additional_args (list | ndarray)

  • T (float)

  • dt (float)

Returns:

Trajectory array for integrated coefficients.

Return type:

ndarray

stop_integrate_check(t, y)
Stop the inspiral when close to the separatrix (forwards integration)

or when close to the outer grid boundary (backwards integration).

Parameters:

  • t (float) – Time.

  • y (ndarray) – Current position of integrator.

Returns:

True if integration should be stopped. False otherwise.

Return type:

bool

inner_func_forward(t_step)

Evaluates the distance from the inner boundary at t=t_step. Also caches the state of the system as self._y_inner_cache.

inner_func_backward(t_step)

Evaluates the distance from the outer boundary at t=t_step. Also caches the state of the system as self._y_inner_cache.

finishing_function(t, y)

This function is called when the integrator has finished due to reaching a stopping condition. The function identifies the stopping condition and places the final point at the stopping boundary accordingly.

Parameters:
abstractmethod get_pex(y)

Handles integrator-specific conversion from y to pex and returns the result.

Parameters:

y (ndarray)

Return type:

ndarray

class few.trajectory.integrate.APEXIntegrate(func, integrate_constants_of_motion=False, downsample=None, buffer_length=10000, enforce_schwarz_sep=False, max_iter=None, **kwargs)

Bases: Integrate

An subclass of Integrate for integrating in pex coordinates.

Parameters:

  • func (Type[ODEBase])

  • integrate_constants_of_motion (bool)

  • buffer_length (int)

  • enforce_schwarz_sep (bool)

  • max_iter (Optional[int])

get_pex(y)

Handles integrator-specific conversion from y to pex and returns the result.

action_function(t, y)

Act on the integrator.

Parameters:

  • t (float) – Time.

  • y (ndarray) – Current position of integrator.

Return type:

None

property background

Spacetime background (spin or no spin)

property circular

Circular orbit.

property convert_Y

Convert to Y if needed.

property dense_stepping: bool

If True, trajectory is using fixed stepping.

distance_to_outer_boundary(y)

Distance to outer boundary of interpolation grid.

property equatorial

Orbit limited to equatorial plane.

eval_integrator_derivative_spline(t_new, order=1)

Evaluate integration derivatives at new time array.

Parameters:

  • t_new (ndarray) – New time array.

  • order (int) – Order of derivative to evaluate. Defaults to 1.

Returns:

New trajectory derivatives.

eval_integrator_spline(t_new)

Evaluate integration at new time array.

Parameters:

t_new (ndarray) – New time array.

Returns:

New trajectory.

Return type:

ndarray

finishing_function(t, y)

This function is called when the integrator has finished due to reaching a stopping condition. The function identifies the stopping condition and places the final point at the stopping boundary accordingly.

Parameters:
initialize_integrator(err=1e-11, DENSE_STEPPING=False, integrate_backwards=False, max_step_size=None, **kwargs)

Setup the integrator.

Parameters:

  • err (float) – Absolute tolerance of the integrator.

  • DENSE_STEPPING (bool) – If True, use fixed stepping.

  • integrate_backwards (bool) – If True, integrate backwards.

  • **kwargs – For future interoperability.

  • max_step_size (float | None)

inner_func_backward(t_step)

Evaluates the distance from the outer boundary at t=t_step. Also caches the state of the system as self._y_inner_cache.

inner_func_forward(t_step)

Evaluates the distance from the inner boundary at t=t_step. Also caches the state of the system as self._y_inner_cache.

integrate(t0, y0)

Integrate from (t0, y0).

Parameters:

  • t0 (float) – Initial time.

  • y0 (ndarray) – Initial position.

Returns:

(Time array, Position array).

Return type:

Tuple[ndarray, ndarray]

property integrate_constants_of_motion

Integrating in ELQ.

property integrator_spline_coeff: ndarray

Spline coefficients from dopr backend.

property integrator_t_cache: ndarray

Spline coefficient time array.

log_failed_step(lastError)

Add to and check failed step count.

Raises:

ValueError – Too many failed tries.

Parameters:

lastError (Exception)

property nparams: int

Dimensionality of ODE

property npoints: int

Number of points in the trajectory.

property num_add_args: int

Number of additional arguments to the ODE

run_inspiral(m1, m2, a, y0, additional_args, T=1.0, dt=10.0, **kwargs)

Run inspiral integration.

Parameters:

  • m1 (float) – Larger mass in solar masses.

  • m2 (float) – Small mass in solar masses.

  • a (float) – Dimensionless spin of central black hole.

  • y0 (ndarray) – Initial set of coordinates for integration.

  • additional_args (list | ndarray)

  • T (float)

  • dt (float)

Returns:

Trajectory array for integrated coefficients.

Return type:

ndarray

save_point(t, y, spline_output=None)

Save point in trajectory array.

Parameters:

  • t (float) – Time.

  • y (ndarray) – Current position of integrator.

  • spline_output (ndarray | None) – Spline coefficients from the backend.

Return type:

None

property separatrix_buffer_dist

Buffer distance from separatrix.

stop_integrate_check(t, y)
Stop the inspiral when close to the separatrix (forwards integration)

or when close to the outer grid boundary (backwards integration).

Parameters:

  • t (float) – Time.

  • y (ndarray) – Current position of integrator.

Returns:

True if integration should be stopped. False otherwise.

Return type:

bool

take_step(t, h, y)

Take a step of the integrator.

Parameters:

  • t (float) – Time.

  • h (float) – Step size.

  • y (ndarray) – Current position of integrator.

Returns:

(new time, new step size, new position)

Return type:

Tuple[float, float, ndarray]

property tolerance: float

Absolute tolerance of the integrator.

property trajectory: ndarray

Trajectory array.

tune_initial_step_size(t0, y0, h_max)

Tune the initial step size of the integrator if using adaptive stepping.

Adapted from scipy.integrate._ivp.common.select_initial_step

Parameters:

  • t0 (float) – Initial time.

  • y0 (ndarray) – Initial state.

  • h_max (float) – Largest step size permitted.

Returns:

Initial step size.

Return type:

float

class few.trajectory.integrate.AELQIntegrate(func, integrate_constants_of_motion=False, downsample=None, buffer_length=10000, enforce_schwarz_sep=False, max_iter=None, **kwargs)

Bases: Integrate

An subclass of Integrate for integrating in ELQ coordinates.

Parameters:

  • func (Type[ODEBase])

  • integrate_constants_of_motion (bool)

  • buffer_length (int)

  • enforce_schwarz_sep (bool)

  • max_iter (Optional[int])

action_function(t, y)

Act on the integrator.

Parameters:

  • t (float) – Time.

  • y (ndarray) – Current position of integrator.

Return type:

None

property background

Spacetime background (spin or no spin)

property circular

Circular orbit.

property convert_Y

Convert to Y if needed.

property dense_stepping: bool

If True, trajectory is using fixed stepping.

distance_to_outer_boundary(y)

Distance to outer boundary of interpolation grid.

property equatorial

Orbit limited to equatorial plane.

eval_integrator_derivative_spline(t_new, order=1)

Evaluate integration derivatives at new time array.

Parameters:

  • t_new (ndarray) – New time array.

  • order (int) – Order of derivative to evaluate. Defaults to 1.

Returns:

New trajectory derivatives.

eval_integrator_spline(t_new)

Evaluate integration at new time array.

Parameters:

t_new (ndarray) – New time array.

Returns:

New trajectory.

Return type:

ndarray

finishing_function(t, y)

This function is called when the integrator has finished due to reaching a stopping condition. The function identifies the stopping condition and places the final point at the stopping boundary accordingly.

Parameters:
get_pex(y)

Handles integrator-specific conversion from y to pex and returns the result.

initialize_integrator(err=1e-11, DENSE_STEPPING=False, integrate_backwards=False, max_step_size=None, **kwargs)

Setup the integrator.

Parameters:

  • err (float) – Absolute tolerance of the integrator.

  • DENSE_STEPPING (bool) – If True, use fixed stepping.

  • integrate_backwards (bool) – If True, integrate backwards.

  • **kwargs – For future interoperability.

  • max_step_size (float | None)

inner_func_backward(t_step)

Evaluates the distance from the outer boundary at t=t_step. Also caches the state of the system as self._y_inner_cache.

inner_func_forward(t_step)

Evaluates the distance from the inner boundary at t=t_step. Also caches the state of the system as self._y_inner_cache.

integrate(t0, y0)

Integrate from (t0, y0).

Parameters:

  • t0 (float) – Initial time.

  • y0 (ndarray) – Initial position.

Returns:

(Time array, Position array).

Return type:

Tuple[ndarray, ndarray]

property integrate_constants_of_motion

Integrating in ELQ.

property integrator_spline_coeff: ndarray

Spline coefficients from dopr backend.

property integrator_t_cache: ndarray

Spline coefficient time array.

log_failed_step(lastError)

Add to and check failed step count.

Raises:

ValueError – Too many failed tries.

Parameters:

lastError (Exception)

property nparams: int

Dimensionality of ODE

property npoints: int

Number of points in the trajectory.

property num_add_args: int

Number of additional arguments to the ODE

run_inspiral(m1, m2, a, y0, additional_args, T=1.0, dt=10.0, **kwargs)

Run inspiral integration.

Parameters:

  • m1 (float) – Larger mass in solar masses.

  • m2 (float) – Small mass in solar masses.

  • a (float) – Dimensionless spin of central black hole.

  • y0 (ndarray) – Initial set of coordinates for integration.

  • additional_args (list | ndarray)

  • T (float)

  • dt (float)

Returns:

Trajectory array for integrated coefficients.

Return type:

ndarray

save_point(t, y, spline_output=None)

Save point in trajectory array.

Parameters:

  • t (float) – Time.

  • y (ndarray) – Current position of integrator.

  • spline_output (ndarray | None) – Spline coefficients from the backend.

Return type:

None

property separatrix_buffer_dist

Buffer distance from separatrix.

stop_integrate_check(t, y)
Stop the inspiral when close to the separatrix (forwards integration)

or when close to the outer grid boundary (backwards integration).

Parameters:

  • t (float) – Time.

  • y (ndarray) – Current position of integrator.

Returns:

True if integration should be stopped. False otherwise.

Return type:

bool

take_step(t, h, y)

Take a step of the integrator.

Parameters:

  • t (float) – Time.

  • h (float) – Step size.

  • y (ndarray) – Current position of integrator.

Returns:

(new time, new step size, new position)

Return type:

Tuple[float, float, ndarray]

property tolerance: float

Absolute tolerance of the integrator.

property trajectory: ndarray

Trajectory array.

tune_initial_step_size(t0, y0, h_max)

Tune the initial step size of the integrator if using adaptive stepping.

Adapted from scipy.integrate._ivp.common.select_initial_step

Parameters:

  • t0 (float) – Initial time.

  • y0 (ndarray) – Initial state.

  • h_max (float) – Largest step size permitted.

Returns:

Initial step size.

Return type:

float

ODE Classes

Contains the ODEBase baseclass that handles evaluating the ODE

class few.trajectory.ode.base.ODEBase(*args, use_ELQ=False, integrate_backwards=False, **kwargs)

Bases: Citable

A baseclass for handling the evaluation of ODE derivatives in the trajectory module.

To define a new trajectory function, subclass this function and define evaluate_rhs. Make sure to update the relevant class attributes as well. See the documentation for examples on how to do this.

use_ELQ

If True, the ODE will take as input (and output derivatives of) the integrals of motion (E, L, Q). Defaults to False.

Type:

bool

num_add_args

Number of additional arguments being passed to the ODE function.

Type:

int

integrate_backwards

If True, the ODE corresponds to integrating backwards in time. Defaults to False.

Type:

bool

property convert_Y

If True, the inclination coordinate is assumed to be Y and is converted accordingly. Defaults to False.

property equatorial

If True, the inclination coordinate is assumed to be +/- 1. Defaults to False.

property circular

If True, the eccentricity coordinate is assumed to be 0. Defaults to False.

property supports_ELQ

If True, this ODE can take as input (and output derivatives of) the integrals of motion (E, L, Q) if initialised with use_ELQ=True. Defaults to False.

property background

A string describing the background spacetime. Either “Kerr” or “Schwarzschild”. Defaults to “Kerr”.

property separatrix_buffer_dist

A float describing the value of “p” at which the trajectory should terminate at, with respect to the separatrix. A value of 0 would mean that the trajectory terminates at the separatrix. Defaults to 0.05

property nparams

An integer describing the number of parameters this ODE will integrate. Defaults to 6 (three orbital elements, three orbital phases).

property required_add_args

A list of names of the additional arguments the user must supply for this model (beyond the orbital elements and phases). If fewer than these are provided, the trajectory raises a clear ValueError naming the missing argument(s), rather than silently defaulting them to zero. Defaults to an empty list (no required additional arguments).

property flux_output_convention

A string describing the coordinate convention of the fluxes for this model, as output by evaluate_rhs. These are either “pex” or “ELQ”. If “ELQ”, a Jacobian transformation is performed if the output derivatives of the model are expected to be in the “pex” convention. Defaults to “pex”.

For models that do not perform interpolation to generate fluxes (e.g., PN5), this property is not accessed.

evaluate_rhs(y, **kwargs)

This function evaluates the right-hand side of the ODE at the point y. An ODE model can be defined by subclassing the ODEBase class and implementing this method.

Parameters:

y (ndarray)

Return type:

NotImplementedError

modify_rhs_before_Jacobian(ydot, y, **kwargs)

This function allows the user to modify the right-hand side of the ODE before the Jacobian transform is applied. This is particularly useful if the user wishes to apply modifications to fluxes of the integrals of motion “ELQ”, but integrate a trajectory in terms of “pex”.

By default, this function returns the input right-hand side unchanged.

Parameters:

  • ydot (ndarray)

  • y (ndarray)

Return type:

ndarray

min_p(e, x=1, a=0, separatrix_buffer=None)

Computes the minimum value of the radial coordinate p for a given eccentricity and inclination for this model. Trajectory models implementing their own interpolants should override this function to return the minimum value corresponding to the precomputed grid boundaries.

By default, this function assumes things are rectilinear and returns p_sep + self.separatrix_buffer_dist.

Parameters:
Return type:

float | ndarray

max_p(e, x=1.0, a=0.0)

Computes the maximum value of the semilatus rectum p for a given eccentricity and inclination for this model. Trajectory models implementing their own interpolants should override this function to return the maximum value corresponding to the precomputed grid boundaries.

By default, this function returns np.inf (assumes no bound on p).

Parameters:
Return type:

float | ndarray

min_e(p, x=1.0, a=0.0)

Computes the minimum value of the eccentricity e for a given semilatus rectum and inclination for this model. Trajectory models implementing their own interpolants should override this function to return the minimum value corresponding to the precomputed grid boundaries.

By default, this function assumes minimal eccentricity corresponds to circular orbits and returns 0.

Parameters:
Return type:

float | ndarray

max_e(p, x=1, a=0)

Computes the maximum value of the eccentricity e for a given semilatus rectum and inclination for this model. Trajectory models implementing their own interpolants should override this function to return the minimum value corresponding to the precomputed grid boundaries.

By default, this function assumes no orbital bounds on eccentricity and returns np.inf.

Parameters:
Return type:

float | ndarray

modify_rhs(ydot, y, **kwargs)

This function allows the user to modify the right-hand side of the ODE after any required Jacobian transforms have been applied.

By default, this function returns the input right-hand side unchanged.

Parameters:

  • ydot (ndarray)

  • y (ndarray)

Return type:

ndarray

interpolate_flux_grids(*args)

This function handles the interpolation of the fluxes from the precomputed grids, including parameter transformations. Each stock model implements this function to handle the specifics of their own interpolants. To easily incorporate interpolated fluxes into their own models, users can subclass the stock models; the interpolated fluxes can then be accessed from evaluate_rhs by calling this function.

This method should also handle checking that the input parameters are within the bounds of the precomputed grids. Failure to do so may result in erroneous behaviour during trajectory evaluation in which fluxes are extrapolated.

Return type:

NotImplementedError

distance_to_outer_boundary(y)

This function returns the distance to the outer boundary of the interpolation grid. This is necessary for backwards integration, which performs root-finding to ensure that the trajectory ends on this outer boundary. Root-finding is initiated when this function returns a negative value.

Each stock model implements this function to handle the specifics of their own interpolants. For models that do not use interpolation, this function returns a positive constant (such that root-finding never occurs) and does not need to be implemented.

Parameters:

y (ndarray)

Return type:

float

get_pex(y)

This function converts the integrals of motion (E, L, Q) to the orbital elements (p, e, x), if required.

Parameters:

y (ndarray)

Return type:

tuple[float]

classmethod all_citations()

Return all the citations from the registry as printable BibTeX string

Return type:

str

cache_values_and_check_bounds(y)

This function checks the input points to ensure they are within the physical bounds of the inspiral parameter space. These checks include ensuring that the separatrix has not been crossed, and that the eccentricity is within bounds.

Returns a boolean indicating whether the input was in bounds.

Parameters:

y (ndarray)

Return type:

bool

classmethod citation()

Return the module references as a printable BibTeX string.

Return type:

str

classmethod module_references()

Method implemented by each class to define its list of references

Return type:

Iterable[REFERENCE | str]

__call__(y, out=None, **kwargs)

Call self as a function.

Parameters:

  • y (list | ndarray)

  • out (ndarray | None)

  • kwargs (dict | None)

Return type:

ndarray

class few.trajectory.ode.flux.SchwarzEccFlux(*args, use_ELQ=False, **kwargs)

Bases: ODEBase

Schwarzschild eccentric flux ODE.

Parameters:

use_ELQ (bool) – If True, the ODE will output derivatives of the orbital elements of (E, L, Q). Defaults to False.

property supports_ELQ

If True, this ODE can take as input (and output derivatives of) the integrals of motion (E, L, Q) if initialised with use_ELQ=True. Defaults to False.

min_p(e, x=1, a=0, separatrix_buffer=None)

Computes the minimum value of the radial coordinate p for a given eccentricity and inclination for this model. Trajectory models implementing their own interpolants should override this function to return the minimum value corresponding to the precomputed grid boundaries.

By default, this function assumes things are rectilinear and returns p_sep + self.separatrix_buffer_dist.

max_p(e, x=1, a=0)

Computes the maximum value of the semilatus rectum p for a given eccentricity and inclination for this model. Trajectory models implementing their own interpolants should override this function to return the maximum value corresponding to the precomputed grid boundaries.

By default, this function returns np.inf (assumes no bound on p).

max_e(p, x=1, a=0)

Computes the maximum value of the eccentricity e for a given semilatus rectum and inclination for this model. Trajectory models implementing their own interpolants should override this function to return the minimum value corresponding to the precomputed grid boundaries.

By default, this function assumes no orbital bounds on eccentricity and returns np.inf.

distance_to_outer_boundary(y)

This function returns the distance to the outer boundary of the interpolation grid. This is necessary for backwards integration, which performs root-finding to ensure that the trajectory ends on this outer boundary. Root-finding is initiated when this function returns a negative value.

Each stock model implements this function to handle the specifics of their own interpolants. For models that do not use interpolation, this function returns a positive constant (such that root-finding never occurs) and does not need to be implemented.

interpolate_flux_grids(p, e, Omega_phi, pLSO=None)

This function handles the interpolation of the fluxes from the precomputed grids, including parameter transformations. Each stock model implements this function to handle the specifics of their own interpolants. To easily incorporate interpolated fluxes into their own models, users can subclass the stock models; the interpolated fluxes can then be accessed from evaluate_rhs by calling this function.

This method should also handle checking that the input parameters are within the bounds of the precomputed grids. Failure to do so may result in erroneous behaviour during trajectory evaluation in which fluxes are extrapolated.

Parameters:
Return type:

tuple[float]

evaluate_rhs(y)

This function evaluates the right-hand side of the ODE at the point y. An ODE model can be defined by subclassing the ODEBase class and implementing this method.

Parameters:

y (list[float] | ndarray)

Return type:

list[float | ndarray]

class few.trajectory.ode.flux.KerrEccEqFlux(*args, use_ELQ=False, downsample=None, flux_output_convention='pex', **kwargs)

Bases: ODEBase

Kerr eccentric equatorial flux ODE.

Parameters:

  • use_ELQ (bool) – If True, the ODE will output derivatives of the orbital elements of (E, L, Q). Defaults to False.

  • downsample – List of two 3-tuples of integers to downsample the flux grid in u, w, z. The first list element refers to the inner grid, the second to the outer. Useful for testing error convergence. Defaults to None (no downsampling).

  • flux_output_convention – String specifying the output convention for the fluxes. Can be either “pex” (default) or “ELQ”.

property separatrix_buffer_dist_grid

The distance from the separatrix to the minimum p grid value.

Type:

float

property supports_ELQ

Always True for Kerr eccentric equatorial flux.

Type:

bool

isvalid_x(x)

Check if x is valid. For Kerr equatorial flux, x must be either 1 or -1.

Parameters:

x (float or np.ndarray) – Cosine of the inclination angle.

Raises:

ValueError – If any value of x is not 1 or -1.

isvalid_e(e, e_buffer=[0, 0])

Check if e is valid. For Kerr equatorial flux, e must be between emin and emax, where emin is the lower bound and emax is the upper bound adjusted by the buffer.

Parameters:

  • e (float or np.ndarray) – Eccentricity.

  • e_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid e range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

Raises:

ValueError – If any value of e is out of bounds.

isvalid_p(p, p_buffer=[0, 0])

Check if p is valid. For Kerr equatorial flux, p must be between pmin and pmax, where pmin is the lower bound adjusted by the separatrix buffer and pmax is the upper bound adjusted by the buffer.

Parameters:

  • p (float or np.ndarray) – Semi-latus rectum.

  • p_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid p range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

Raises:

ValueError – If any value of p is out of bounds.

isvalid_a(a, a_buffer=[0, 0])

Check if a is valid. For Kerr equatorial flux, a must be between amin and amax, where amin is the lower bound adjusted by the buffer and amax is the upper bound adjusted by the buffer.

Parameters:

  • a (float or np.ndarray) – Primary spin parameter.

  • a_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid a range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

Raises:

ValueError – If any value of a is out of bounds.

min_p(e=0, x=1, a=0, separatrix_buffer=None, max_e_buffer=0)

Return the minimum value of p for given e, x, and a, considering the separatrix buffer and maximum e buffer.

Parameters:

  • e (float) – Eccentricity.

  • x (float) – Inclination parameter.

  • a (float) – Primary spin parameter.

  • separatrix_buffer (float, optional) – Buffer distance from the separatrix. Defaults to None, which uses the class’s separatrix_buffer_dist_grid.

  • max_e_buffer (float, optional) – Buffer for the maximum e value. Defaults to 0.

max_p(e=0, x=1, a=0)

Return the maximum value of p for given e, x, and a.

Parameters:

  • e (float) – Eccentricity.

  • x (float) – Inclination parameter.

  • a (float) – Primary spin parameter.

min_e(p=20, x=1, a=0)

Return the minimum value of e for given p, x, and a.

Parameters:

  • p (float) – Semi-latus rectum.

  • x (float) – Inclination parameter.

  • a (float) – Primary spin parameter.

max_e(p=20, x=1, a=0, separatrix_buffer=None, max_e_buffer=0)

Return the maximum value of e for given p, x, and a.

Parameters:

  • p (float) – Semi-latus rectum.

  • x (float) – Inclination parameter.

  • a (float) – Primary spin parameter.

min_a(p=20, e=0, x=1)

Return the minimum value of a for given p, e, and x.

Parameters:

  • p (float) – Semi-latus rectum.

  • e (float) – Eccentricity.

  • x (float) – Inclination parameter.

max_a(p=20, e=0, x=1)

Return the maximum value of a for given p, e, and x.

Parameters:

  • p (float) – Semi-latus rectum.

  • e (float) – Eccentricity.

  • x (float) – Inclination parameter.

bounds_a(p=20, e=0, x=1, a_buffer=[0, 0])

Return the minimum and maximum bounds of a for given p, e, and x, considering the buffer.

Parameters:

  • p (float) – Semi-latus rectum.

  • e (float) – Eccentricity.

  • x (float) – Inclination parameter.

  • a_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid a range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

bounds_p(e=0, x=1, a=0, p_buffer=[0, 0], separatrix_buffer=None, max_e_buffer=0)

Return the minimum and maximum bounds of p for given e, x, and a, considering the buffer.

Parameters:

  • e (float) – Eccentricity.

  • x (float) – Inclination parameter.

  • a (float) – Primary spin parameter.

  • p_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid p range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

  • separatrix_buffer (float, optional) – Buffer distance from the separatrix. Defaults to None, which uses the class’s separatrix_buffer_dist_grid.

  • max_e_buffer (float, optional) – Buffer for the maximum e value. Defaults to 0.

bounds_e(p=20, x=1, a=0, e_buffer=[0, 0], separatrix_buffer=None, max_e_buffer=0)

Return the minimum and maximum bounds of e for given p, x, and a, considering the buffer.

Parameters:

  • p (float) – Semi-latus rectum.

  • x (float) – Inclination parameter.

  • a (float) – Primary spin parameter.

  • e_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid e range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

  • separatrix_buffer (float, optional) – Buffer distance from the separatrix. Defaults to None, which uses the class’s separatrix_buffer_dist_grid.

  • max_e_buffer (float, optional) – Buffer for the maximum e value. Defaults to 0.

isvalid_pex(p=20, e=0, x=1, a=0, p_buffer=[0, 0], e_buffer=[0, 0], a_buffer=[0, 0], separatrix_buffer=None, max_e_buffer=0)

Raise an error if the given (p, e, x, a) values are not valid for interpolation, considering the specified buffers.

Parameters:

  • p (float) – Semi-latus rectum.

  • e (float) – Eccentricity.

  • x (float) – Inclination parameter.

  • a (float) – Primary spin parameter.

  • p_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid p range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

  • e_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid e range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

  • a_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid a range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

  • separatrix_buffer (float, optional) – Buffer distance from the separatrix. Defaults to None, which uses the class’s separatrix_buffer_dist_grid.

  • max_e_buffer (float, optional) – Buffer for the maximum e value. Defaults to 0.

Raises:

AssertionError – If p is outside the valid range given e, x, and a.

distance_to_outer_boundary(y)

Return the distance to the outer boundary of the valid (p, e) region for interpolation, given the state vector y.

interpolate_flux_grids(p, e, x=1, a=0, pLSO=None)

Return the interpolated fluxes (Edot, Ldot) or (pdot, edot) for given (p, e, x, a), using the precomputed flux grids.

Parameters:

  • p (float) – Semi-latus rectum.

  • e (float) – Eccentricity.

  • x (float) – Inclination parameter. Must be either 1 or -1. Defaults to 1.

  • a (float) – Primary spin parameter. Defaults to 0.

  • pLSO (float, optional) – The value of p at the last stable orbit (LSO) for the given (a, e, x). If not provided, it will be computed using get_separatrix(a, e, x).

Return type:

tuple[float]

evaluate_rhs(y)

Returns the right-hand side of the ODE for the Kerr eccentric equatorial flux, given the state vector y. The state vector can be either (p, e, x, Phi_phi, Phi_theta, Phi_r) or (E, L, Q, Phi_phi, Phi_theta, Phi_r) depending on the use_ELQ flag.

Parameters:

y (list[float] or np.ndarray) – The state vector, either (p, e, x, Phi_phi, Phi_theta, Phi_r) or (E, L, Q, Phi_phi, Phi_theta, Phi_r).

Return type:

list[float | ndarray]

class few.trajectory.ode.circ1pat1r.TrajectoryCirc1PAT1R(*args, downsample=None, evolve_primary=True, **kwargs)

Bases: ODEBase

Trajectory of the first-post-adiabatic, slowly-spinning circular inspiral model 1PAT1R of Ref. 2510.16113.

Circular equatorial inspiral of an arbitrary spinning secondary (chit2) BH around a slowly-spinning primary (chit1), at first post-adiabatic order. The model also includes 1PA corrections due to the evolution of the primary’s mass and spin.

The state vector has 8 components:

y = [p, e, xI, Phi_phi, Phi_theta, Phi_r, deltaM, delta_chit1]

where e = 0 and xI = 1 (circular, equatorial motion). The extra slots deltaM and delta_chit1 are the deviations of the total mass and reduced primary spin from their initial values.

The flux interpolants are read from the HDF5 data file Trajectory1PAT1R.h5.

Parameters:

  • *args

  • downsample (dict, optional) – Dictionary of integer downsampling factors for each flux grid, keyed by grid name. Useful for testing convergence. Defaults to None (no downsampling). Get the default dictionary of downsampling factors with class attribute default_downsample.

  • evolve_primary (bool, optional) – Whether to evolve the primary’s mass and spin. Defaults to True. If False, the primary’s mass and spin are held fixed at their initial values.

  • **kwargs

Raises:

  • ValueError – If any downsampling value is not an integer.

  • ValueError – If a downsampling key is not in the set of grid names.

  • ValueError – If the loaded flux data contains NaN values.

classmethod module_references()

Return citations related to this module

Return type:

list[REFERENCE]

property circular

Always True; this model is restricted to circular orbits (e = 0).

Type:

bool

property background

Spacetime background for this trajectory.

Returns ‘Kerr’ so that the primary spin (chi1) passed as ‘a’ is not zeroed out by the inspiral integrator.

Type:

str

property nparams

An integer describing the number of parameters this ODE will integrate.

property required_add_args

The secondary spin chi2 must always be supplied as an additional argument.

isvalid_x(x, **kwargs)

Raise ValueError if the cosine of the inclination angle is outside the valid range.

Parameters:

x (float) – Cosine of the inclination angle.

min_e(p, x=1.0, a=0.0, **kwargs)

Computes the minimum value of the eccentricity e for a given semilatus rectum and inclination for this model. Trajectory models implementing their own interpolants should override this function to return the minimum value corresponding to the precomputed grid boundaries.

By default, this function assumes minimal eccentricity corresponds to circular orbits and returns 0.

Parameters:

  • p (float or np.ndarray) – Semilatus rectum.

  • x (float or np.ndarray) – Cosine of the inclination angle.

  • a (float or np.ndarray) – Primary spin.

Return type:

float | ndarray

max_e(p, x=1, a=0, **kwargs)

Computes the maximum value of the eccentricity e for a given semilatus rectum and inclination for this model. Trajectory models implementing their own interpolants should override this function to return the minimum value corresponding to the precomputed grid boundaries.

By default, this function assumes no orbital bounds on eccentricity and returns np.inf.

Parameters:

  • p (float or np.ndarray) – Semilatus rectum.

  • x (float or np.ndarray) – Cosine of the inclination angle.

  • a (float or np.ndarray) – Primary spin.

Return type:

float | ndarray

isvalid_e(e, **kwargs)

Raise ValueError if the eccentricity is outside the valid range.

Parameters:

e (float) – Eccentricity.

Raises:

ValueError – If e is outside [min_e, max_e]

isvalid_a(a, **kwargs)

Raise ValueError if the primary spin is outside the valid range.

Parameters:

a (float) – Primary dimensionless spin.

Raises:

ValueError – If a is outside [_min_chi1, _max_chi1].

min_p(e=0.0, x=1.0, a=0.0, separatrix_buffer=None, **kwargs)

Return the minumum valid p for this model (data grid upper boundary).

Parameters:

  • e (float) – Eccentricity (ignored since this model is circular).

  • x (float) – Cosine of the inclination angle (ignored since this model is equatorial

  • a (float) – Primary spin

  • separatrix_buffer (float, optional) – Buffer distance to add to the minimum valid p. Defaults to None, which is treated as 0.0.

max_p(e=0.0, x=1.0, a=0.0, **kwargs)

Return the maximum valid p for this model (data grid upper boundary).

Parameters:

  • e (float) – Eccentricity (ignored since this model is circular).

  • x (float) – Cosine of the inclination angle (ignored since this model is equatorial

  • a (float) – Primary spin

bounds_a(a_buffer=[0, 0], **kwargs)

Return the valid a range for this model with optional buffer.

Parameters:

a_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid a range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

bounds_p(e=0.0, x=1.0, a=0.0, p_buffer=[0.0, 0.0], separatrix_buffer=None, **kwargs)

Return the valid p range for this model, accounting for the separatrix and optional buffer.

Parameters:

  • e (float) – Eccentricity (ignored since this model is circular).

  • x (float) – Cosine of the inclination angle (ignored since this model is equatorial).

  • a (float) – Primary spin.

  • p_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid p range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

  • separatrix_buffer (float, optional) – Buffer distance to add to the minimum valid p. Defaults to None, which is treated as 0.0.

isvalid_p(p, e=0.0, x=1.0, a=0.0, separatrix_buffer=None, **kwargs)

Raise ValueError if p is outside the flux-data grid.

Parameters:

  • p (float) – Dimensionless semi-latus rectum p=Omega_phi^(-2/3).

  • e (float) – Eccentricity (ignored since this model is circular).

  • x (float) – Cosine of the inclination angle (ignored since this model is equatorial).

  • a (float) – Primary spin.

  • separatrix_buffer (float, optional) – Buffer distance to add to the minimum valid p. Defaults to None, which is treated as 0.0.

isvalid_pex(p=20, e=0, x=1, a=0, p_buffer=[0, 0], separatrix_buffer=None, **kwargs)

Check whether p, e, x, and primary spin are within the valid ranges.

Parameters:

  • p (float) – Dimensionless semi-latus rectum p=Omega_phi^(-2/3).

  • e (float) – Eccentricity (ignored since this model is circular).

  • x (float) – Cosine of the inclination angle (ignored since this model is equatorial).

  • a (float) – Dimensionless primary spin.

  • p_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid p range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

  • separatrix_buffer (float, optional) – Buffer distance to add to the minimum valid p. Defaults to None, which is treated as 0.0.

Raises:

AssertionError – If p is outside the valid range given e, x, and a.

add_fixed_parameters(m1, m2, chi1, additional_args)

Store the symmetric mass ratio nu, the spin parameters, and the reduced spins chit1, chit2.

The additional_args array must contain, in order:

  1. chi2 — dimensionless spin of the secondary compact object.

Parameters:

  • m1 (float) – Mass of the primary in solar masses.

  • m2 (float) – Mass of the secondary in solar masses.

  • chi1 (float) – Dimensionless spin of the primary.

  • additional_args (array-like) – Extra arguments; additional_args[0] is the secondary spin chi2.

Raises:

ValueError – If nu, chi1, or chi2 fall outside their respective valid ranges.

balance_law_numerator(nu, p, chit1, chit2, deltaM)
Total 1PA numerator of the flux-balance law, Eq. (84) of Ref. 2510.16113.

plus an additional correction to the delta M term that is missing in Eq. (84).

Returns the numerator of the flux-balance law at 1PA order.

Parameters:

  • nu (float) – Symmetric mass ratio nu = m1*m2/(m1+m2)^2.

  • p (float) – Dimensionless semi-latus rectum p=Omega_phi^(-2/3).

  • chit1 (float) – Primary reduced spin.

  • chit2 (float) – Secondary reduced spin.

  • deltaM (float) – Deviation of the total mass from its initial value.

Returns:

Numerator of the flux-balance law.

Return type:

float

dEdOmega(nu, p, chit1, chit2, deltaM)

1PA derivative of the binding energy with respect to the orbital frequency, Eq. (108) of Ref. 2510.16113.

Parameters:

  • nu (float) – Symmetric mass ratio nu = m1*m2/(m1+m2)^2.

  • p (float) – Dimensionless semi-latus rectum p=Omega_phi^(-2/3).

  • chit1 (float) – Primary reduced spin.

  • chit2 (float) – Secondary reduced spin.

  • deltaM (float) – Deviation of the total mass from its initial value.

Returns:

dE/dOmega_phi at 1PA order.

Return type:

float

evaluate_rhs(y)

Evaluate the right-hand side of the 1PAT1R ODE.

Computes the time derivatives of all state-vector components at the current state y. The state vector layout is:

y = [p, e, xI, nu*Phi_phi, nu*Phi_theta, nu*Phi_r, deltaM, delta_chit1]

The evolution equations are derived from the 1PA flux-balance law, Eq. (84) of Ref. 2510.16113 [with an additional correction to the delta M term that is missing in Eq. (84)].

Parameters:

y (array-like) – Current state vector of length 8: [p, e, xI, nu*Phi_phi, nu*Phi_theta, nu*Phi_r, deltaM, delta_chit1].

Returns:

Derivatives

[pdot, 0, 0, Omega_phi, 0, 0, deltaMdot, chit1dot].

Return type:

list[float]

Raises:

ValueError – If p is outside the grid bounds.

__call__(y, out=None, **kwargs)

Evaluate the ODE and return the derivative array.

Overrides __call__() to pre-allocate an 8-element output buffer (instead of the base-class default of 6) and to skip the ELQ Jacobian transform, which is not used in this model.

Parameters:

  • y (array-like) – Current state vector of length 8.

  • out (np.ndarray, optional) – Pre-allocated output array of length 8. If None, a new array is allocated.

  • **kwargs

Returns:

Derivative array of length 8.

Return type:

np.ndarray

class few.trajectory.ode.pn5.PN5(*args, use_ELQ=False, integrate_backwards=False, **kwargs)

Bases: ODEBase

property supports_ELQ

Whether the ODE supports evolution of E, Lz, and Q.

Type:

bool

max_p(e, x, a)

Returns the maximum valid value of p for given e, x, and a. For this ODE, there is no upper limit on p, so it returns infinity.

Parameters:

  • e (float) – Eccentricity.

  • x (float) – Cosine of the inclination angle.

  • a (float) – Primary spin parameter.

min_p(e, x, a, separatrix_buffer=None)

Returns the minimum valid value of p for given e, x, and a, which is the separatrix value plus a buffer distance.

Parameters:

  • e (float) – Eccentricity.

  • x (float) – Cosine of the inclination angle.

  • a (float) – Primary spin parameter.

  • separatrix_buffer (float, optional) – The buffer distance from the separatrix. If not provided, it will be computed using self.separatrix_buffer_dist.

Returns:

The minimum valid value of p.

Return type:

float

isvalid_x(x, x_buffer=[0, 0])

Raises a ValueError if x is outside the valid range of [-1 + x_buffer[0], 1 - x_buffer[1]].

Parameters:

  • x (float) – Cosine of the inclination angle.

  • x_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid x range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

Raises:

ValueError – If x is outside the valid range.

isvalid_e(e, e_buffer=[0, 0])

Raises a ValueError if e is outside the valid range of [0 + e_buffer[0], 1 - e_buffer[1]].

Parameters:

  • e (float) – Eccentricity.

  • e_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid e range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

Raises:

ValueError – If e is outside the valid range.

isvalid_p(p, p_buffer=[0, 0], separatrix_buffer=None)

Raises a ValueError if p is less than the minimum valid value of p.

Parameters:

  • p (float) – Semi-latus rectum.

  • p_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid p range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

  • separatrix_buffer (float, optional) – The buffer distance from the separatrix. If not provided, it will be computed using self.separatrix_buffer_dist.

Raises:

ValueError – If p is less than the minimum valid value of p.

isvalid_a(a, a_buffer=[0, 0])

Raises a ValueError if a is outside the valid range of [-1 + a_buffer[0], 1 - a_buffer[1]].

Parameters:

  • a (float) – Primary spin parameter.

  • a_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid a range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

Raises:

ValueError – If a is outside the valid range.

bounds_p(e, x=1, a=0, p_buffer=[0, 0], separatrix_buffer=None)

Returns the minimum and maximum valid values of p for given e, x, and a, taking into account the specified buffers.

Parameters:

  • e (float) – Eccentricity.

  • x (float) – Cosine of the inclination angle. Defaults to 1.

  • a (float) – Primary spin parameter. Defaults to 0.

  • p_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid p range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

  • separatrix_buffer (float, optional) – The buffer distance from the separatrix. If not provided, it will be computed using self.separatrix_buffer_dist.

max_e(p, x, a)

Returns the maximum valid value of e for given p, x, and a. For this ODE, there is no upper limit on e, so it returns 1.

isvalid_pex(p=20, e=0, x=1, a=0, p_buffer=[0, 0], e_buffer=[0, 0], x_buffer=[0, 0], a_buffer=[0, 0], separatrix_buffer=None, **kwargs)

Checks if the given values of p, e, x, and a are valid according to the specified buffers and raises an AssertionError if any of them are invalid.

Parameters:

  • p (float) – Semi-latus rectum.

  • e (float) – Eccentricity.

  • x (float) – Cosine of the inclination angle.

  • a (float) – Primary spin parameter.

  • p_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid p range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

  • e_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid e range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

  • x_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid x range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

  • a_buffer (list[float, float], optional) – Two-element list of non-negative floats specifying how much to shrink the valid a range from below and above, respectively. Defaults to [0.0, 0.0] (no buffer).

  • separatrix_buffer (float, optional) – The buffer distance from the separatrix. If not provided, it will be computed using self.separatrix_buffer_dist.

Raises:

AssertionError – If p is outside the valid range given e, x, and a.

evaluate_rhs(y)

Returns the right-hand side of the ODE system for the given state vector y, which contains the parameters (p, e, Y, Phi_phi, Phi_theta, Phi_r).

Parameters:

y (list[float] | ndarray)

Return type:

list[float | ndarray]