Skip to content

Commit

Permalink
Merge pull request #34 from BlackHolePerturbationToolkit/ExposeDeltas
Browse files Browse the repository at this point in the history 
Expose deltas
  • Loading branch information
@barrywardell
barrywardell authored Jul 29, 2022
2 parents 764df0b + 8c33407 commit ee02610
Showing 1 changed file with 77 additions and 37 deletions.
114 changes: 77 additions & 37 deletions Kernel/KerrGeoOrbit.m
View file
Original file line number Diff line number Diff line change
Expand Up @@ -151,7 +151,8 @@
"Inclination" -> x,
"Energy" -> En,
"AngularMomentum" -> L,
"CarterConstant" -> Q
"CarterConstant" -> Q,
"InitialPhases" -> initPhases
];

KerrGeoOrbitFunction[a, p, e, x, assoc]
Expand Down Expand Up @@ -302,7 +303,8 @@
"a" -> a,
"p" -> p,
"e" -> e,
"Inclination" -> x
"Inclination" -> x,
"InitialPhases" -> initPhases
];

KerrGeoOrbitFunction[a,p,e,x,assoc]
Expand Down Expand Up @@ -391,7 +393,8 @@
"a" -> a,
"p" -> p,
"e" -> e,
"Inclination" -> x
"Inclination" -> x,
"InitialPhases" -> initPhases
];

KerrGeoOrbitFunction[a,p,e,x,assoc]
Expand Down Expand Up @@ -431,13 +434,14 @@
"PolarFrequency" -> \[CapitalUpsilon]\[Theta],
"AzimuthalFrequency" -> \[CapitalUpsilon]\[Phi],
"RadialRoots" -> {r1,r2,r3,r4},
"Frequencies" -> <|"\!\(\*SubscriptBox[\(\[CapitalUpsilon]\), \(r\)]\)" -> \[CapitalUpsilon]r, "\!\(\*SubscriptBox[\(\[CapitalUpsilon]\), \(\[Theta]\)]\)" -> \[CapitalUpsilon]\[Theta], "\!\(\*SubscriptBox[\(\[CapitalUpsilon]\), \(\[Phi]\)]\)" -> \[CapitalUpsilon]\[Phi] |>,
"Frequencies" -> <|"\!\(\*SubscriptBox[\(\[CapitalUpsilon]\), \(r\)]\)" -> \[CapitalUpsilon]r, "\!\(\*SubscriptBox[\(\[CapitalUpsilon]\), \(\[Theta]\)]\)" -> \[CapitalUpsilon]\[Theta], "\!\(\*SubscriptBox[\(\[CapitalUpsilon]\), \(\[Phi]\)]\)" -> \[CapitalUpsilon]\[Phi], "\!\(\*SubscriptBox[\(\[CapitalUpsilon]\), \(t\)]\)" -> \[CapitalUpsilon]t |>,
"Trajectory" -> {t,r,\[Theta],\[Phi]},
"FourVelocity" -> velocity,
"a" -> a,
"p" -> p,
"e" -> e,
"Inclination" -> x
"Inclination" -> x,
"InitialPhases" -> initPhases
];

KerrGeoOrbitFunction[a, p, e, x, assoc]
Expand All @@ -450,7 +454,7 @@
(*Generic (Mino)*)


