Chapter 3
Astrodynamics

3.1 Ellipse main parameters

Figure 3.1: Ellipse

3.2 Table of common equations

The following table contains the common relations to use for elliptic motion. Equation of ellipse is \(\frac {x^{2}}{a^{2}}+\frac {y^{2}}{b^{2}}=1\)

term to find	relation
conversion between \(E\) and \(\theta \)	\(\begin {aligned} \tan \left (\frac {\theta }{2}\right ) &=\sqrt {\frac {1+e}{1-e}}\tan \left ( \frac {E}{2}\right ) \\ \cos E &= \frac {e+\cos \theta }{1+e\cos \theta }\\ \cos \theta &= \frac {e-\cos E}{e\cos E-1} \end {aligned} \)
position of satellite at time \(t\) Solve for \(E\), then find \(\theta \). \(\tau \) here is time at \(perigee\) and \(n\) is mean satellite speed.	\(E-e\sin E=n\left ( t-\tau \right )\)
eccentricity \(e\)	\(e=\frac {c}{a}=\frac {r_{a}-r_{p}}{r_{a}+r_{p}} =\sqrt {1+\frac {2\mathcal {E}h^{2}}{\mu ^{2}}}\)
Major axes \(a\)	\(\begin {aligned} a &= \frac {r_p (1+e)}{1-e^2} \\ &= \frac {r_a (1-e)}{1-e^2} \\ &= -\frac {\mu }{2\mathcal {E}}\\ &= \sqrt {b^{2}+c^{2}}\\ &= \frac {p}{1-e^{2}} \end {aligned} \)
Minor axes \(b\)	\(\begin {aligned} b &=a\sqrt {1-e^{2}} \end {aligned} \)
\(r_p\)	\(\begin {aligned} r_p &= \frac {a (1-e^2)}{1+e}\\ &= a(1-e) \end {aligned} \)
\(r_a\)	\(\begin {aligned} r_a &= \frac {a (1-e^2)}{1-e}\\ &= a\left ( 1+e\right ) \\ &=\frac {h^{2}}{\mu }\frac {1}{1-e} \end {aligned} \) \(\frac {r_{p}}{r_{a}} =\frac {1-e}{1+e}\)
\(p\)	\(p=a\left ( 1-e^{2}\right ) =\frac {h^{2}}{\mu }=r_{p}\left ( 1+e\right ) =r_{a}\left ( 1-e\right ) \)
specific angular momentum \(h\)	\(\begin {aligned} h &= r_{p}v_{p}=r_{a}v_{a}=\vec {r}\times \vec {v}=\sqrt {p\mu }\\ h &=\sqrt {\mu r} \quad \text {(circular orbit)} \end {aligned} \)
Total Energy \(\mathcal {E}\)	\(\mathcal {E}=\frac {v^{2}}{2}-\frac {\mu }{r} =-\frac {\mu }{2a}\)
velocity \(v\)	\(\begin {aligned} v &= \sqrt {\mu \left ( \frac {2}{r}-\frac {1}{a}\right ) } \quad \text {(vis-viva)}\\ v_{escape} &=\sqrt {\frac {2\mu }{r}} \quad \text {(escape velocity for parabola)}\\ v_{radial} &=\sqrt {\frac {\mu }{p}}e\sin \theta \\ v_{normal} &=\sqrt {\frac {\mu }{p}}\left ( 1+e\cos \theta \right ) \end {aligned} \)
\(v_{perigee}\) (closest)	\(\begin {aligned} v_p &= \sqrt {\frac {\mu }{a}\left ( \frac {1+e}{1-e}\right ) }\\ &= \sqrt {\frac {\mu }{p}}(1+e)\\ &= \sqrt {\mu \left (\frac {2}{r_p} - \frac {1}{a} \right )} \end {aligned} \)
\(v_{apogee}\) (furthest)	\( \begin {aligned} v_{a} &= \sqrt {\frac {\mu }{a}\left ( \frac {1-e}{1+e}\right ) }\\ &= \sqrt {\frac {\mu }{p}}(1-e) \\ &= \sqrt {\mu \left (\frac {2}{r_a} - \frac {1}{a} \right )} \end {aligned} \)
magnitude of \(\vec {r}\)	\(\begin {aligned} r &= \frac {a\left ( 1-e^{2}\right ) }{1+e\cos \theta }=\frac {h^{2}}{\mu }\frac {1}{1+e\cos \theta }\\ r\cos \theta &= a\left ( \cos E-e\right ) \\ r &= a\left ( 1-e\cos E\right ) \quad \text {(eq 4.2-14 Bate book)} \end {aligned} \)
period \(T\)	\(T = \frac {2}{h}\pi ab=2\pi \sqrt {\frac {a^{3}}{\mu }}\)
mean satellite speed \(n\)	\(n=\frac {2\pi }{T}=\sqrt {\frac {\mu }{a^{3}}}\)
eccentric anomaly \(E\)	\(\tan \frac {\theta }{2}=\sqrt {\frac {1+e}{1-e}}\tan \frac {E}{2}\)
area sweep rate	\(\frac {dA}{dt}=\frac {h}{2}\)
equation of motion	\(\ddot {\vec {r}}+\frac {\mu }{r^{3}}\vec {r}=0\)
spherical coordinates relation	\(\cos (i) = \sin (A_z) \cos (\phi )\) where \(i\) is the inclination and \(A_z\) is the azimuth and \(\phi \) is latitude 1

