Low-thrust spiraling to escape, is the flight path angle (gamma) at C3=0 always 39 degrees?

@MarkAdler's comment led me to ask Why would a slow spiral from a C3 of zero take about 2.4 times as much ΔV as an impulsive maneuver? which resulted in this tidy and efficient @MarkAdler answer which points to another thoughtful answer about slowly spiraling out from a circular orbit to escape in the limit of very weak prograde propulsion, which (at first counterintuitively) slows you down while raising your orbit.

Below that answer is yet another easter-egg-like comment gem.

Always aligned with the velocity vector. That is the most efficient use of the thrust in order to increase the specific energy. The final γ is 31°.

In this answer @Julio provides a diagram showing definitions for both $\beta$ and $\gamma$ angles which measure the angle between the instantaneous velocity vector and the radial and the tangential directions, respectively.

In this answer @TomSpilker elaborates on these angles, and in this answer I give a little more information on how to calculate them.

Now I've gone back and calculated an outwardly spiraling orbit under low thrust using various conditions. Invariably I end up with a final angle $\gamma$ (gamma) of about 39 degrees when checking the moment where C3 = 0, not 31 degrees.

I'm doing a unitless calculation where GM = 1.0 and the period of an r=1.0 orbit is $2 \pi$. In this case C3 = v^2 - 2/r.

note: For this calculation, thrust is always in the same direction as velocity $\mathbf{v}$, rather than in the tangential direction (perpendicular to $\mathbf{r}$) and I'm beginning to wonder if herein lies the difference between 31 and 39 degrees.

Question: Is this ~39 degrees at C3=0 correct, and is it expected to be invariant like this?

      starting conditions                              at C3 = 0
-------------------------------     ------------------------------------------
rstart  vstart    C3    thrust      time   delta-v  gamma(deg)    r       v        C3
1.0     1.0    -1.0    0.01        74.5    0.745     38.9       8.78    0.477   0.000
1.0     1.0    -1.0    0.001       856.3   0.856     39.2      27.80    0.268   0.000
1.0     1.0    -1.0    0.0001      9192.1  0.919     39.2      87.91    0.151   0.000
4.0     0.5    -0.25   0.0001      4192.1  0.419     39.1      87.90    0.151   0.000


def deriv(X, t):
x, v  = X.reshape(2, -1)
vnorm = v / np.sqrt((v**2).sum())
acc_g = -x * ((x**2).sum())**-1.5
acc_t = thrust * vnorm
return np.hstack((v, acc_g + acc_t))

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint as ODEint

halfpi, pi, twopi = [f*np.pi for f in (0.5, 1, 2)]

T    = 16 * twopi        # or 160, 1600

ntot = 20001
time = np.linspace(0, T, ntot)

rstart = 1.0             # or 4.0
vstart = np.sqrt(1./rstart)

X0     = np.array([rstart, 0, 0, vstart])

thrust = 0.01            # or 0.001, 0.0001

answer, info = ODEint(deriv, X0, time, full_output= True)

xx, vv = answer.T.reshape(2, 2, -1)

r   = np.sqrt((xx**2).sum(axis=0))
vsq =         (vv**2).sum(axis=0)
C3 = vsq - 2./r

nstop = np.argmax(C3>0) + 1

dotted     = (xx*vv).sum(axis=0)
rabs, vabs = [np.sqrt((thing**2).sum(axis=0)) for thing in (xx, vv)]
gamma      = np.arcsin(dotted/(rabs*vabs))   # Per Tom Spilker's answer Eq. 3

print 'C3 min, max: ', C3.min(), C3.max()
print 'nstop, ntot: ', nstop, ntot
if True:
plt.figure()

plt.subplot(1, 2, 1)
plt.plot(xx[0, :nstop], xx[1, :nstop])

plt.subplot(3, 2, 2)
plt.plot(time[:nstop], r[:nstop])
plt.ylabel('r')

plt.subplot(3, 2, 4)
plt.plot(time[:nstop], C3[:nstop])
plt.plot(time[:nstop], np.zeros_like(C3)[:nstop], '-k')
plt.ylabel('C3')

plt.subplot(3, 2, 6)
plt.plot(time[:nstop], degs*gamma[:nstop])
plt.ylabel('gamma (deg)')

plt.suptitle('thrust = 0.0001, start at r=4, time=4192.1, gamma=39.12 deg, r=87.90', fontsize=16)

plt.show()


Sorry, must have been a typo in the comment. I went back to the original notebook in which I made the plots, and indeed the final $$\gamma$$ for the 0.001 acceleration case was 39.2°

It is not always 39.2°, but it does go to that asymptotically as the acceleration becomes smaller. Here is a plot of the $$\gamma$$ in degrees at $$C_3=0$$ as a function of the relative acceleration:

I am not aware of a way to determine that $$\gamma$$ analytically.

Below is the same plot for when accelerating tangentially, as opposed to in the velocity direction. It looks identical, except for the y-axis, where here it converges to 32.3°.

Though you wouldn't do that, since accelerating in the velocity direction is more efficient.

• This is really beautiful! Thanks for the snazzy plot. It must be nice to own "grown up" software. It turns out that purely tangential acceleration (perpendicular to the instantaneous radius vector) wiggles around 32 degrees while acceleration parallel to the instantaneous velocity wiggles around 39 degrees (at C3=0). Is there any chance you have enough remaining electrons to run that as well?
– uhoh
Jul 17 '18 at 5:21
• @uhoh that looks like a Mathematica plot, if interested there ARE open source variations on Mathematica's software. Though, admittedly, the paid Mathematica version is amazing. Jul 17 '18 at 18:47
• Yep. Mathematica. Easily worth the \$\$. (Not a shill.) Jul 17 '18 at 18:59
• fyi I've just asked Ratio of low-thrust slow spiral to Hohmann transfer delta-v? I thought there was something about that around already, but can't find it.
– uhoh
Feb 10 '19 at 4:34
• – uhoh
Oct 8 '19 at 4:19