KerrGeoOrbitMino[a_,p_,e_,x_,initPhases:{_,_,_,_}:{0,0,0,0}]:=Module[{M=1,consts,En,L,Q,assoc,\[CapitalUpsilon]r,\[CapitalUpsilon]\[Theta],\[CapitalUpsilon]\[Phi],\[CapitalUpsilon]t,r1,r2,r3,r4,zp,zm,kr,k\[Theta],rp,rm,hr,hp,hm,rq,zq,\[Psi]r,tr,\[Phi]f,\[Psi]z,tz,\[Phi]z,qt0,qr0,qz0,q\[Phi]0,t,r,\[Theta],\[Phi],\[Phi]t,\[Phi]r,Ct,C\[Phi],velocity},
KerrGeoOrbitPhases[a_,p_,e_,x_]:=Module[{M=1,consts,En,L,Q,assoc,\[CapitalUpsilon]r,\[CapitalUpsilon]\[Theta],\[CapitalUpsilon]\[Phi],\[CapitalUpsilon]t,r1,r2,r3,r4,zp,zm,kr,k\[Theta],rp,rm,hr,hp,hm,rq,zq,\[Psi]r,tr,\[Phi]f,\[Psi]z,tz,\[Phi]z,t,r,\[Theta],\[Phi],\[Phi]t,\[Phi]r,Ct,C\[Phi],velocity,\[CapitalDelta]tr,\[CapitalDelta]t\[Theta],\[CapitalDelta]\[Phi]r,\[CapitalDelta]\[Phi]\[Theta]},
consts = KerrGeoConstantsOfMotion[a,p,e,x];
{En,L,Q} = Values[consts];
{\[CapitalUpsilon]r,\[CapitalUpsilon]\[Theta],\[CapitalUpsilon]\[Phi],\[CapitalUpsilon]t} = Values[KerrGeodesics`OrbitalFrequencies`Private`KerrGeoMinoFrequencies[a,p,e,x]];
Expand Down Expand Up @@ -487,47 +491,80 @@
tz[qz_]:= 1/(1-En^2) En zp ( EllipticE[k\[Theta]]2((qz+\[Pi]/2)/\[Pi])-EllipticE[\[Psi]z[qz],k\[Theta]]);
\[Phi]z[qz_]:= -1/zp L ( EllipticPi[zm^2,k\[Theta]]2((qz+\[Pi]/2)/\[Pi])-EllipticPi[zm^2,\[Psi]z[qz],k\[Theta]]);

{qt0, qr0, qz0, q\[Phi]0} = {initPhases[[1]], initPhases[[2]], initPhases[[3]], initPhases[[4]]};
(*Calculate normalization constants so that t=0 and \[Phi]=0 at \[Lambda]=0 when qt0=0 and q\[Phi]0=0 *)
Ct=tr[qr0]+tz[qz0]/.i_/;i==0:>0;
C\[Phi]=\[Phi]r[qr0]+\[Phi]z[qz0]/.i_/;i==0:>0;

t=Function[{Global`\[Lambda]}, Evaluate[ qt0 + \[CapitalUpsilon]t Global`\[Lambda] + tr[\[CapitalUpsilon]r Global`\[Lambda] + qr0] + tz[\[CapitalUpsilon]\[Theta] Global`\[Lambda] + qz0]-Ct ], Listable];
r=Function[{Global`\[Lambda]}, Evaluate[ rq[\[CapitalUpsilon]r Global`\[Lambda]+ qr0] ], Listable];
\[Theta]=Function[{Global`\[Lambda]}, Evaluate[ ArcCos[zq[\[CapitalUpsilon]\[Theta] Global`\[Lambda] + qz0]] ], Listable];
\[Phi]=Function[{Global`\[Lambda]}, Evaluate[ q\[Phi]0 + \[CapitalUpsilon]\[Phi] Global`\[Lambda] + \[Phi]r[\[CapitalUpsilon]r Global`\[Lambda]+ qr0] + \[Phi]z[\[CapitalUpsilon]\[Theta] Global`\[Lambda] + qz0]-C\[Phi] ], Listable];

velocity = Values[KerrGeoFourVelocity[a,p,e,x,{initPhases[[2]],initPhases[[3]]}]];
t=Function[{Global`qt,Global`qr,Global`q\[Theta]}, Evaluate[ Global`qt+ tr[Global`qr] + tz[Global`q\[Theta]]], Listable];
r=Function[{Global`qr}, Evaluate[ rq[Global`qr] ], Listable];
\[Theta]=Function[{Global`q\[Theta]}, Evaluate[ ArcCos[zq[Global`q\[Theta]]] ], Listable];
\[Phi]=Function[{Global`q\[Phi],Global`qr,Global`q\[Theta]}, Evaluate[ Global`q\[Phi] + \[Phi]r[Global`qr] + \[Phi]z[Global`q\[Theta]]], Listable];

\[CapitalDelta]tr=Function[{Global`qr}, Evaluate[ tr[Global`qr] ], Listable];
\[CapitalDelta]t\[Theta]=Function[{Global`q\[Theta]}, Evaluate[ tz[Global`q\[Theta]] ], Listable];
\[CapitalDelta]\[Phi]r=Function[{Global`qr}, Evaluate[ \[Phi]r[Global`qr] ], Listable];
\[CapitalDelta]\[Phi]\[Theta]=Function[{Global`q\[Theta]}, Evaluate[ \[Phi]z[Global`q\[Theta]] ], Listable];

assoc = Association[
"a" -> a,
"p" -> p,
"e" -> e,
"Inclination" -> x,
"Parametrization"->"Mino",
"Parametrization"->"Phases",
"Energy" -> En,
"AngularMomentum" -> L,
"CarterConstant" -> Q,
"ConstantsOfMotion" -> consts,
"RadialFrequency" -> \[CapitalUpsilon]r,
"PolarFrequency" -> \[CapitalUpsilon]\[Theta],
"AzimuthalFrequency" -> \[CapitalUpsilon]\[Phi],
"Frequencies" -> <|"\!\(\*SubscriptBox[\(\[CapitalUpsilon]\), \(r\)]\)" -> \[CapitalUpsilon]r, "\!\(\*SubscriptBox[\(\[CapitalUpsilon]\), \(\[Theta]\)]\)" -> \[CapitalUpsilon]\[Theta], "\!\(\*SubscriptBox[\(\[CapitalUpsilon]\), \(\[Phi]\)]\)" -> \[CapitalUpsilon]\[Phi] |> ,
"Frequencies" -> <|"\!\(\*SubscriptBox[\(\[CapitalUpsilon]\), \(r\)]\)" -> \[CapitalUpsilon]r, "\!\(\*SubscriptBox[\(\[CapitalUpsilon]\), \(\[Theta]\)]\)" -> \[CapitalUpsilon]\[Theta], "\!\(\*SubscriptBox[\(\[CapitalUpsilon]\), \(\[Phi]\)]\)" -> \[CapitalUpsilon]\[Phi] ,"\!\(\*SubscriptBox[\(\[CapitalUpsilon]\), \(t\)]\)" -> \[CapitalUpsilon]t |> ,
"Trajectory" -> {t,r,\[Theta],\[Phi]},
"FourVelocity"-> velocity,
"TrajectoryDeltas" -> <|"\[CapitalDelta]tr"-> \[CapitalDelta]tr,"\[CapitalDelta]t\[Theta]"-> \[CapitalDelta]t\[Theta],"\[CapitalDelta]\[Phi]r"-> \[CapitalDelta]\[Phi]r,"\[CapitalDelta]\[Phi]\[Theta]"-> \[CapitalDelta]\[Phi]\[Theta]|>,
"RadialRoots" -> {r1,r2,r3,r4},
"a" -> a,
"p" -> p,
"e" -> e,
"Inclination" -> x
"Inclination" -> x,
"InitialPhases" -> {0, 0, 0, 0}
];

KerrGeoOrbitFunction[a,p,e,x,assoc]

]


KerrGeoOrbitMino[a_,p_,e_,x_,initPhases:{_,_,_,_}:{0,0,0,0}]:=Module[{M=1,assoc,qt0,qr0,qz0,q\[Phi]0,tr,tz,t,r,\[Theta],\[Phi],\[Phi]z,\[Phi]r,Ct,C\[Phi],velocity,\[CapitalUpsilon]r,\[CapitalUpsilon]\[Theta],\[CapitalUpsilon]\[Phi],\[CapitalUpsilon]t},
assoc= Last@KerrGeoOrbitPhases[a,p,e,x];

{qt0, qr0, qz0, q\[Phi]0} = {initPhases[[1]], initPhases[[2]], initPhases[[3]], initPhases[[4]]};

{\[CapitalUpsilon]r,\[CapitalUpsilon]\[Theta],\[CapitalUpsilon]\[Phi],\[CapitalUpsilon]t}=Values@assoc["Frequencies"];
{t,r,\[Theta],\[Phi]}=assoc["Trajectory"];

tr=assoc["TrajectoryDeltas"]["\[CapitalDelta]tr"];
tz=assoc["TrajectoryDeltas"]["\[CapitalDelta]t\[Theta]"];
\[Phi]r=assoc["TrajectoryDeltas"]["\[CapitalDelta]\[Phi]r"];
\[Phi]z=assoc["TrajectoryDeltas"]["\[CapitalDelta]\[Phi]\[Theta]"];

(*Calculate normalization constants so that t=0 and \[Phi]=0 at \[Lambda]=0 when qt0=0 and q\[Phi]0=0 *)
Ct=tr[qr0]+tz[qz0]/.i_/;i==0:>0;
C\[Phi]=\[Phi]r[qr0]+\[Phi]z[qz0]/.i_/;i==0:>0;

t=Function[{Global`\[Lambda]}, Evaluate[ t[qt0 + \[CapitalUpsilon]t Global`\[Lambda], \[CapitalUpsilon]r Global`\[Lambda] + qr0, \[CapitalUpsilon]\[Theta] Global`\[Lambda] + qz0]-Ct ], Listable];
r=Function[{Global`\[Lambda]}, Evaluate[ r[\[CapitalUpsilon]r Global`\[Lambda]+ qr0] ], Listable];
\[Theta]=Function[{Global`\[Lambda]}, Evaluate[ \[Theta][\[CapitalUpsilon]\[Theta] Global`\[Lambda] + qz0] ], Listable];
\[Phi]=Function[{Global`\[Lambda]}, Evaluate[ \[Phi][q\[Phi]0 + \[CapitalUpsilon]\[Phi] Global`\[Lambda] ,\[CapitalUpsilon]r Global`\[Lambda]+ qr0,\[CapitalUpsilon]\[Theta] Global`\[Lambda] + qz0]-C\[Phi] ], Listable];

velocity = Values[KerrGeoFourVelocity[a,p,e,x,{qr0,qz0}]];

assoc["Parametrization"] = "Mino";
assoc["Trajectory"] = {t,r,\[Theta],\[Phi]};
assoc["FourVelocity"] = velocity;
assoc["InitialPhases"] = {qt0, qr0, qz0, q\[Phi]0};

KerrGeoOrbitFunction[a,p,e,x,assoc]

]


(* ::Subsection::Closed:: *)
(*Generic (Fast Spec - Mino)*)

Expand Down Expand Up @@ -599,7 +636,7 @@


(* ::Subsubsection::Closed:: *)
(*Subroutines for calculating \[CapitalDelta]\[Phi]r(\[Lambda]), \[CapitalDelta]\[Phi]\[Theta](\[Lambda]), \[CapitalDelta]tr(\[Lambda]), \[CapitalDelta]t\[Theta](\[Lambda])*)
(*Subroutines for calculating Subscript[\[CapitalDelta]\[Phi], r](Subscript[q, r]), Subscript[\[CapitalDelta]\[Phi], \[Theta]](Subscript[q, \[Theta]]), Subscript[\[CapitalDelta]t, r](Subscript[q, r]), Subscript[\[CapitalDelta]t, \[Theta]](Subscript[q, \[Theta]])*)


PhiOfMinoFastSpecR[a_,p_,e_,x_,{\[CapitalUpsilon]r_,minoSampleR_}]:=
Expand Down Expand Up @@ -890,7 +927,7 @@
f[n_]:=fList[[n+1]];

(* Construct integrated series solution *)
integratedF[mino_]:=2*Sum[f[n]/n Sin[n freq mino],{n,1,sampleMax}];
integratedF[phase_]:=2*Sum[f[n]/n Sin[n phase],{n,1,sampleMax}];
(* Allow function to evaluate lists by threading over them *)

integratedF
Expand All @@ -902,7 +939,7 @@


Clear[KerrGeoOrbitFastSpec];
Options[KerrGeoOrbitFastSpec]={InitialPosition->{0,0,0,0}};
Options[KerrGeoOrbitFastSpec]={"InitialPosition"->{0,0,0,0}};
KerrGeoOrbitFastSpec[a_,p_,e_,x_,initPhases:{_,_,_,_}:{0,0,0,0},opts:OptionsPattern[]]:=
Module[{M=1,consts,En,L,Q,\[CapitalUpsilon]r,\[CapitalUpsilon]\[Theta],\[CapitalUpsilon]\[Phi],\[CapitalUpsilon]t,r1,r2,r3,r4,p3,p4,\[Alpha],\[Beta],zp,zm,assoc,var,\[Chi]0,\[Psi]0,
r0,\[Theta]0,qt0,qr0,q\[Theta]0,q\[Phi]0,\[Lambda]t0,\[Lambda]r0,\[Lambda]\[Theta]0,\[Lambda]\[Phi]0,t,r,\[Theta],\[Phi],\[Psi],\[Chi],\[CapitalDelta]\[Lambda]r,\[Lambda]r,\[CapitalDelta]\[Lambda]\[Theta],\[Lambda]\[Theta],rC,\[Theta]C,\[CapitalDelta]r,\[CapitalDelta]\[Theta],
Expand Down Expand Up @@ -953,35 +990,35 @@
using spectral integration*)
If[e>0,
\[Psi]=PhaseOfMinoFastSpec[\[CapitalUpsilon]r,\[Lambda]rSample],
\[Psi][\[Lambda]_]:=0
\[Psi][qr_]:=0
];
If[x^2<1,
\[Chi]=PhaseOfMinoFastSpec[\[CapitalUpsilon]\[Theta],\[Lambda]\[Theta]Sample],
\[Chi][\[Lambda]_]:=0
\[Chi][q\[Theta]_]:=0
];

(*Spectral integration of t and \[Phi] as functions of \[Lambda]*)
If[e>0,
\[CapitalDelta]tr=TimeOfMinoFastSpecR[a,p,e,x,{\[CapitalUpsilon]r,\[Lambda]rSample}];
\[CapitalDelta]\[Phi]r=PhiOfMinoFastSpecR[a,p,e,x,{\[CapitalUpsilon]r,\[Lambda]rSample}],
\[CapitalDelta]tr[\[Lambda]_]:=0;
\[CapitalDelta]\[Phi]r[\[Lambda]_]:=0;
\[CapitalDelta]tr[qr_]:=0;
\[CapitalDelta]\[Phi]r[qr_]:=0;
];
If[x^2<1,
(* Calculate theta dependence for non-equatorial orbits *)
\[CapitalDelta]t\[Theta]=TimeOfMinoFastSpecTheta[a,p,e,x,{\[CapitalUpsilon]\[Theta],\[Lambda]\[Theta]Sample}];
\[CapitalDelta]\[Phi]\[Theta]=PhiOfMinoFastSpecTheta[a,p,e,x,{\[CapitalUpsilon]\[Theta],\[Lambda]\[Theta]Sample}],
(* No theta dependence for equatorial orbits *)
\[CapitalDelta]t\[Theta][\[Lambda]_]:=0;
\[CapitalDelta]\[Phi]\[Theta][\[Lambda]_]:=0;
\[CapitalDelta]t\[Theta][q\[Theta]_]:=0;
\[CapitalDelta]\[Phi]\[Theta][q\[Theta]_]:=0;
];

(*Collect initial Mino time phases*)
{qt0, qr0, q\[Theta]0, q\[Phi]0} = {initPhases[[1]], initPhases[[2]], initPhases[[3]], initPhases[[4]]};

(* If the user specifies a valid set of initial coordinate positions,
find phases that give these initial positions *)
{tInit,rInit,\[Theta]Init,\[Phi]Init}=OptionValue[InitialPosition];
{tInit,rInit,\[Theta]Init,\[Phi]Init}=OptionValue["InitialPosition"];
If[{tInit,rInit,\[Theta]Init,\[Phi]Init}!={0,0,0,0},
If[rInit<=r1&&rInit>=r2&&\[Theta]Init>=zRoots[[1]]&&\[Theta]Init<=zRoots[[2]],
If[e==0,\[Psi]0=0,\[Psi]0=ArcCos[p M/(e rInit)-1/e]];
Expand All @@ -997,13 +1034,13 @@
If[\[CapitalUpsilon]\[Theta]==0,\[Lambda]\[Theta]0=0,\[Lambda]\[Theta]0=q\[Theta]0/\[CapitalUpsilon]\[Theta]];

(* Find integration constants for t and \[Phi], so that t(\[Lambda]=0)=qt0 and \[Phi](\[Lambda]=0)=q\[Phi]0 *)
\[Phi]C=\[CapitalDelta]\[Phi]r[\[Lambda]r0]+\[CapitalDelta]\[Phi]\[Theta][\[Lambda]\[Theta]0];
tC=\[CapitalDelta]tr[\[Lambda]r0]+\[CapitalDelta]t\[Theta][\[Lambda]\[Theta]0];
\[Phi]C=\[CapitalDelta]\[Phi]r[qr0]+\[CapitalDelta]\[Phi]\[Theta][q\[Theta]0];
tC=\[CapitalDelta]tr[qr0]+\[CapitalDelta]t\[Theta][q\[Theta]0];

t=Function[{Global`\[Lambda]}, Evaluate[ \[CapitalDelta]tr[Global`\[Lambda]+\[Lambda]r0]+\[CapitalDelta]t\[Theta][Global`\[Lambda]+\[Lambda]\[Theta]0]+\[CapitalUpsilon]t Global`\[Lambda]+qt0-tC ], Listable];
r=Function[{Global`\[Lambda]}, Evaluate[ r0[\[Psi][Global`\[Lambda]+\[Lambda]r0]+\[CapitalUpsilon]r Global`\[Lambda]+qr0] ], Listable];
\[Theta]=Function[{Global`\[Lambda]}, Evaluate[ \[Theta]0[\[Chi][Global`\[Lambda]+\[Lambda]\[Theta]0]+\[CapitalUpsilon]\[Theta] Global`\[Lambda]+q\[Theta]0] ], Listable];
\[Phi]=Function[{Global`\[Lambda]}, Evaluate[ \[CapitalDelta]\[Phi]r[Global`\[Lambda]+\[Lambda]r0]+\[CapitalDelta]\[Phi]\[Theta][Global`\[Lambda]+\[Lambda]\[Theta]0]+\[CapitalUpsilon]\[Phi] Global`\[Lambda]+q\[Phi]0-\[Phi]C ], Listable];
t=Function[{Global`\[Lambda]}, Evaluate[ \[CapitalDelta]tr[\[CapitalUpsilon]r Global`\[Lambda] + qr0]+\[CapitalDelta]t\[Theta][\[CapitalUpsilon]\[Theta] Global`\[Lambda] + q\[Theta]0]+\[CapitalUpsilon]t Global`\[Lambda]+qt0-tC ], Listable];
r=Function[{Global`\[Lambda]}, Evaluate[ r0[\[Psi][\[CapitalUpsilon]r Global`\[Lambda] + qr0]+\[CapitalUpsilon]r Global`\[Lambda]+qr0] ], Listable];
\[Theta]=Function[{Global`\[Lambda]}, Evaluate[ \[Theta]0[\[Chi][\[CapitalUpsilon]\[Theta] Global`\[Lambda]+ q\[Theta]0]+\[CapitalUpsilon]\[Theta] Global`\[Lambda]+q\[Theta]0] ], Listable];
\[Phi]=Function[{Global`\[Lambda]}, Evaluate[ \[CapitalDelta]\[Phi]r[\[CapitalUpsilon]r Global`\[Lambda] + qr0]+\[CapitalDelta]\[Phi]\[Theta][\[CapitalUpsilon]\[Theta] Global`\[Lambda] + q\[Theta]0]+\[CapitalUpsilon]\[Phi] Global`\[Lambda]+q\[Phi]0-\[Phi]C ], Listable];

velocity = Values[KerrGeoFourVelocity[a,p,e,x,{initPhases[[2]],initPhases[[3]]}]];

Expand All @@ -1018,13 +1055,15 @@
"AzimuthalFrequency" -> \[CapitalUpsilon]\[Phi],
"Frequencies" -> <|"\!\(\*SubscriptBox[\(\[CapitalUpsilon]\), \(r\)]\)" -> \[CapitalUpsilon]r, "\!\(\*SubscriptBox[\(\[CapitalUpsilon]\), \(\[Theta]\)]\)" -> \[CapitalUpsilon]\[Theta], "\!\(\*SubscriptBox[\(\[CapitalUpsilon]\), \(\[Phi]\)]\)" -> \[CapitalUpsilon]\[Phi] |> ,
"Trajectory" -> {t,r,\[Theta],\[Phi]},
"TrajectoryDeltas" -> <|"\[CapitalDelta]tr"-> \[CapitalDelta]tr,"\[CapitalDelta]t\[Theta]"-> \[CapitalDelta]t\[Theta],"\[CapitalDelta]\[Phi]r"-> \[CapitalDelta]\[Phi]r,"\[CapitalDelta]\[Phi]\[Theta]"-> \[CapitalDelta]\[Phi]\[Theta]|>,
"FourVelocity"-> velocity,
"RadialRoots" -> {r1,r2,r3,r4},
"PolarRoots" -> zRoots,
"a" -> a,
"p" -> p,
"e" -> e,
"Inclination" -> x
"Inclination" -> x,
"InitialPhases" -> {qt0, qr0, q\[Theta]0, q\[Phi]0}
];

KerrGeoOrbitFunction[a,p,e,x,assoc]
Expand Down Expand Up @@ -1062,6 +1101,7 @@
If[method == "Analytic",
(*Changed "KerrGeoOrbitDarwin" to "KerrGeoOrbitEquatorialDarwin"*)
If[param == "Mino", Return[KerrGeoOrbitMino[a, p, e, x, initPhases]]];
If[param == "Phases", Return[KerrGeoOrbitPhases[a, p, e, x]]];
If[param == "Darwin",
If[PossibleZeroQ[a], Return[KerrGeoOrbitSchwarzDarwin[p, e]], Return[KerrGeoOrbitEquatorialDarwin[a,p,e,x,initPhases]]]
];
Expand Down

0 comments on commit ee02610

Please sign in to comment.