Skip to content

integrate_backwards=True in EMRIInspiral returns inverted physical phases #183

@perturber

Description

@perturber

From FEWv2.0.0 tag:

When generating a backward trajectory using the integrate_backwards=True flag, the initial phase signs are flipped to use the same integrator module (see few.trajectory.inspiral.EMRIInspiral.get_inspiral Line 242).

However, when returning the trajectory in Line 269, the phases are not flipped back, which results in $-1 \Phi_{\alpha}(t)$ as the output. Here, $\Phi_{\alpha}(t)$ are the time evolutions of the true initial physical phases with $\alpha = \phi, \theta, r$.

Reproduce the bug:

import numpy as np
from few.trajectory.inspiral import EMRIInspiral
from few.trajectory.ode import KerrEccEqFlux
from few.utils.utility import get_p_at_t

kerr_traj = EMRIInspiral(func=KerrEccEqFlux)
T = 1.0
m1, m2, a, e0, xI0 = 1e6, 1e1, 0.9, 0.2, 1.0
Phi_phi0, Phi_theta0, Phi_r0 = 1.0, 0.0, 1.0
p0 = 12.0

# 1. Forward integration
traj = kerr_traj(m1, m2, a, p0, e0, xI0, 
                 Phi_phi0=Phi_phi0, Phi_theta0=Phi_theta0, Phi_r0=Phi_r0, 
                 T=T, err_tol=1e-11)

# 2. Backward integration from the final state
traj_back = kerr_traj(m1, m2, a, traj[1][-1], traj[2][-1], traj[3][-1], 
                      Phi_phi0=-traj[4][-1], # <= user must supply the negative of the final phases
                      Phi_theta0=-traj[5][-1], 
                      Phi_r0=-traj[6][-1], 
                      T=T, upsample=False, integrate_backwards=True, err_tol=1e-11)

# The raw returned phase is the negative of the true physical phase
print(f"Input Initial Phase: {Phi_phi0}")
print(f"Raw Returned t=0 Phase: {traj_back[4][0]}")
print(f"Negated Returned t=0 Phase: {-traj_back[4][0]}")

Output:

Input Initial Phase: 1.0
Raw Returned t=0 Phase: -0.9999967547482811
Negated Returned t=0 Phase: 0.9999967547482811

Proposed fix:
Add a condition to negate the phase columns of the out array right before return:

if temp_kwargs["integrate_backwards"]:
    out[:, 4] *= -1
    out[:, 5] *= -1
    out[:, 6] *= -1

t, p, e, x, Phi_phi, Phi_theta, Phi_r = out.T.copy()
return t, p, e, x, Phi_phi, Phi_theta, Phi_r

Metadata

Metadata

Assignees

No one assigned

    Labels

    bugSomething isn't working

    Type

    No fields configured for Bug.

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions