Ecuación para la velocidad y la distancia desde el Sol de una nave espacial impulsada por velas solares

Las matemáticas se vuelven muy complicadas muy rápidamente. Antes incluso de considerar el empuje, recuerde que la solución de la órbita de dos cuerpos de Kepler en sí misma no es una función explícita del tiempo: el radio se da como una función del ángulo (anomalía verdadera), pero el tiempo es proporcional a la anomalía media, que debe convertirse en excéntrica. anomalía y luego anomalía verdadera por solución numérica de la ecuación trascendental de Kepler .

Incluso para los casos más simples, estás hablando de series o expansiones en fracciones continuas de integrales elípticas. Sin embargo, si tiene el libro de Battin sobre astrodinámica, él muestra cómo hacerlo a mano, de la forma en que lo hicieron personas como Gauss y Legendre, y que se esperaba que todos sus estudiantes en el MIT lucharan, en las páginas 408 a 418. Dicho esto, comienza la sección con

Una de las posibles aceleraciones perturbadoras que afectan a una nave espacial es la aceleración de empuje producida por los motores del vehículo. La determinación de la trayectoria en estas circunstancias, como en casi todos los problemas de movimiento perturbado, generalmente requiere la aplicación de técnicas numéricas especiales. Sin embargo, hay varios ejemplos de algún interés práctico para los cuales se puede hacer un análisis considerable. Específicamente, si la aceleración del empuje es de magnitud constante y está dirigida radial, tangencial o circunferencialmente, entonces es posible obtener, al menos parcialmente, algunos resultados matemáticos bastante interesantes.

Una encuesta de algunos métodos más, que tienen en cuenta exactamente cómo planeas salir en espiral, que en general encuentro más accesible y útil en este caso, es Petropoulos y Sims 2002 , "A Review of Some Exact Solutions to the Planar Equations de movimiento de una nave espacial de empuje". Tenga en cuenta que están más interesados ​​en el uso de paneles solares para alimentar un propulsor de iones, lo que le permite empujar en la dirección más ventajosa, ¡que está muy lejos de ser radial!

Solo respuesta parcial, no sé qué estoy haciendo exactamente, pero lo intentaré.

Actualización: los cohetes funcionan porque sus fallas espectaculares se analizaron a la luz del día en lugar de ocultarse. Si bien este comentario expresa la preocupación de que las matemáticas a continuación puedan confundir negativamente a los lectores, confío en que para los lectores aquí quedará claro que no da la respuesta correcta y que la técnica de comparar una solución matemática con una solución numérica sirve como un buen ejemplo de cómo se puede verificar una respuesta.

Matemáticas

Dada la aceleración inicial a 1 AU es$8.17 \times 10^{-4} m/s^2$.

Inclínelo a 45 grados para que el empuje sea tangencial, divídalo por$\sqrt{2}$ya que ahora es oblicuo al Sol, y tenga en cuenta la caída con la distancia desde el Sol:

en la dirección progresiva (la misma dirección que la velocidad actual).

ecuación vis-viva :

Ahora

dónde$GM$es el parámetro gravitatorio estándar del Sol $1.327 \times 10^{20} m^3/s^2$. El signo menos aparece porque sabemos que, contrariamente al primer instinto, cuando tenemos una fuerza de aceleración en la dirección prograda, desaceleramos de manera contraria a la intuición en la misma cantidad. Esto también se cita en varias otras publicaciones aquí, buscaré otras respuestas para citar ...

Reescribe y resuelve:

si establecemos$t(v_0) = 0$en otras palabras, en el tiempo cero nos estamos moviendo a velocidad orbital$v_0 = \sqrt{GM/AU}$ obtenemos 0.408 años.

Numéricamente

Curiosamente, cuando intento simular lo mismo numéricamente, ¡ obtengo un tiempo de escape de 0,515 años! Esto me sorprende por dos razones:

1. La velocidad nunca llega a cero, lo que significa que las matemáticas anteriores son incorrectas.
2. ¡Y sin embargo está bastante cerca!

Conclusión

Tardará aproximadamente medio año en escapar, pero aún no tengo la ecuación completa.

Otras lecturas:

• ¿Es posible llegar al Sol sin gastar combustible/masa de reacción?
• ¿Cuál es el ángulo óptimo para una vela solar desorbitada hacia el Sol cuando se incluye el empuje radial? (actualmente sin respuesta)

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

AU = 150E+06 * 1000   # meters
GM = 1.327E+20   # m^3/s^2
a0 = 8.17E-04 / np.sqrt(2)  # dv/dt at 1 AU
year = 365.2564 * 24 * 3600 

coef = ( GM**2 / (4 * a0 * AU**2) )

v0 = np.sqrt(GM/AU)
print('v0: ', v0)
print('coef: ', coef)

t0 = coef / v0**3
print('t0 / year: ', t0 / year)

def deriv(X, t):
    x, v = X.reshape(2, -1)
    vhat = v / np.sqrt((v**2).sum())
    rsq = (x**2).sum()
    acc_thrust = vhat * a0 * rsq / AU**2
    acc = acc_thrust - GM * x * rsq**-1.5
    return np.hstack((v, acc))

X0 = np.array([AU, 0, 0, v0])

time = np.linspace(0, 0.6, 10001) * year  # half year

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

print(answer.shape)

x, y, vx, vy = answer.T
r, speed = [np.sqrt((thing**2).sum(axis=0))
            for thing in answer.T.reshape(2, 2, -1)]

E = 0.5 * speed**2 - GM/r
E_norm = np.abs(E[0])

i_esc = np.argmax(E>=0)
things = time, r, x, y, speed, E
t_esc, r_esc, x_esc, y_esc, s_esc, E_esc = [thing[i_esc]
                                            for thing in things]

print('t_esc / year: ', t_esc / year)

if True:
    plt.figure()
    plt.subplot(2, 2, 1)
    plt.plot(time/year, r/AU)
    plt.plot([t_esc/year], [r_esc/AU], 'ok')
    plt.ylabel('r/AU')
    plt.xlabel('time (years)')
    
    plt.subplot(2, 2, 2)
    plt.plot(time/year, speed/1000)
    plt.plot([t_esc/year], [s_esc/1000], 'ok')
    plt.ylabel('speed (km/s)')
    plt.xlabel('time (years)')

    plt.subplot(2, 2, 3)
    plt.plot(time/year, E/E_norm)
    plt.plot([t_esc/year], [E_esc/E_norm], 'ok')
    plt.plot(time/year, np.zeros_like(time), '-k')
    plt.ylabel('Energy (norm)')
    plt.xlabel('time (years)')

    plt.subplot(2, 2, 4)
    plt.plot(x/AU, y/AU)
    plt.plot([x_esc/AU], [y_esc/AU], 'ok')
    plt.plot([0], [0], 'oy')
    th = np.linspace(0, 2*np.pi, 201)
    plt.plot(np.cos(th), np.sin(th), '-r', linewidth=0.5)
    plt.ylim(-1, 1.5)
    plt.gca().set_aspect('equal')
    plt.xlabel('AU')
    plt.ylabel('AU')
    plt.show()