Notice in the above, that the period \(T\) of satellite depends only on \(a\) (for same \(\mu \))

In the above, \(\mu =GM\) where \(M\) is the mass of the body at the focus of the ellipse and \(G\) is the gravitational constant. \(h\) is the specific mass angular momentum (moment of linear momentum) of the satellite. Hence the units of \(\frac {h^{2}}{\mu }\) is length.

To draw the locus of the satellite (the small body moving around the ellipse, all what we need is the eccentricity \(e\) and \(a\), the major axes length. Then by changing the angle \(\theta \) the path of the satellite is drawn. I have a demo on this here

See http://nssdc.gsfc.nasa.gov/planetary/factsheet/earthfact.html for earth facts

This table below is from my class EMA 550 handouts (astrodynamics, spring 2014)

Figure 3.2: Astrodynamics constants

3.3 Flight path angle for ellipse \(\gamma \)

Figure 3.3: Flight path angle

\begin {align*} v & =\left \vert \boldsymbol {\vec {v}}\right \vert =\sqrt {\left \vert \boldsymbol {\vec {v}}_{r}\right \vert ^{2}+\left \vert \boldsymbol {\vec {v}}_{n}\right \vert ^{2}}\\ & =\sqrt {\mu \left ( \frac {2}{r}-\frac {1}{a}\right ) } \end {align*}

To find \(\gamma \), if \(r\) is given then use

\[ \cos \gamma =\sqrt {\frac {ap}{r\left ( 2a-r\right ) }}=\sqrt {\frac {a^{2}\left ( 1-e^{2}\right ) }{r\left ( 2a-r\right ) }} \]

If \(\theta \) is given, then use

\[ \tan \gamma =\frac {e\sin \theta }{1+e\cos \theta } \]

3.4 Parabolic trajectory

This diagram shows the parabolic trajectory

Figure 3.4: parabolic trajectory

3.5 Hyperbolic trajectory

This diagram shows the hyperbolic trajectory

Figure 3.5: hyperbolic trajectory

This diagram below from Orbital mechanics for Engineering students, second edition, by Howard D. Curtis, page 109

Figure 3.6: diagram below from Orbital mechanics textbook

3.6 showing that energy is constant

Showing that energy \(\mathcal {E}=\frac {v^{2}}{2}-\frac {\mu }{r}\) is constant.

Most of such relations starts from the same place. The equation of motion of satellite under the assumption that its mass is much smaller than the mass of the large body (say earth) it is rotating around. Hence we can use \(\nu =GM \) and the equation of motion reduces to \[ \ddot {\vec {r}} + \frac {\mu }{r^{3}} \vec {r} = 0 \] In the above equation, the vector \(\vec {r}\) is the relative vector from the center of the earth to the center of the satellite. The reason the center of earth is used as the origin of the inertial frame of reference is due to the assumption that \(M\gg {m}\) where \(M\) is the mass of earth (or the body at the focal of the ellipse) and \(m\) is the mass of the satellite. Hence the median center of mass between the earth and the satellite is taken to be the center of earth. This is an approximation, but a very good approximation.

The first step is to dot product the above equation with \(\dot {\vec {r}}\) giving

\begin {align} \dot {\vec {r}} \cdot \ddot {\vec {r}} + \dot {\vec {r}} \cdot \frac {\mu }{r^{3}} \vec {r} & = 0 \tag {1} \end {align}

And there is the main trick. We look ahead and see that \(\dot {\vec {r}}\cdot \ddot {\vec {r}}=\dot {r}\ddot {r}\) but \(\dot {r}\ddot {r}=\frac {d}{dt}\left ( \frac {\dot {r}^{2}}{2}\right ) \) and we also see that \(\dot {\vec {r}}\cdot \frac {\mu }{r^{3}}\vec {r}=\mu \frac {\dot {r}}{r^{2}}\) but \(\mu \frac {\dot {r}}{r^{2}}=\frac {d}{dt}\left ( \frac {-\mu }{r}\right ) \) Hence equation 1 above can be written as \[ \frac {d}{dt}\left ( \frac {v^{2}}{2}-\frac {r}{\mu }\right ) =0 \] Hence \[ \mathcal {E}=\frac {v^{2}}{2}-\frac {r}{\mu } \] Where \(\mathcal {E}\) is a constant, which is the total energy of the satellite.

3.7 Earth satellite Transfer orbits

3.7.1 Hohmann transfer

This diagram shows the Hohmann transfer

Figure 3.7: Hohmann transfer

3.7.2 Bi-Elliptic transfer orbit

Figure 3.8: Bi-Elliptic transfer orbit

3.7.3 semi-tangential elliptical transfer

Figure 3.9: semi-tangential elliptical transfer

3.8 Rocket engines, Hohmann transfer, plane change at equator

two cases: Hohmann transfer, 2 burns, or semi-tangential. All burns at equator.

Figure 3.10: Hohmann transfer cases

3.8.0.1 Rocket equation with plane change not at equator

Figure 3.11: Rocket equation with plane change not at equator

3.9 Spherical coordinates

Figure 3.12: Spherical coordinates

3.10 interplanetary transfer orbits

3.10.0.1 interplanetary hohmann transfer orbit, case one

Figure 3.13: interplanetary hohmann transfer orbit, case one

The following are the steps to accomplish the above. The first stage is getting into the Hohmann orbit from planet 1, then reaching the sphere of influence of the second planet. Then we either do a fly-by or do a parking orbit around the second planet. These steps below show how to reach the second planet and do a parking orbit around it.

The input is the following.

1. \(\mu _{1}\) planet one standard gravitational parameter
2. \(\mu _{2}\) planet two standard gravitational parameter
3. \(\mu _{sun}\) standard gravitational parameter for the sun \(1.327\times 10^{8}\) km
4. \(r_{1}\) planet one radius
5. \(r_{2}\) planet two radius
6. \(alt_{1}\) original satellite altitude above planet one. For example, for LEO use 300 km
7. \(alt_{2}\) satellite altitude above second planet. (since goal is to send satellite for circular orbit around second planet)
8. \(R_{1}\) mean distance of center of first planet from the sun. For earth use \(AU=1.495978\times 10^{8}\) km
9. \(R_{2}\) mean distance of center of second planet from the sun. For Mars use \(1.524\) \(AU\)
10. \(SOI_{1}\) sphere of influence for first planet. For earth use \(9.24\times 10^{8}\) km
11. \(SOI_{2}\) sphere of influence for second planet.

Given the above input, there are the steps to achieve the above maneuver

1. Find the burn out distance of the satellite \(r_{bo}=r_{1}+alt_{1}\)
2. Find satellite speed around planet earth (relative to planet) \(V_{sat}=\sqrt {\frac {\mu _{1}}{r_{bo}}}\)
3. Find Hohmann ellipse \(a=\frac {R_{1}+R_{2}}{2}\)
4. Find speed of satellite at perigee relative to sun \(V_{perigee}=\sqrt {\mu _{sun}\left ( \frac {2}{R_{1}}-\frac {1}{a}\right ) }\)
5. Find speed of earth (first planet) relative to sun \(V_{1}=\sqrt {\frac {\mu _{sun}}{R_{1}}}\)
6. Find escape velocity from first planet \(V_{\infty ,out}=V_{perigee}-V_{1}\)
7. Find burn out speed at first planet by solving the energy equation \(\frac {V_{bo}^{2}}{2}-\frac {\mu _{1}}{r_{bo}}=\frac {V_{\infty ,out}^{2}}{2}-\frac {\mu _{1}}{SOI_{1}}\)for \(V_{bo}\)
8. Find \(\Delta V_{1}\) needed at planet one \(\Delta V_{1}=V_{bo}-V_{sat}\)
9. Find \(e\) the eccentricity of the escape hyperbola \(e=\sqrt {1+\frac {V_{\infty }^{2}V_{bo}^{2}r_{bo}^{2}}{\mu _{1}^{2}}}\)
10. Find the angle with the path of planet one velocity vector \(\eta =\arccos \left ( -\frac {1}{e}\right ) \)
11. Find the dusk-line angle \(\theta =180^{0}-\eta \)

The above completes the first stage, now the satellite is in the Hohmann transfer orbit. Assuming it reached the orbit of the second planet ahead of it as shown in the diagram above. Now we start the second stage to land the satellite on a parking orbit around the second planet at altitude \(alt_{2}\) above the surface of the second planet. These are the steps needed.

1. Find the apogee speed of the satellite \(V_{apogree}=\sqrt {\mu _{sun}\left ( \frac {2}{R_{2}}-\frac {1}{a}\right ) }\)
2. Find speed of second planet \(V_{2}=\sqrt {\frac {\mu _{sun}}{R_{2}}}\)
3. Find \(V_{\infty }\) entering the second planet sphere of influence \(V_{\infty ,in}=V_{2}-V_{apogree}\)
4. Find burn in radius where the satellite will be closest to the second planet. \(r_{bo}=r_{1}+alt_{2}\)
5. Find burn out speed at second planet by solving the energy equation \(\frac {V_{bo}^{2}}{2}-\frac {\mu _{2}}{r_{bo}}=\frac {V_{\infty ,in}^{2}}{2}-\frac {\mu _{2}}{SOI_{2}}\)for \(V_{bo}\)
6. Find impact parameter \(b\) on entry to second planet SOI \(b=\frac {r_{bo}V_{bo}}{V_{\infty ,in}}\)
7. Find the required satellite speed around the second planet \(V_{sat}=\sqrt {\frac {\mu _{2}}{r_{bo}}}\)
8. Find \(\Delta V_{2}\) needed at planet two \(\Delta V_{2}=V_{sat}-V_{bo}\)
9. Find \(e\) the eccentricity of the approaching hyperbola on second planet \(e=\sqrt {1+\frac {V_{\infty }^{2}V_{bo}^{2}r_{bo}^{2}}{\mu _{2}^{2}}}\)
10. Find the angle with the path of planet two velocity vector \(\eta =\arccos \left ( -\frac {1}{e}\right ) \)
11. Find the dusk-line angle for second planet \(\theta =2\eta -180^{0}\)

3.10.1 rendezvous orbits

3.10.1.1 Two satellite, walking rendezvous using Hohmann transfer

Figure 3.14: Two satellite, walking rendezvous using Hohmann transfer

Algorithm 1 Hohmann Walking Rendezvous Orbit, case 1

1:   function hohmann_walking_rendezvous(\(\theta _0,r,altitude,\mu \))
2:    \(\theta _0 := \theta _0 \frac {\pi }{180}\) \(\triangleright \) convert from degrees to radian
3:    \(r_a := r+\text {altitude}\)
4:    \(T := 2\pi \sqrt {\frac {r_a^3}{\mu }}\) \(\triangleright \) period of circular orbit
5:    \(N:=1\)
6:    done:=false
7:    while not(done) do
8:    TOF := \(\left (N-\frac {\theta _0}{2\pi }\right )T\)
9:    \(a := \textbf {solve}\, \left (TOF = 2\pi \sqrt {\frac {a^3}{\mu }}\right ) \text {for}\, a\)
10:    \(r_p:=2 a-r_a\)
11:    if \(rp<r\) then
12:    \(N=N+1\)
13:    else
14:    done:=true
15:    end if
16:    end while
17:    \(V_{befor}:=\sqrt {\frac {\mu }{h}}\)
18:    \(V_{after}:=\sqrt {\mu \left ( \frac {2}{h} - \frac {1}{a} \right )}\)
19:    \(\Delta V:= 2(V_{after}-V_{before})\)
20:    return \((\text {TOF},\Delta V)\)
21:   end function

An example implementation is below

hohmannRendezvousSameOrbit[\[Theta]00_, r_, alt_, mu_] := 
Module[{\[Theta]0 = \[Theta]00*Pi/180, n = 1, delT, v1, v2, period, a, 
rp, ra, done = False, vBefore, vAfter}, 
ra     = r + alt; 
period = 2 Pi Sqrt[ra^3/mu]; 
While[Not[done], 
delT = (n - \[Theta]0 /(2 Pi)) period ; 
a = First@Select[a /. NSolve[delT == 2 Pi Sqrt[a^3/mu], a], 
Element[#, Reals] &]; 
rp = 2 a - ra; 
If[rp < r,(*we hit the earth, try again*) 
n = n + 1, 
done = True 
] 
]; 
vBefore = Sqrt[mu/h]; 
vAfter  = Sqrt[mu (2/h - 1/a)]; 
{delT, 2 (vAfter - vBefore)}  (*return value*) 
]
 

And calling the above

mu = 324859; 
alt = 1475.776; 
r = 6052; 
\[Theta]0 = 3.80562; (*degree*) 
hohmannRendezvousSameOrbit[\[Theta]0, r, alt, mu]
 

gives

{7123.89, -0.0467913}
 

3.10.1.2 Two satellite, separate orbits, rendezvous using Hohmann transfer, coplaner

Figure 3.15: Two satellite, separate orbits, rendezvous using Hohmann transfer, coplaner

Algorithm 2 Hohmann rendezvous algorithm, case 1

1:   function hohmann_rendezvous_1(\(\theta _0,r_a,r_b,\mu \))
2:    \(\theta _0 := \theta _0 \frac {\pi }{180}\) \(\triangleright \) convert from degrees to radian
3:    \(a := \frac {r_a+r_b}{2}\) \(\triangleright \) Hohmann orbit semi-major axes
4:    \(\text {TOF} := \pi \sqrt {\frac {a^3}{\mu }}\) \(\triangleright \) time of flight on Hohmann orbit
5:    \(\theta _H := \pi \left ( 1- \left ( \frac {r_a+r_b}{2 r_b} \right )^{3/2} \right )\) \(\triangleright \) required phase angle before starting Hohmann transfer
6:    \(\omega _a := \sqrt {\frac {\mu }{r_a^3}}\) \(\triangleright \) angular speed of lower rad/sec
7:    \(\omega _b := \sqrt {\frac {\mu }{r_b^3}}\) \(\triangleright \) angular speed of higher satellite rad/sec
8:    if \(\theta _0\leq \theta _H\) then \(\triangleright \) adjust initial angle if needed
9:    \(\theta _0 := \theta _0+2 \pi \)
10:    end if
11:    \(\text {wait\_time} := \frac {\theta _0-\theta _H}{\omega _a-\omega _b}\) \(\triangleright \) how long to wait before starting Hohmann transfer
12:    \(\text {wait\_time} := \text {wait\_time} + \text {TOF}\) \(\triangleright \) now ready to go, add Hohmann transfer time
13:    return wait_time
14:   end function

An example implementation is below (in Maple)

hohmann_rendezvous_1:= proc({ 
theta::numeric:=0, 
r1::numeric:=0, 
r2::numeric:=0, 
mu::numeric:=3.986*10^5}) 
local theta0,thetaH,TOF,a,omega1,omega2,wait_time; 
theta0  := evalf(theta*Pi/180); 
a       := (r1+r2)/2; 
TOF     := Pi*(sqrt(a^3/mu)); 
omega1  := sqrt(mu/r1^3); 
omega2  := sqrt(mu/r2^3); 
thetaH  := evalf(Pi*(1-((r1+r2)/(2*r2))^(3/2))); 
if theta0 <= thetaH then 
theta0 := theta0+2*Pi; 
fi; 
wait_time := TOF+(theta0-thetaH)/(omega1-omega2); 
eval(wait_time); 
end proc:
 

And calling the above for two different cases gives (times in hrs)

TOF:=hohmann_rendezvous_1(r1=6678,r2=6878,theta=0): 
evalf(TOF/(60*60)); 
35.23480353
 

And

TOF:=hohmann_rendezvous_1(r1=6678,r2=6878,theta=280): 
evalf(TOF/(60*60)); 
27.49212919
 

3.10.1.3 Two satellite, separate orbits, rendezvous using bi-elliptic transfer, coplaner

Figure 3.16: Two satellite, separate orbits, rendezvous using bi-elliptic transfer, coplaner

In this transfer, the lower (fast satellite) does not have to wait for phase lock as in the case with Hohmann transfer. The transfer can starts immediately. There is a free parameter \(N\) that one select depending on fuel cost requiments or any limitiation on the first transfer orbit semi-major axes distance required. One can start with \(N=0\) and adjust as needed.

Algorithm 3 Hohmann rendezvous algorithm, case 2

1:   function hohmann_rendezvous_2(\(\theta _0,r_1,r_2,N,\mu \))
2:    \(\theta _0 := \theta _0 \frac {\pi }{180}\) \(\triangleright \) convert from degrees to radian
3:    \(\theta _H := \pi \left ( 1- \left ( \frac {r_1+r_2}{2 r_2} \right )^{3/2} \right )\) \(\triangleright \) Find Hohmann ideal phase angle before transfer
4:    if \(\theta _0 = \theta _H\, \textbf {and}\, N=0\) then \(\triangleright \) adjust for special case
5:    \(a := \frac {r_2+r_1}{2}\)
6:    \(\text {TOF} := \pi \left ( \sqrt {\frac {a^3}{\mu } } \right )\)
7:    else
8:    \(\omega _2 := \sqrt {\frac {\mu }{r_2^3}}\) \(\triangleright \) angular speed of slower satellite in rad/sec
9:    \(t_2 :=\frac {(2 \pi -\theta _0)+ 2\pi N}{\omega _2}\) \(\triangleright \) find time of light of the slower satellite
10:    \(a_1 := \frac {r_t+r_1}{2}\)
11:    \(a_2 := \frac {r_t+r_2}{2}\)
12:    \(\text {TOF} := \pi \left ( \sqrt {\frac {a_1^3}{\mu } + \frac {a_2^3}{\mu }} \right )\) \(\triangleright \) time of flight for the fast satellite
13:    \(r_t := \textbf {solve}\,\left (t_2 = TOF \right ) \text {for}\, r_t\) \(\triangleright \) Solve numerically for \(r_t\)
14:    end if
15:    return TOF
16:   end function

An example implementation is below in Maple

hohmann_rendezvous_2:= proc({ 
theta::numeric:=0, 
r1::numeric:=0, 
r2::numeric:=0, 
N::nonnegint:=0, 
mu::numeric:=3.986*10^5}) 
local theta0,thetaH,TOF; 
theta0  := theta*Pi/180; 
thetaH  := Pi*(1-((r1+r2)/(2*r2))^(3/2)); 
if theta0 = thetaH and N = 0 then 
proc() 
local a:=(r1+r2)/2; 
TOF:= Pi*(sqrt(a^3/mu)); 
end proc() 
else 
proc() 
local t2,a1,a2,rt,omega2; 
omega2 := sqrt(mu/r2^3); 
t2     := ((2*Pi-theta0)+2*Pi*N)/omega2; 
a1     := (rt+r1)/2; 
a2     := (rt+r2)/2; 
TOF    := Pi*(sqrt(a1^3/mu)+sqrt(a2^3/mu)); 
rt     := op(select(is, [solve(t2=TOF,rt)], real)); 
end proc() 
fi; 
eval(TOF); 
end proc:
 

And calling the above for two different cases gives

TOF:=hohmann_rendezvous_2(theta=0,r1=6678,r2=6878,N=0): 
evalf(TOF/(60*60)); #in hrs 
1.576892101 
TOF:=hohmann_rendezvous_2(theta=160,r1=6678,r2=6878,N=1): 
evalf(TOF/(60*60));  #in hrs 
2.452943266
 

3.10.2 Semi-tangential transfers, elliptical, parabolic and hyperbolic

Figure 3.17: Semi-tangential transfers, elliptical, parabolic and hyperbolic

Figure 3.18: Semi-tangential transfers, elliptical, parabolic and hyperbolic (2)

3.10.3 Lagrange points

Figure 3.19: Lagrange points

3.10.4 Orbit changing by low contiuous thrust

Figure 3.20: Orbit changing by low contiuous thrust

3.10.5 References

1. Orbital mechanics for Engineering students, second edition, by Howard D. Curtis
2. Orbital Mechanics, Vladimir Chobotov, second edition, AIAA
3. Fundamentals of Astrodynamics, Bates, Muller and White. Dover1971

¹Image is from my class notes, page 7-9, EMA 550, Univ. Of Wisconsin, Madison by professor S. Sandrik