Falcon 9 descent attempt

#! /usr/bin/python
import numpy as np
from scipy.integrate import odeint

# Approximating the powered descent of a falcon 9 with a single Merlin 1D over the last 30 seconds of flight.
# Also an attempt to learn odeint, and how to use it to deal with a non-linear 2nd order ODE.

rho = 1.2   # sealevel atmo density, ish
d = 3.7     # Falcon 9 diamter
Ve = 2730.0 # Effective exhaust velocity
F = 934e3   # Thrust
m0 = F/25.0 # final mass is unclear, so assumed to allow an acceleration of 25 m/s/s

y0 = [0.0, 0.0] # y = 0, y' = 0
t = np.arange(0.0, 30.0, 0.1)

# time is running backwards, hence the increasing mass and upwards drag
def rocket(y, t, F, Ve, m0, rho, d): 
    alt, spd = y 
    rho1 = np.exp(-alt/7000.0)              # air density altitude
    m1 = m0 + t*F/(Ve*9.8)                  # Mass at any given time
    dvdt = [spd, F/m1 + 5.26*rho1*d*spd*spd/m1 - 9.8]   # Engine thrust, drag, gravity. Blame ProjectThoth for the odd drag formulation.
    return dvdt

solution = odeint(rocket, y0, t, args=(F, Ve, m0, rho, d)) 

print(solution)