Rocket Landing Guidance via Convex MPC

3-DOF simulation of powered-descent guidance for a reusable booster. A second-order cone program (G-FOLD lineage: log-mass substitution + lossless convexification, plus successive convexification for drag) is used two ways: solved once and played back open-loop, and re-solved from the true state every 0.5 s as an MPC. The landing burn is also chained into a full return profile (post-MECO coast, boostback burn, coast, entry burn, coast, landing burn), flown from the real separation state of a Falcon 9 RTLS mission and compared against that mission's public telemetry.

Everything is written from scratch: dynamics, atmosphere, solver setup, plots and the web page. No flight software, no proprietary data. See NOTICE.md.

Python 3.11+, CVXPY + Clarabel, pure Python/NumPy otherwise.

Run

python -m pip install -e ".[dev]"
pytest                                     # 11 tests, seconds

python src/closed_loop_sim.py              # open-loop vs MPC, single landing burn
python src/plot_results.py                 # -> results/comparison.png

python src/full_profile_sim.py             # full return profile
python src/plot_full_profile.py            # -> results/full_profile.png

python src/validate_drag_model.py          # drag model vs real coast telemetry
python src/drag_holdout_validation.py      # leave-one-flight-out variant
python src/monte_carlo.py --trials 20      # seeded dispersions

python src/export_web_data.py && python web/build.py   # -> web/index.html
python src/plot_3d_playback.py --scenario full         # -> results/*.mp4

Scenario constants live in src/params.py. The full-profile block holds real CRS-11/CRS-12 telemetry values, so changing them means deliberately diverging from the real flight.

results/ holds the generated artifacts of a full run, checked in so the numbers below are inspectable without re-running anything.

Model

• Dynamics (dynamics.py): point mass; states are position, velocity, mass. Thrust, constant gravity, aerodynamic drag. RK4.
• Atmosphere (atmosphere.py): US Standard Atmosphere 1976 density, three-segment Mach-dependent Cd (flat subsonic, linear transonic rise, flat hypersonic), reference area 10.75 m2 (pi/4 * 3.7^2).
• Guidance (gfold_mpc.py): SOCP with the z = ln(m) mass substitution (exact, not an approximation), the ||u|| <= Gamma lossless relaxation, throttle bounds, thrust-pointing cone, glideslope cone, velocity bound, and objective min sum(Gamma)*dt. Drag is folded in by successive convexification: each pass evaluates drag along the previous pass's trajectory, where it is a per-step constant. The iteration is damped and seeded from an unpowered drag propagation, not from a vacuum solution.
• MPC (closed_loop_sim.py): shrinking horizon, re-solved every 0.5 s from the true state. If the fixed-time terminal equality becomes infeasible near touchdown, a fallback solve keeps every state and actuator constraint hard and softens only terminal position/velocity; its use is recorded in the diagnostics.
• Flight-time search (time_search.py): ascending scan for the earliest feasible duration, plus a margin. Feasibility is not monotonic in flight time, so bisection is not valid here (see Findings).
• Deliberate simplification: the physical thrust bounds T_min/m(t) <= Gamma(t) <= T_max/m(t) are frozen at the current measured mass for a whole horizon instead of being relaxed along a reference mass trajectory as in the original papers. Inside an MPC that re-solves every 0.5 s, the error self-corrects at the next tick.

Diagnostics written with every run: SCvx fixed-point residual, SOC relaxation gap, terminal residuals, physical thrust range, dry-mass margin, prediction-vs-truth burn endpoint error.

Results

Standalone landing burn, identical 1.5-5.5 s wind gust applied to the truth model only, never to either controller's internal prediction:

	miss distance	residual speed	propellant
Open-loop playback	174.6 m	8.03 m/s	2367.4 kg
Closed-loop MPC	0.001 m	0.004 m/s	2379.0 kg

Feedback buys disturbance rejection, not a cheaper trajectory (+0.5 % propellant). The throttle profile comes out max -> min -> max ("hoverslam"), which is a property of the fuel-optimal solution, not a tuned input.

