Utility Functions
Fundamental Frequencies
Get dimensionless fundamental frequencies from Schmidt 2002:
[1]:
import numpy as np
[2]:
from few.utils.utility import get_fundamental_frequencies
a = 0.5
num = 10
p = np.linspace(8.0, 10.0, num)
e = np.linspace(0.2, 0.4, num)
x = np.linspace(0.3, 0.4, num)
OmegaPhi, OmegaTheta, OmegaR = get_fundamental_frequencies(a, p, e, x)
out = np.array([p, e, x, OmegaPhi, OmegaTheta, OmegaR]).T
print("p\t\te\t\tx\tOmegaPhi\tOmegaTheta\tOmegaR")
print(out)
p e x OmegaPhi OmegaTheta OmegaR
[[ 8. 0.2 0.3 0.04344785 0.04162395 0.02349858]
[ 8.22222222 0.22222222 0.31111111 0.0412017 0.03953685 0.02299122]
[ 8.44444444 0.24444444 0.32222222 0.03905666 0.03753599 0.02241441]
[ 8.66666667 0.26666667 0.33333333 0.03700405 0.03561442 0.02178189]
[ 8.88888889 0.28888889 0.34444444 0.03503633 0.03376608 0.02110469]
[ 9.11111111 0.31111111 0.35555556 0.0331469 0.03198567 0.02039181]
[ 9.33333333 0.33333333 0.36666667 0.03133003 0.03026855 0.01965064]
[ 9.55555556 0.35555556 0.37777778 0.02958068 0.02861065 0.01888732]
[ 9.77777778 0.37777778 0.38888889 0.02789444 0.02700838 0.01810703]
[10. 0.4 0.4 0.02626744 0.02545863 0.01731413]]
Separatrix in Generic Kerr
Get the separatrix in generic Kerr from Stein & Warburton 2020:
[3]:
from few.utils.utility import get_separatrix
a = 0.5
num = 10
e = np.linspace(0.2, 0.4, num)
x = np.linspace(0.3, 0.4, num)
p_sep = get_separatrix(a, e, x)
out = np.array([e, x, p_sep]).T
print("e\t\tx\tseparatrix")
print(out)
e x separatrix
[[0.2 0.3 5.69505166]
[0.22222222 0.31111111 5.71503602]
[0.24444444 0.32222222 5.73483554]
[0.26666667 0.33333333 5.75445904]
[0.28888889 0.34444444 5.7739152 ]
[0.31111111 0.35555556 5.79321253]
[0.33333333 0.36666667 5.81235943]
[0.35555556 0.37777778 5.83136416]
[0.37777778 0.38888889 5.85023483]
[0.4 0.4 5.86897945]]
Converting between the integrals of motion \((E, L, Q)\) and the orbital elements \((p, e, x)\)
We can convert to the three constants of motion \((E,L,Q)\) in generic Kerr spacetime, and go back to \((p, e, x)\) via Hughes 2024:
[4]:
from few.utils.utility import get_kerr_geo_constants_of_motion, ELQ_to_pex
a = 0.5
num = 11
p = np.linspace(8.0, 10.0, num)
e = np.linspace(0.2, 0.4, num)
x = np.linspace(0.3, 0.4, num)
print("To ELQ:")
E, L, Q = get_kerr_geo_constants_of_motion(a, p, e, x)
out = np.array([p, e, x, E, L, Q]).T
print("p\t\te\t\tx\tE\tL\tQ")
print(out)
print("To pex:")
p_new, e_new, x_new = ELQ_to_pex(a, E, L, Q)
out = np.array([p_new, e_new, x_new, E, L, Q]).T
print("p\t\te\t\tx\tE\tL\tQ")
print(out)
To ELQ:
p e x E L Q
[[ 8. 0.2 0.3 0.94869805 1.0488073 11.14493312]
[ 8.2 0.22 0.31 0.94995592 1.09042731 11.20587846]
[ 8.4 0.24 0.32 0.9512292 1.1326752 11.26724045]
[ 8.6 0.26 0.33 0.95251428 1.17554328 11.32839477]
[ 8.8 0.28 0.34 0.95380826 1.21902438 11.38878189]
[ 9. 0.3 0.35 0.95510877 1.26311183 11.44789629]
[ 9.2 0.32 0.36 0.95641391 1.30779937 11.50527773]
[ 9.4 0.34 0.37 0.9577221 1.35308111 11.5605041 ]
[ 9.6 0.36 0.38 0.9590321 1.3989515 11.61318563]
[ 9.8 0.38 0.39 0.96034289 1.44540531 11.66296002]
[10. 0.4 0.4 0.96165365 1.49243756 11.70948848]]
To pex:
p e x E L Q
[[ 8. 0.2 0.3 0.94869805 1.0488073 11.14493312]
[ 8.2 0.22 0.31 0.94995592 1.09042731 11.20587846]
[ 8.4 0.24 0.32 0.9512292 1.1326752 11.26724045]
[ 8.6 0.26 0.33 0.95251428 1.17554328 11.32839477]
[ 8.8 0.28 0.34 0.95380826 1.21902438 11.38878189]
[ 9. 0.3 0.35 0.95510877 1.26311183 11.44789629]
[ 9.2 0.32 0.36 0.95641391 1.30779937 11.50527773]
[ 9.4 0.34 0.37 0.9577221 1.35308111 11.5605041 ]
[ 9.6 0.36 0.38 0.9590321 1.3989515 11.61318563]
[ 9.8 0.38 0.39 0.96034289 1.44540531 11.66296002]
[10. 0.4 0.4 0.96165365 1.49243756 11.70948848]]
Convert between \(x_I\) and \(Y\)
\(Y\equiv\cos{\iota}=L/\sqrt{L^2 + Q}\) is different than \(x_I\equiv \cos{I}\), which is accepted for relativistic waveforms and in the generic waveform interface discussed above. \(I\) is the inclination angle of the orbital plane to the to equatorial plane.
[5]:
from few.utils.pn_map import xI_to_Y, Y_to_xI
a = 0.5
num = 11
p = np.linspace(8.0, 10.0, num)
e = np.linspace(0.2, 0.4, num)
x = np.linspace(0.3, 0.4, num)
Y = xI_to_Y(a, p, e, x)
out = np.array([p, e, x, Y]).T
print("To Y:\np\t\te\t\tx\tY")
print(out)
x_new = Y_to_xI(a, p, e, Y)
out = np.array([p, e, x_new, Y]).T
print("To x:\np\t\te\t\tx\tY")
print(out)
To Y:
p e x Y
[[ 8. 0.2 0.3 0.29972126]
[ 8.2 0.22 0.31 0.30972412]
[ 8.4 0.24 0.32 0.31972764]
[ 8.6 0.26 0.33 0.32973176]
[ 8.8 0.28 0.34 0.33973642]
[ 9. 0.3 0.35 0.34974158]
[ 9.2 0.32 0.36 0.35974718]
[ 9.4 0.34 0.37 0.36975319]
[ 9.6 0.36 0.38 0.37975956]
[ 9.8 0.38 0.39 0.38976626]
[10. 0.4 0.4 0.39977325]]
To x:
p e x Y
[[ 8. 0.2 0.3 0.29972126]
[ 8.2 0.22 0.31 0.30972412]
[ 8.4 0.24 0.32 0.31972764]
[ 8.6 0.26 0.33 0.32973176]
[ 8.8 0.28 0.34 0.33973642]
[ 9. 0.3 0.35 0.34974158]
[ 9.2 0.32 0.36 0.35974718]
[ 9.4 0.34 0.37 0.36975319]
[ 9.6 0.36 0.38 0.37975956]
[ 9.8 0.38 0.39 0.38976626]
[10. 0.4 0.4 0.39977325]]
Get \(\mu\) based on desired duration of waveform
If you have a desired length of waveform to analyze, this function will give you the value of \(\mu\) that corresponds to the proper evolution time.
[6]:
from few.trajectory.inspiral import EMRIInspiral
from few.trajectory.ode import SchwarzEccFlux
from few.utils.utility import get_mu_at_t
traj_module = EMRIInspiral(func=SchwarzEccFlux)
# set initial parameters
M = 1e6
p0 = 11.0
e0 = 0.7
traj_args = [M, 0.0, p0, e0, 1.0]
traj_kwargs = {}
index_of_mu = 1
t_out = 1.5
# run trajectory
mu_new = get_mu_at_t(
traj_module,
t_out,
traj_args,
index_of_mu=index_of_mu,
traj_kwargs=traj_kwargs,
xtol=2e-12,
rtol=8.881784197001252e-16,
bounds=None,
)
print(
"mu = {} will create a waveform that is {} years long, given the other input parameters.".format(
mu_new, t_out
)
)
mu = 18.804307109481528 will create a waveform that is 1.5 years long, given the other input parameters.
Get \(p_0\) based on desired duration of waveform
If you have a desired length of waveform to analyze, this function will give you the value of \(p_0\) that corresponds to the proper evolution time.
[7]:
from few.utils.utility import get_p_at_t
traj_module = EMRIInspiral(func=SchwarzEccFlux)
# set initial parameters
M = 1e6
mu = 5e1
e0 = 0.7
traj_args = [M, mu, 0.0, e0, 1.0]
traj_kwargs = {}
index_of_p = 3
t_out = 1.5
# run trajectory
p_new = get_p_at_t(
traj_module,
t_out,
traj_args,
index_of_p=3,
index_of_a=2,
index_of_e=4,
index_of_x=5,
traj_kwargs={},
xtol=2e-12,
rtol=8.881784197001252e-16,
bounds=None,
)
print(
"p0 = {} will create a waveform that is {} years long, given the other input parameters.".format(
p_new, t_out
)
)
p0 = 13.059789978787903 will create a waveform that is 1.5 years long, given the other input parameters.