% Kerr space-time
% First integrals and orbits
% Continets configuration
% Time-like geodesics
syms E delta1 e r M Q a;
fir = e*r^3 +2*M*delta1*r^2 -Q*r +2*M*Q;
Delta = r^2 -2*M*r +a^2;
Dr = r*Delta*fir;
L_plus = (2*M*a*E*r +sqrt(Dr))/r/(2*M-r);
L_minus = (2*M*a*E*r -sqrt(Dr))/r/(2*M-r);
E = solve( E^2 -delta1-e, E); E = E(1); % E = (e+delta1)^(1/2)
M = 1; a = sqrt(.84); e = .4; Q = 1; delta1 = 1;
%r = 0:0.01:2*M;
subs(L_plus);
L_plus