Monte Carlo: a seeded three-trial full-budget run (initial state, mass, drag and wind dispersions) gave 0/3 open-loop and 3/3 MPC landings against a 5 m / 2 m/s criterion. That is a pipeline check, not a statistically meaningful campaign; raise --trials for a real distribution. Samples and quantiles are in results/monte_carlo_summary.json.

Full profile vs. CRS-11

Chained from the real separation state through two burns and three unpowered arcs. Clock is separation-relative in both columns.

milestone	CRS-11 telemetry	this sim	delta
Separation (MECO)	t+0 s, 62.4 km, 1662 m/s	identical	sim input
Boostback ignition	t+16 s, 79.2 km, 1579 m/s	t+17 s, 79.7 km, 1555 m/s	+0.6 % alt, -1.5 % speed
Boostback burnout	t+66 s, 113.5 km, 588 m/s	t+65 s, 110.8 km, 582 m/s	-2.4 % alt, -1.1 % speed
Entry-burn ignition	t+227 s, 47.5 km, 1260 m/s	t+220 s, 48.6 km, 1246 m/s	+2.5 % alt, -1.2 % speed
Entry-burn cutoff	t+239 s, 35.0 km, 954 m/s	t+233 s, 35.5 km, 908 m/s	+1.4 % alt, -4.8 % speed
Landing-burn ignition	t+284 s, 4.74 km, 320 m/s	t+276 s, 4.70 km, 326.6 m/s	-0.9 % alt, +2.1 % speed
Touchdown	t+317 s, ~0 m/s	t+310 s, 0.001 m miss, 0.005 m/s	-2.2 % mission time

Landing-burn ignition speed (326.6 vs 320 m/s) is the informative number: it is an output of ~280 s of chained simulation, not a fitted input. Total propellant 32,773 kg of 37,000 kg available.

The sim reaches every event ~2 % early, consistently. Same sign as the drag validation bias below: the model falls slightly too readily.

The boostback and entry solves reuse the core dynamics constraints but drop the glideslope and thrust-pointing cones. At those phases the vehicle is flying sideways and would have to flip attitude, which a 3-DOF point mass cannot represent honestly. That is the gap the 6-DOF item under Next steps addresses.

Drag validation

Between entry-burn cutoff (~35 km, ~950 m/s) and landing-burn ignition (~4.7 km, ~320 m/s) the engines are off. The booster falls 30 km and loses 630 m/s. Gravity is known, so drag is the only unknown force over that segment: propagate the model from the real cutoff state for the real duration and compare. Two independent targets (distance fallen, final speed), one free coefficient, nothing tuned afterwards.

Rejected by that test: a single constant Cd. Cd*A ~ 10.8 m2 matches final speed to ~2 % but overshoots distance fallen by ~11 %; roughly doubling it matches distance and misses speed by ~20 %. The two metrics probe different flow regimes.

With the three-segment Cd curve, the vehicle's real base area, and Cd_hypersonic = 1.45 fitted on the coast:

	distance fallen	final speed
CRS-11	+4.5 %	-5.7 %
CRS-12	+4.8 %	-7.5 %
vacuum model, for scale	-	+335 %

That table is in-sample: the default coefficient used both missions. drag_holdout_validation.py fits on one mission and reports errors on the untouched other, then swaps. The folds fit Cd_hypersonic = 1.530 and 1.538; holdout speed RMSE 11.8-18.6 m/s, distance error +3.5-3.9 %, final-speed error -6.5 % to -8.1 %.

Cd_subsonic is left at a standard 0.9 rather than fitted; the coast is barely subsonic before the landing burn, so this data does not constrain it (sweeping 0.5-0.9 moves the fit by <0.1 %).

Solver-side caveat: each SCvx pass is optimal only for its own frozen drag estimate, so the result is a fixed-point iterate, not a certified optimum of the true nonconvex problem. The residual and whether the iteration budget converged are in every summary. The truth model applies exact nonlinear drag, so whatever the iteration misses becomes real controller model error.

