In this example we solve the universal Kepler’s equation to determine universal anomaly after dt time with a given initial radius r0, velocity v0 and semimajor axis of a spacecraft.
clc; clear all; mu = 398600; % Earth’s gravitational parameter [km^3/s^2] % Initial conditions r0 = 12000; %[km] Initial radius vr0 = 3; %[km/s] Initial radial speed dt = 7200; %[s] Time interval a = -20000; %[km] Semimajor axis alfa = 1/a; %[km^-1] Reciprocal of the semimajor axis; X0 = mu^0.5*abs(alfa)*dt; %[km^0.5]Initial estimate of X0 Xi = X0; tol = 1E-10; % Tolerance while(1) zi = alfa*Xi^2; [ Cz,Sz] = Stumpff( zi ); fX = r0*vr0/(mu)^0.5*Xi^2*Cz + (1 - alfa*r0)*Xi^3*Sz + r0*Xi -(mu)^0.5*dt; fdX = r0*vr0/(mu)^0.5*Xi*(1 - alfa*Xi^2*Sz) + (1 - alfa*r0)*Xi^2*Cz + r0; eps = fX/fdX; Xi = Xi - eps; if(abs(eps) < tol ) break end end fprintf('Universal anomaly X = %4.3f [km^0.5] \n',Xi) % Equations from the book % Orbital Mechanics for Engineering Students,2nd Edition,Aerospace Engineering
Universal anomaly X = 173.324 [km^0.5]
Small Satellites 