Data

Real telemetry: CRS-11 (2017-06-03) and CRS-12 (2017-08-14), both Falcon 9 RTLS landings at LZ-1, stage-1 recovery telemetry included. It is a community reconstruction digitized from public launch webcasts, public domain (Unlicense), from github.com/shahar603/Telemetry-Data. Vendored in data/real_flight_telemetry/; see SOURCE.md there. CRS-11 is the plotted reference curve, i.e. one actual flight, not an average.

An earlier version of this project used point values from general spaceflight commentary instead (MECO quoted as ~80 km / ~2778 m/s). The reconstructed telemetry puts real MECO at 61-62 km and 1655-1684 m/s. The secondary sources were not imprecise, they were wrong by ~30 % in altitude.

Findings

Things that were not obvious before building it:

• Feasibility is not monotonic in flight time. Too short: the vehicle cannot decelerate in time. Too long: minimum throttle still burns propellant every second, so a lazy descent can hit the dry-mass floor before touchdown. Hence an ascending scan rather than a bisection.
• Deep throttling is a feasibility constraint, not a realism detail. At ~40 % minimum throttle the chosen boundary conditions are geometrically infeasible: minimum sustained braking stops the vehicle above the ground and it cannot idle down further to fall the rest of the way. 20 % minimum (a deep-throttle test article, not a specific engine) fixes it.
• The constraint you expect to bind, doesn't. The landing burn was infeasible at every horizon. Relaxing one constraint at a time: the glideslope cone does not bind anywhere between 5 and 30 degrees (identical results), while the thrust-pointing cone is the sole blocker - 20 degrees infeasible everywhere, 45 degrees solves. Check the sign of each knob first: a larger glideslope angle is a tighter constraint, so the first sweep "relaxed" it by tightening it.
• Drag is load-bearing, not a refinement. Same sweep: with the pointing cone at 45 degrees, the landing burn is feasible with drag and infeasible without it. The atmosphere does most of the braking, ~630 m/s of the whole recovery.
• Seeding SCvx at zero drag fails exactly when drag matters. It assumes drag is a correction to a feasible vacuum trajectory, but the vehicle arrives at 4.7 km doing ~330 m/s, more deceleration than one engine can produce. Pass 0 is infeasible and the loop aborts before drag ever enters. Seed from an unpowered drag propagation instead.
• A 17 s "do nothing" phase was worth ~20 km. The real mission coasts 16-19 s after MECO before boostback ignition. Skipping it reverses the trajectory 20 km too early and lands ~18 km past the pad, at a 14 degree elevation angle, i.e. outside the 25 degree glideslope cone before the landing burn even starts. Unpowered phases feel free because nothing is commanded; the trajectory disagrees.
• A minimum-fuel reference can leave the controller infeasible with zero disturbance. The chained hand-off state uses 75 % of available propellant against 35 % in the standalone demo. At tf_margin=1.12 the MPC went infeasible partway through touchdown with the gust switched off - not from a disturbance, but because the frozen-Gamma-bounds approximation adds a small error at every re-solve and there was no slack left after ~15 of them. The chained landing phase uses tf_margin=1.35. Fuel margin is also headroom for the controller's own approximations.
• Two entry points quietly shared output filenames. full_profile_sim.py wrote its landing phase to the same open_loop.npz / closed_loop_mpc.npz the standalone demo uses, so plot_results.py paired one run's trajectories with the other run's summary and nothing errored. Chained logs now go to landing_*.npz.
• Plant and controller need separate checks when adding a force. After wiring drag through the solver and the coast propagation, the landing phase's truth integrator was still being called without it, so the MPC was planning against drag the simulated world never applied. Nothing errored; the numbers were just wrong. Grep every integrator call site.
• The real flight's own q is ~0 from MECO through entry-burn ignition (47-50 km), so vacuum physics is a fair approximation there. q only becomes large between entry-burn cutoff and landing-burn ignition, exactly the phase the drag model has to get right.

Web visualisation

web/index.html is a single self-contained file: no CDN, no external fonts, no charting or 3D library. It contains an interactive 3D playback of the return trajectory with the real CRS-11 path as a ghost line, the milestone table, charts for altitude, velocity, dynamic pressure and drag validation, and a 1x real-time replay of the 34 s landing burn with a thrust-vector display (pointing-cone dial against the 45 degree limit, throttle-scaled plume, tilt and throttle sparklines). The telemetry console shows 18 channels computed in export_web_data.py from the simulation's own state history - thrust from the mass history via T = mdot*Isp*g0, sensed acceleration by differentiating inertial velocity and removing gravity.

web/index_template.html is the source; build.py only injects the data, so edit the template, not the generated file. ?t=<seconds> opens the page paused at a given mission time.

Caveat stated on the page: 3-DOF does not separate vehicle attitude from engine gimbal, so the displayed angle is the commanded thrust direction, which a real vehicle achieves with both together.

Layout

src/
  params.py                   vehicle / environment / guidance constants
  dynamics.py                 nonlinear truth model + RK4
  atmosphere.py               US76 density + Mach-dependent drag
  gfold_mpc.py                SOCP solves: landing MPC + boostback/entry burns
  coast.py                    unpowered propagation with altitude trigger
  time_search.py              feasible / fuel-optimal flight-time search
  closed_loop_sim.py          open-loop vs closed-loop MPC experiment
  full_profile_sim.py         chained six-phase return profile
  real_flight_reference.py    loads the vendored CRS-11/CRS-12 telemetry
  validate_drag_model.py      drag check against real coast telemetry
  drag_holdout_validation.py  leave-one-flight-out Cd fit
  monte_carlo.py              seeded dispersions
  plot_results.py             results/comparison.png
  plot_full_profile.py        results/full_profile.png
  plot_3d_playback.py         3D MP4/GIF (ffmpeg, else Pillow GIF)
  export_web_data.py          results/web_data.json
web/
  index_template.html         page source (edit this one)
  build.py                    inlines web_data.json -> web/index.html
  index.html                  generated, standalone
data/real_flight_telemetry/   vendored CRS-11/CRS-12 telemetry + SOURCE.md
results/                      generated artifacts of a full run
tests/                        physics, telemetry and diagnostic regression tests

Next steps

• 6-DOF: attitude, rotational dynamics, gimbal actuator lag, grid fins. The largest remaining gap. It would remove the boostback/entry pointing-constraint omission and replace the 45 degree cone that currently stands in for how far the vehicle can pitch. The SCvx machinery drag already required covers most of what rotational dynamics would need.
• Better drag: a transonic peak rather than a monotone rise, grid-fin area during entry, angle-of-attack dependence. The current model is biased ~4.5 % long on the coast, which propagates into the ~2 % early profile.
• Sensor noise and state estimation. The MPC currently receives the exact truth state every tick.
• Robust or stochastic MPC. Monte Carlo measures the dispersion; the controller is still nominal. Tube or chance constraints and a designed safety controller would make solver failure an explicit part of the control design.
• Swapping the physics for Basilisk if higher-fidelity gravity, actuator and sensor models matter more than having every equation visible.

Background

Lossless convexification for powered-descent guidance: Acikmese & Ploen 2007; Acikmese, Carson & Blackmore 2013 (G-FOLD). Lars Blackmore, a co-author on that work, is SpaceX's Principal Landing GNC engineer and has stated publicly that this convex-optimization approach is conceptually the basis of Falcon 9's landing guidance. This project uses the published method and public data only.

License

MIT, see LICENSE. The vendored telemetry is public domain (Unlicense); attribution, dependency licenses and trademark notes are in NOTICE.md.

Not affiliated with, endorsed by, or sponsored by SpaceX. "SpaceX" and "Falcon 9" are trademarks of Space Exploration Technologies Corp., used here descriptively to name the real vehicle and the real missions the public data set refers to.