An introduction to using deep reinforcement learning to autonomously optimize fuel consumption across multi-planetary missions. You can visit mazzola.dev​​​​‌‍​‍​‍‌‍‌​‍‌‍‍‌‌‍‌‌‍‍‌‌‍‍​‍​‍​‍‍​‍​‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍​‍​‍​​‍​‍‌‍‍​‌​‍‌‍‌‌‌‍‌‍​‍​‍​‍‍​‍​‍‌‍‍​‌‌​‌‌​‌​​‌​​‍‍​‍​‍‌‍​‍‌‍​‌‌​‌‍‍‌‌​​‍‍‌‍​‌‌‍‌​‌‍‌‌‍‍‌‌‍‍​‍‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​‌‍‌​‍​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌‌‌‌​‌‍‌‌​​‌‍‌‌‌​​‍‌‌‍​​‍‌‌​‌‌‌‌‌‌‌‌​‌‌‌​​‌​‍‌‌​‌‌​‌​‍‌​‌​‌‌‍‍‌​‌​​‌‌‌‌‌​‌‌‍‍‍‌‍‌‌‌‌‌‌‌​‌‌​‌​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌​​‍‌‍​‍‌‌‌​‌‍‍‌​​‌‌​‌‌‍​‌‌‌​‌​‍‌‌‌​‌‌​​‌‍‌‍‌​‌​​‍‌‍‍‌‌​​​‍‌‌‌​‌‌​‌​‍‌‌​‍‌​‍‍‌‌​‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​​‍​‌‍‌‍‌‍‍‌‌‍‌‌‌‍​‌‍‌​‌‌​​‌‍​‌‌‌​‌‍‍​​‌‌‍​‍‌‍​‌‍‌‍‌‌‌‌​​‌‍‌​‌‌​‌​​‍‍​​​​‍‍‌‍​‌‌​‍‌‌​‌‍‍‌‌‍​‌‍​‌‍‌‌​‍​‍‌‌ for nicer formatting.

01 — The Problem with Getting Somewhere

Space is expensive because of the fuel. Every kilogram you want to move requires more fuel, adding more weight, and finally requiring even more fuel to lift. This is described by the Tsiolkovsky rocket equation and it makes the cost of space-travel increase exponentially. Delivering one kilogram to low Earth orbit costs roughly $7,342. Delivering that same kilogram to the surface of the Moon costs nearly $1.2 million. Push past Mars, and the numbers become staggering.

Destination | Cost per kg

Low Earth Orbit | ~$7,342

Moon Surface | ~$1.2 million

NASA Mars Sample Return (total mission) | $8–11 billion

There is, however, a trick that spacecraft engineers have long exploited: the gravity assist. By going past a planet at the right angle and speed, spacecraft can “steal” energy from a planet’s orbital motion around the Sun, dramatically reducing the fuel required for deep-space missions. The Voyager probes, launched in 1977, used a rare planetary alignment to achieve a “Grand Tour” of the outer solar system that would have been physically impossible without gravity assists. Cassini used three separate gravity assists (Venus twice, then Earth) just to build up enough speed to reach Saturn.

Left: Voyager 2's trajectory through the outer solar system, using gravity assists at Jupiter, Saturn, Uranus, and Neptune.
Right: How a gravity assist works in two reference frames: in the planet's frame the spacecraft's speed is unchanged, but in the Sun's frame it gains velocity from the planet's orbital motion.

The problem is that planning these trajectories is extraordinarily complex. The search space for multi-gravity-assist routes is enormous as every possible sequence and timing represents a different option, and traditional optimization methods struggle with this “combinatorial explosion.” These methods are also slow as they require extensive computation up front and cannot adapt when reality differs from the model. As a result, they need ground operators to intervene every time a correction is required. For missions to the outer planets, where communication delays can exceed an hour each way, waiting for a human to correct course is simply not an option.

The Core Problem: Traditional methods find near-optimal trajectories but are rigid and computationally expensive. A spacecraft that could learn its own guidance policy that adapts in real-time to sensor noise and unexpected perturbations would be a fundamentally different kind of mission capability.

This is the problem that motivates my senior project. I want to explore whether modern deep reinforcement learning, specifically an algorithm called Proximal Policy Optimization (PPO), can learn fuel-optimal trajectory policies for multi-planetary missions while remaining robust to real-world uncertainties while still running onboard in milliseconds.

02 — Week 1: Into the Literature

This first week was devoted entirely to reading to build a strong theoretical foundation. My reading list covered three major sources, each addressing a different layer of the problem: the mathematics of learning, the algorithm I plan to use, and the state of the art in applying that algorithm to spacecraft.

Source | Authors | Focus

Reinforcement Learning: An Introduction (2nd ed.) — Ch. 1, 3, 5, 13 | Sutton & Barto, 2018 | MDPs, Monte Carlo, Policy Gradients

Proximal Policy Optimization Algorithms | Schulman, Wolski, Dhariwal, Radford & Klimov, 2017 | Core algorithm for the project

Reinforcement Learning for Low-Thrust Trajectory Design of Interplanetary Missions | Zavoli & Federici, 2021 | Most directly related prior work

What follows is a walkthrough of the key ideas I pulled from each source and how they connect to what I’m trying to build.

03 — Learning by Doing: The RL Framework

Sutton and Barto’s textbook is essentially the foundational reference for the entire field of reinforcement learning. Rather than learning from labeled examples (supervised learning) or by finding structure in data (unsupervised learning), RL is about learning from interaction. An agent takes actions in an environment, observes what happens, and gradually learns which actions lead to the most reward. As Sutton and Barto describe it in Chapter 1, the RL agent is always trying to maximize a cumulative reward signal, and crucially, nobody tells it how to do this. It must discover this on its own through trial and error.

Chapter 3 formalizes this with the concept of the Markov Decision Process (MDP), which is the mathematical framework used to model trajectory problems. In an MDP, at each time step t, an agent observes a state s, chooses an action a, receives a reward r, and transitions to a new state s’. The agent’s goal is to find a policy (a mapping from states to actions) maximizing its expected cumulative future reward.

Spacecraft trajectory optimization surprisingly fits naturally into this framework. The spacecraft’s position, velocity, and remaining fuel mass define the state. The thrust vector (magnitude and direction) defines the action. The negative of fuel consumed is a simple reward. And the optimal policy is the one that gets the spacecraft to its destination using the least propellant. Zavoli and Federici (2021) explicitly formalize a three-dimensional Earth-Mars rendezvous mission this way, defining the spacecraft state at time step k as a seven-dimensional vector comprising position, velocity, and mass.

Chapter 5 introduced me to Monte Carlo methods, which I’ll need for validating my trained policies. The key idea in MC methods is that you can estimate the expected return of a policy by simply running it many times and averaging the results. I’ll run the same policy across hundreds or thousands of randomly varied scenarios, adding different levels of sensor noise and unexpected perturbations, and measure how often it succeeds. Sutton and Barto frame this as estimating value functions directly from sampled episodes, which makes MC methods especially powerful in environments where we don’t have a perfect model of the dynamics (which is always true for real spacecraft).

Chapter 13 on Policy Gradient Methods was perhaps the most directly relevant to my project. Rather than learning the value of states and then deriving a policy, policy gradient methods directly optimize the parameters of the policy itself, following the gradient of expected return. For continuous action spaces, the same as spacecraft thrust commands, Sutton and Barto recommend parameterizing the policy as a Gaussian distribution, with the neural network outputting the mean (and sometimes standard deviation) of that distribution. This allows the agent to output continuous-valued thrust vectors rather than being restricted to discrete choices.

Key Insight from Sutton & Barto: The critical tradeoff in RL is between exploration (trying new actions to learn what is better) and exploitation (using what has already been learned). For spacecraft guidance, this has a real cost: exploration means occasionally wasting fuel on suboptimal maneuvers during training. The design of the reward function and the training environment must balance these pressures carefully.

04 — PPO: Stable Learning at Scale

The Schulman et al. (2017) paper introduced Proximal Policy Optimization, which has since become one of the most widely used deep RL algorithms. When you update a neural network policy using gradient ascent, it’s easy to accidentally take a step that drastically changes the policy’s behavior, destroying everything it has learned. Earlier approaches like Trust Region Policy Optimization (TRPO) addressed this with an explicit constraint on how much the policy could change per update, but TRPO required second-order optimization (expensive and complex to implement).

PPO’s insight is to instead use a clipped surrogate objective, which modifies the loss function to discourage large policy updates without requiring any expensive constraint machinery. Let r(θ) be the probability ratio between the new and old policy for a given action. The clipped objective is:

# PPO Clipped Objective (Schulman et al., 2017, Eq. 7)
r_t(θ) = π_θ(a_t | s_t) / π_θ_old(a_t | s_t)
L^CLIP(θ) = E[ min( r_t(θ)·Â_t, clip(r_t(θ), 1-ε, 1+ε)·Â_t ) ]
# ε = 0.2 is the standard hyperparameter
# Â_t = advantage estimate at timestep t
# Full Objective (Eq. 9)
L^total(θ) = E[ L^CLIP(θ) - c₁·L^VF(θ) + c₂·S[π_θ](s_t) ]
# L^VF = value function MSE loss (Critic)
# S[π] = entropy bonus (encourages exploration)

The clipping means that if an update would move the policy too far in a beneficial direction (i.e., when rt(θ) > 1 + ε and Ât > 0), that excess gain is simply ignored and the algorithm doesn’t get credit for going too far. This creates a natural “pessimistic bound” that discourages aggressive updates. Schulman et al. found that using ε = 0.2 performed best in their benchmarks, outperforming TRPO on almost all continuous control tasks while being dramatically simpler to implement.

The entropy bonus is particularly important for my application. Without it, the policy could converge prematurely to a deterministic and suboptimal solution. Adding entropy to the objective pushes the policy to remain somewhat stochastic during training, maintaining exploration. Schulman et al. validated PPO across 49 Atari games and 7 MuJoCo continuous control benchmarks, showing it struck the best overall balance between complexity and final performance.

Why PPO for Spacecraft? Spacecraft trajectory optimization involves a high-dimensional, continuous action space (thrust magnitude and direction) with complex, nonlinear dynamics. PPO’s strengths of stable training without second-order methods and empirically strong performance on continuous control make it a natural fit.

05 — Prior Works

Zavoli and Federici's (2021) paper is the work I keep returning to as the closest precedent
for my project. They apply PPO to a three-dimensional, time-fixed Earth-Mars rendezvous mission. It's a minimum-fuel problem where a spacecraft must match Mars's position and velocity at a specified arrival time. What makes the paper especially relevant is that they explicitly introduce four distinct types of uncertainty: state noise (unmodeled dynamics), observation noise (sensor error), control actuation errors (thruster imprecision), and missed thrust events (complete loss of thrust for one or more consecutive time steps).

The MDP formulation they use maps directly onto the Sutton & Barto framework I studied this week. The spacecraft state at each time step is a seven-dimensional vector of position, velocity, and mass. The commanded action is an impulsive delta-V (a change in velocity), subject to a maximum magnitude constraint derived from the Tsiolkovsky equation and the Sims-Flanagan trajectory model. The reward function penalizes propellant consumption, constraint violations on the maximum delta-V, and terminal state errors (missing Mars at the end).

In the deterministic (no-noise) case, the RL policy trained by PPO achieved a final spacecraft mass of 600.23 kg which is nearly identical to the 599.98 kg obtained by a classical indirect optimization method that took seconds to compute (versus 2–3 hours for the RL training). This validates that PPO can, in principle, match optimal solutions. The real power emerged in the stochastic scenarios: when the policy was trained with state uncertainty and then tested across 500 Monte Carlo episodes, it succeeded 80% of the time. By contrast, the deterministic policy (trained without noise) succeeded essentially 0% of the time in those same noisy conditions as it had no ability to adapt.

The Gap This Paper Leaves Open: Zavoli and Federici focus on a single fixed mission: Earth to Mars in 358.79 days with one trajectory arc. They explicitly note that the high computational cost “discourages the use of a model-free RL algorithm” for deterministic problems. Multi-gravity-assist trajectories, which involve sequential planetary encounters and vastly larger search spaces, remain largely unexplored in the RL literature. This is the gap I aim to address.

Zavoli and Federici used an ε-constraint relaxation technique, where the terminal constraint tolerance starts loose (10−2) and tightens to the final value (10−3) halfway through training. This gives the policy room to explore freely at first before being held to strict accuracy standards. I anticipate needing a similar curriculum-style approach in my own work, especially as mission complexity increases with multi-gravity-assist scenarios.

06 — The Question I’m Trying to Answer

All three readings this week converge on a single, sharp inquiry: Can deep reinforcement learning, specifically PPO, enable a spacecraft to autonomously optimize fuel consumption across multi-gravity-assist trajectories while correcting for real-world uncertainties in real time?

This breaks down into three sub-questions that will drive my experimental design:

Sub-Question 1: What RL architecture (feedforward MLP, LSTM, or transformer) most effectively learns optimal control policies across the different phases of a multi-gravity-assist mission?

Sub-Question 2: How robust are the learned policies to uncertainties like state noise, observation noise, control errors, and missed thrust events as quantified by Monte Carlo validation campaigns?

Sub-Question 3: Can a trained RL policy achieve fuel efficiency comparable to or better than traditional methods (indirect optimization, heuristic controllers) while delivering real-time inference in milliseconds?

The third question is the most practically significant. Sutton and Barto describe RL policies as fundamentally different from traditional planners: once trained, they are just a single forward pass through a neural network, which would take milliseconds regardless of mission complexity. An indirect optimizer for a complex gravity-assist sequence might take hours or days. For a spacecraft already en route to Jupiter with communication delays of 35+ minutes, the difference between a 2-millisecond onboard decision and a 70-minute round-trip consultation with Earth could be the difference between mission success and failure.

07 — What I Believe I’ll Find

I believe that a PPO-trained policy, when given a well-designed reward function and integrated with a high-fidelity astrodynamics simulator (Basilisk), will successfully learn fuel-optimal guidance for at least simplified gravity-assist trajectories, and that this policy will generalize robustly to uncertain environments in a way that deterministic trajectories simply cannot.

Specifically, I expect the RL policy to trade a small amount of propellant efficiency (perhaps 1–5% worse than the theoretical optimum) in exchange for dramatically improved robustness, or the ability to correct mid-flight for sensor noise, thruster errors, and missed burns. This is precisely the tradeoff Zavoli and Federici observed: their robust policies used slightly more fuel in ideal conditions but maintained mission success rates above 70% in the presence of severe perturbations, while the deterministic policy failed entirely.

I also believe the architecture choice will matter. Gravity-assist trajectories have distinct mission phases (cruise, approach, flyby, post-flyby correction), and I suspect that architectures with memory like LSTMs may outperform simple feedforward networks by better representing the temporal structure of the trajectory. This remains an open question in the literature that I hope to shed some light on.

The deeper question I’m exploring is whether RL can connect theoretically optimal trajectories and practical real-time implementation. Traditional methods often produce trajectories that assume perfect knowledge and a perfectly performing spacecraft. The real universe is messier. Reinforcement learning, as Sutton and Barto describe it, is fundamentally a framework for decision-making under uncertainty. An RL agent learns a policy, which maps any possible situation to the best action in that situation. That is exactly what a spacecraft navigating an unpredictable deep-space environment needs.

08 — Looking Ahead

Next week, I’ll be shifting from theory to orbital mechanics. I’ll be reading Schaub and Junkins on the analytical mechanics of space systems, and Sarli and Taheri on interplanetary gravity-assist trajectory optimization. I’ll also be setting up my development environment, which means installing the Basilisk astrodynamics simulation framework from the University of Colorado Boulder and verifying that its Python API works as expected.

One thing that’s been on my mind after this first week of reading: the Zavoli and Federici paper used a relatively simple 40-step mission model (Sims-Flanagan with 40 impulses) and still took 10–12 hours of training on a modern CPU. My project aims to incorporate more complex dynamics through Basilisk’s higher-fidelity simulation environment. Managing computational cost is going to be a central challenge. I may need to explore techniques like parallelized environment simulation and curriculum learning to make training tractable.

I’m curious whether anyone reading this has experience with integrating high-fidelity simulators into RL training loops. The latency per environment step is the critical bottleneck, and I’d love to hear how others have approached this kind of simulation-in-the-loop RL problem.

Diving into orbital mechanics, Lambert's problem, and what autonomous gravity-assist control research tells us about where reinforcement learning needs to go next. Visit mazzola.dev​​​​‌‍​‍​‍‌‍‌​‍‌‍‍‌‌‍‌‌‍‍‌‌‍‍​‍​‍​‍‍​‍​‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍​‍​‍​​‍​‍‌‍‍​‌​‍‌‍‌‌‌‍‌‍​‍​‍​‍‍​‍​‍‌‍‍​‌‌​‌‌​‌​​‌​​‍‍​‍​‍‌‍​‍‌‍​‌‌​‌‍‍‌‌​​‍‍‌‍​‌‌‍‌​‌‍‌‌‍‍‌‌‍‍​‍‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​‌‍‌​‍​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌‌‌‌​‌‍‌‌​​‌‍‌‌‌​​‍‌‌‍​​‍‌‌​‌‌‌‌‌‌‌‌​‌‌‌​​‌​‍‌‌​‌‌​‌​‍‌​‌​‌‌‍‍‌​‌​​‌‌‌‌‌​‌‌‍‍‍‌‍‌‌‌‌‌‌‌​‌‌​‌​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌​​‍‌‍​‍‌‌‌​‌‍‍‌​​‌‌​‌‌‍​‌‌‌​‌​‍‌‌‌​‌‌​​‌‍‌‍‌​‌​​‍‌‍‍‌‌​​​‍‌‌‌​‌‌​‌​‍‌‌​‍‌​‍‍‌‌​‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​​‍​‌‍‌‍‌‍‍‌‌‍‌‌‌‍​‌‍‌​‌‌​​‌‍​‌‌‌​‌‍‍​​‌‌‍​‍‌‍​‌‍‌‍‌‌‌‌​​‌‍‌​‌‌​‌​​‍‍​​‌​‍‍‌‍​‌‌​‍‌‌​‌‍‍‌‌‍​‌‍​‌‍‌‌​‍​‍‌‌ for nicer formatting.

01 — From Algorithms to Astrodynamics

Welcome back. Last week I researched the learning side of the problem, like the Markov decision processes and policy gradients. This week I pivoted entirely to the physics and mathematics of the trajectories themselves. To build a good environment for an RL agent to learn in, I need to understand the dynamics the environment is supposed to model, so this week was about building that understanding from first principles.

The reading covered orbital mechanics from Schaub and Junkins’ Analytical Mechanics of Space Systems, then shifted to specialized literature on gravity-assist trajectory design, Lambert’s problem, and an investigation of autonomous flyby control by Negri and Prado (2023), which managed to highlight a major challenge for my project. What follows is a walkthrough of the key ideas and how they fit together.

Source | Authors | Focus

Analytical Mechanics of Space Systems (3rd ed.) — Ch. 1, 2, 9 | Schaub & Junkins, 2014 | Two-body dynamics, orbital elements, conics, Lambert’s problem

Revisiting Lambert’s Problem | Izzo, D., 2015 — Celestial Mechanics and Dynamical Astronomy | Efficient algorithm for solving orbital boundary-value problems

Study on Autonomous Gravity-Assists with a Path-Following Control | Negri & Prado, 2023 — arXiv:2307.06185 | Sliding-mode flyby control, Monte Carlo validation, Jovian tours

Interplanetary Gravity-Assist Fuel-Optimal Trajectory Optimization | Sarli & Taheri, 2020 — Celestial Mechanics and Dynamical Astronomy | Optimal control for gravity-assist sequences with perturbations

02 — The Two-Body Problem: Where Everything Starts

Schaub and Junkins ground everything in the two-body problem, which consists of a spacecraft of negligible mass orbiting a single gravitational body. In this idealized case, the equations of motion reduce to a very compact form. The spacecraft’s acceleration is entirely determined by the gravitational parameter μ = GM of the central body and the spacecraft’s distance r from it. The vis-viva equation gives the speed at any point on a conic orbit:

# Two-body equation of motion (Schaub & Junkins, Ch. 2)
r̈ = -(μ/r³)·r
# Vis-viva equation: speed at any point on a conic orbit
v² = μ(2/r - 1/a)
# a = semi-major axis
# a > 0 → ellipse
# a < 0 → hyperbola
# a → ∞ → parabola

The vis-viva equation is the single most useful result in orbital mechanics. It tells you the spacecraft’s speed at any point on any conic orbit, given only the orbit’s semi-major axis and the current distance from the gravitational body. For gravity-assist design, it is essential because the approach speed at a planet, the speed through periapsis of the flyby hyperbola, and the departure speed are all just applications of vis-viva.

The two-body problem’s solutions are conic sections: circles, ellipses, parabolas, and hyperbolas, which can conveniently be described by six numbers called the classical orbital elements: the semi-major axis a, eccentricity e, inclination i, longitude of the ascending node Ω, argument of periapsis ω, and true anomaly ν. These six numbers completely describe where a spacecraft is and where it is going in a two-body gravitational field. Schaub and Junkins spend considerable effort showing how to convert between these elements and the Cartesian position-velocity state vector, which is what a simulation like Basilisk actually propagates.

For my RL environment, the spacecraft’s state needs to encode all information relevant to making good decisions. The Cartesian state vector (position r, velocity v, mass m) is the natural choice because Basilisk propagates in Cartesian coordinates. However, the orbital elements are far more interpretable since an RL agent that can effectively “see” its current orbit shape and orientation relative to its target may learn faster and generalize better than one that only sees raw coordinates. How to represent state is one of the design choices I’ll need to experiment with.

Design Question #1: Should the RL agent observe Cartesian state vectors (r, v, m) or Keplerian orbital elements (a, e, i, Ω, ω, ν, m)? Or both? The orbital elements representation may help the agent learn faster by providing a physically meaningful structure, but it also introduces singularities at zero eccentricity and zero inclination that need to be handled carefully.

03 — Hyperbolic Orbits: The Grammar of the Gravity Assist

Much of orbital mechanics focuses on elliptical orbits, which are the closed paths that planets and satellites follow. But gravity assists happen on hyperbolic orbits, the curves traced by objects that have more than enough energy to escape a gravitational body. Understanding hyperbolic orbits in detail was one of the most important things I did this week because they are the fundamental unit of every gravity-assist maneuver.

On a hyperbolic orbit, the spacecraft arrives from effectively infinite distance with a nonzero speed relative to the planet, called the hyperbolic excess velocity v∞. It swings around the planet along a hyperbola, getting as close as the periapsis radius rp, and then departs in a completely different direction, also with speed v∞ in the planet’s reference frame. Schaub and Junkins derive from the vis-viva equation that the magnitude of v∞ is conserved through the flyby. Hence, the gravity-assist changes the direction of v∞, not its magnitude, in the planet’s frame.

# Hyperbolic orbit (Schaub & Junkins, Ch. 9)
v²_∞ = μ·(-1/a) (a < 0 for hyperbola)
# Half-turning angle: how much the flyby deflects the trajectory
sin δ = μ / (μ + r_p · v²_∞)
# Total turning angle is 2δ; smaller r_p → larger deflection
Δv_∞ = 2·v_∞·sin δ (velocity change in Sun's frame)

The turning angle 2δ determines how much the flyby redirects the spacecraft’s heliocentric velocity. A closer periapsis altitude means a tighter hyperbola, a larger turning angle, and a more powerful assist. But obviously the periapsis must stay above the planet’s surface (and ideally well above the atmosphere). As a result you have to strike a balance between large deflections and the physical constraints of the flyby geometry.

What surprised me in the Sarli and Taheri (2020) paper was how much perturbations complicate this clean picture. They model solar radiation pressure on the spacecraft alongside the planetary gravity assists themselves, and show that the optimal trajectory must account for these forces continuously throughout the transfer and not just at the flyby encounters. The perturbations are relatively small, but over months-long transfers to the outer planets they accumulate into significant trajectory deviations. A trained RL policy that can correct continuously for small errors is therefore far more valuable than a pre-computed optimal trajectory that assumes a clean dynamical model.

Key Insight from Sarli & Taheri: Gravity-assist trajectories are typically designed using the patched-conics approximation, which treats each gravitational body in isolation and stitches together two-body solutions at the boundaries of each planet’s sphere of influence. While this greatly simplifies the math, it introduces errors that must be corrected with additional delta-v. Higher-fidelity N-body simulations reveal that these patched-conics solutions can diverge substantially from real trajectories, especially in the complex gravitational environments of the outer planetary systems.

04 — Lambert’s Problem: Connecting Two Points in Space

The problem is that given two positions in space and a specified travel time, what is the orbit connecting them? This is a boundary-value problem that appears everywhere in trajectory design. Any time you want to know “what velocity must I leave point A with to reach point B in exactly Δt seconds?”, you are solving Lambert’s problem.

The mathematical elegance of Izzo’s (2015) treatment was something I genuinely was not expecting from what I assumed would be a dry algorithmic paper. Izzo builds on Lancaster and Blanchard’s classical formulation and introduces a new variable transformation. All Lambert problems with the same ratio c/s (chord length over semi-perimeter of the triangle formed by the two position vectors) are geometrically equivalent under what Gooding called “L-similarity” as they share the same non-dimensional time-of-flight curves. This collapses what looks like a family of infinitely many problems into a single universal parameter space.

# Lambert's problem: key parameters (Izzo, 2015)
λ² = 1 - c/s, s = (r₁ + r₂ + c) / 2
# Non-dimensional time of flight
T = (1/2)·√(μ/a_m³)·(t₂ - t₁), where a_m = s/2
# T is a function of x (the iteration variable) and λ alone
T(x, λ) = [1/(1-x²)] · [(ψ + Mπ)/√|1-x²| - x + λy]

Izzo’s solver converges very precisely (errors below 10−13) in an average of 2.1 iterations for single-revolution transfers and about 3.3 iterations for multi-revolution cases, using a Householder iterative scheme fed by analytically derived initial guesses. For context, this is faster than the previously state-of-the-art Gooding algorithm, which required about 3 iterations per solve. Since a trajectory design problem may involve thousands or millions of Lambert solves (for example, sweeping over all possible departure and arrival dates to build a “porkchop plot”), the difference matters enormously.

For my project, the Lambert solver is the backbone of the trajectory planning baseline I need to compare against. When I want to check whether my RL policy is doing something sensible, or when I need to construct a reference trajectory for the patched-conics tour design, Lambert solves are the tool. Izzo’s solver is distributed as part of the open-source PyKEP library from ESA, which conveniently interfaces with Python. This is the same library that Zavoli and Federici used in their RL framework, so the integration path is well-tested.

Why Lambert Matters for RL: Lambert’s problem gives me a way to construct “expert” reference trajectories for any desired departure-arrival pair. These can serve as a sanity check for RL policy behavior, as targets for reward shaping, and as comparison baselines when measuring how much efficiency the RL policy sacrifices in exchange for robustness. Izzo’s efficient solver means I can generate these references quickly enough to use them during training.

05 — Autonomous Flyby Control & Limits of Patched Conics

The most immediately useful paper this week was Negri and Prado’s (2023) study on autonomous control of gravity-assist hyperbolic trajectories. Negri and Prado zoom in on what actually happens during the flyby itself. Their central finding is important enough that it reframed how I’m thinking about the scope of my own project.

Their control strategy, called the Robust Keplerian Path-Following Control (RKPFC), is built on sliding mode control theory. Rather than tracking a specific reference trajectory in time (which would require the spacecraft to be at exactly the right place at exactly the right moment), the RKPFC drives the spacecraft to the correct orbital geometry, or the right shape and orientation of the hyperbola, without caring precisely when it arrives at each point. The sliding surface is defined in terms of integrals of motion:

# RKPFC sliding surface (Negri & Prado, 2023, Eq. 20)
s = [ ẽ · (λ_R·r̂ + θ̂) ]
[ h̃ · (ĥ_d · (λ_N·r̂ + θ̂)) ] = 0
# ẽ = e - e_d : eccentricity vector error
# h̃ = h - h_d : angular momentum magnitude error
# When s = 0: spacecraft converges to desired Keplerian geometry

“Path following” rather than “trajectory tracking” is as a result the right control objective for a gravity assist. Once you’re inside a planet’s sphere of influence, what matters is preserving the hyperbolic geometry (the eccentricity, inclination, periapsis) that will produce the desired velocity change. A few minutes early or late at periapsis barely matters, while arriving with the wrong orbit shape destroys the assist entirely.

What does the Monte Carlo validation show? Without control authority, Negri and Prado show that initial position errors of around 50 km (a realistic figure for autonomous guidance) produce wildly divergent outgoing trajectories. In the Titan simulations, the outgoing semi-major axis has a standard deviation of over 2,400 km, and some trajectories interact catastrophically with Titan’s atmosphere. With the RKPFC active, the standard deviation goes to essentially zero and all trajectories converge to the desired geometry.

Negri and Prado also simulated a Callisto-Ganymede-Europa-Io tour designed using a patched-conics model, then executed it in a more realistic N-body simulation of the Jovian system. However, even with their robust flyby controller active, the tour fails. The encounter with Europa occurs early due to chaotic N-body interactions, and the trajectory diverges completely from the design thereafter. A flyby control strategy alone cannot guarantee a tour based on patched-conics guidance in a truly chaotic multi-body environment. The guidance layer deciding the overall trajectory between flybys needs to be high-fidelity and adaptive.

The Critical Gap This Paper Exposed: Negri and Prado show that robust flyby control is a solved problem for individual encounters. The open challenge is the guidance layer that connects encounters into a coherent tour under N-body dynamics. A trained RL policy could potentially learn to account for the chaotic character of the Jovian system directly from simulation experience, treating the chaos as part of the environment rather than a perturbation to be corrected.

06 — The Patched-Conics Approximation and Its Limitations

Most gravity-assist trajectory design relies on the patched-conics approximation. In the zero-SOI version, the encounter with a planet is treated as an instantaneous impulsive maneuver at the planet’s center. In the more refined “finite-SOI” patched-conics model that Negri and Prado describe, the spacecraft’s trajectory inside the sphere of influence is modeled as a genuine two-body hyperbola, and a Lambert solver connects the arcs between successive encounters. The decision variables for the tour become the encounter times tj and periapsis radii rp,j at each flyby body, and the optimization minimizes the sum of all corrective delta-v’s needed to connect the arcs:

# Patched-conics tour optimization objective (Negri & Prado, 2023, Eq. 17)
minimize: Σ_j ( ||ΔV_j⁻|| + ||ΔV_j⁺|| )
# ΔV_j⁻ : impulse before j-th flyby to achieve the desired approach
# ΔV_j⁺ : impulse after j-th flyby to continue on the designed transfer arc
# Decision variables: t_j, r_{p,j}, e_j, i_j, Ω_j, ω_j for each encounter

This framework is clean and computationally manageable, which is why it has been the workhorse of interplanetary mission design for decades. But Negri and Prado’s N-body simulations show that in the outer planets, patched conics can be dramatically wrong. The Jovian moons move fast (Io completes a full orbit in under two days), so even the tens of minutes a spacecraft spends inside a moon’s sphere of influence are enough for the moon to have moved significantly from where the patched-conics model assumed it was. The accumulated error across a multi-flyby tour makes the patched-conics solution essentially a warm-start guess rather than an executable flight plan.

This is where my thinking about the project is sharpening. For the inner solar system (Earth-Venus-Mars gravity assists), patched conics is a reasonable model and RL might simply need to learn to correct small deviations. For outer planet tours, which are the missions where gravity assists matter most, the RL agent will need to cope with fundamentally chaotic N-body dynamics. This could mean training in a high-fidelity simulation from the start, even if it is computationally expensive, rather than first training on a simple patched-conics model and expecting to transfer.

Design Question #2: Should I train the RL policy in a low-fidelity patched-conics environment first (faster training, cleaner reward signal) and then fine-tune in the high-fidelity Basilisk N-body environment, or is the dynamical gap so large that this transfer is harmful rather than helpful? Curriculum learning suggests the former might work, but the Negri & Prado results suggest the gap may be qualitative, not just quantitative.

07 — Building the State Space from First Principles

The spacecraft’s state at any instant lives in a seven-dimensional space: three position components, three velocity components, and one for remaining mass. Basilisk propagates this state forward using high-fidelity numerical integration, accounting for N-body gravitational forces, solar radiation pressure, and atmospheric drag where relevant. The RL agent observes this state and outputs a three-dimensional thrust vector.

But raw Cartesian state is not necessarily the most informative representation for the agent. Inspired by Negri and Prado’s use of the eccentricity vector e and angular momentum vector h as the quantities to be controlled during a flyby, I’m considering augmenting the Cartesian state with derived orbital elements. An agent that can “see” how far its current orbit deviates from the desired flyby hyperbola can form much more targeted commands than one reasoning purely from position and velocity numbers.

However, unlike elliptical orbits, where standard Keplerian elements are always well-defined, the orbital elements of a hyperbolic trajectory inside a planet’s SOI are defined relative to that planet’s gravitational field, not the Sun’s. The RL agent will need to switch reference frames as it moves between gravitational regimes, from heliocentric two-body during the interplanetary cruise, to planetocentric two-body during the flyby, and back again. The patched-conics approximation handles this transition (somewhat crudely) by stitching orbits together at the sphere of influence boundary. Building a smooth and numerically stable version of this for an RL observation function is a concrete implementation challenge ahead of me.

08 — Looking Ahead

Next week I’ll be finishing the reading phase and moving into environment implementation. The first concrete milestone is getting Basilisk installed and running simple test propagations, confirming that a spacecraft in a Keplerian orbit stays on that orbit, that a hyperbolic flyby gives the expected velocity change, and that the sphere-of-influence transitions behave consistently. These are sanity checks, but failing them early is far better than failing them in the middle of RL training.

Two things from this week’s reading are sitting unresolved in my mind. The first is the question of simulation fidelity versus training speed. Negri and Prado’s results suggest that patched-conics guidance fails in the outer solar system because the dynamics are qualitatively different from the model. That means high-fidelity training may not be optional for multi-gravity-assist missions as it is the only setting where the problem is being posed correctly. This will likely force me to be aggressive about computational optimization from the start, possibly using parallelized environment rollouts or cloud computing.

The second unresolved question is about the guidance layer. Negri and Prado’s flyby controller handles the local problem about keeping the spacecraft on the right hyperbola during an encounter but concludes that the global tour planning problem requires something beyond pure control. Reinforcement learning, trained end-to-end on the full multi-flyby problem, could in principle learn to handle both layers simultaneously. That is an ambitious goal, but it is also precisely what makes this research interesting. I’ll be designing the Markov decision process formulation with this dual-layer structure in mind.

One question for anyone reading this: do you know of any open-source implementations of multi-body tour optimization that use something other than patched-conics as the baseline? I’m aware of GMAT and MONTE from JPL, but both are either proprietary or complex to integrate. Any pointers to Python-based tools for generating high-fidelity gravity-assist reference trajectories would be enormously useful.

Understanding Basilisk from the source and running the first real thruster simulations to understand the building blocks that will become my RL training environment. Visit mazzola.dev​​​​‌‍​‍​‍‌‍‌​‍‌‍‍‌‌‍‌‌‍‍‌‌‍‍​‍​‍​‍‍​‍​‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍​‍​‍​​‍​‍‌‍‍​‌​‍‌‍‌‌‌‍‌‍​‍​‍​‍‍​‍​‍‌‍‍​‌‌​‌‌​‌​​‌​​‍‍​‍​‍‌‍​‍‌‍​‌‌​‌‍‍‌‌​​‍‍‌‍​‌‌‍‌​‌‍‌‌‍‍‌‌‍‍​‍‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​‌‍‌​‍​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌‌‌‌​‌‍‌‌​​‌‍‌‌‌​​‍‌‌‍​​‍‌‌​‌‌‌‌‌‌‌‌​‌‌‌​​‌​‍‌‌​‌‌​‌​‍‌​‌​‌‌‍‍‌​‌​​‌‌‌‌‌​‌‌‍‍‍‌‍‌‌‌‌‌‌‌​‌‌​‌​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌​​‍‌‍​‍‌‌‌​‌‍‍‌​​‌‌​‌‌‍​‌‌‌​‌​‍‌‌‌​‌‌​​‌‍‌‍‌​‌​​‍‌‍‍‌‌​​​‍‌‌‌​‌‌​‌​‍‌‌​‍‌​‍‍‌‌​‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​​‍​‌‍‌‍‌‍‍‌‌‍‌‌‌‍​‌‍‌​‌‌​​‌‍​‌‌‌​‌‍‍​​‌‌‍​‍‌‍​‌‍‌‍‌‌‌‌​​‌‍‌​‌‌​‌​​‍‍​​‍​‍‍‌‍​‌‌​‍‌‌​‌‍‍‌‌‍​‌‍​‌‍‌‌​‍​‍‌‌ for nicer formatting

01 — From Pages to Python

The first two weeks of this project were almost entirely about absorbing the theory of reinforcement learning and mastering the orbital mechanics. This week, the balance finally shifted. With the core literature reviewed and consolidated, I installed the Basilisk astrodynamics simulation framework and ran the first real thruster-based simulations.

This post is longer and more technical than the previous two, because what I worked through this week is dense. To understand how I'll take advantage of Basilisk's design I needed to read the example scripts alongside the architecture principles to understand how the pieces connect. I am going to walk through exactly what I found, because the choices I make in the next week when I start designing my Gymnasium wrapper will be direct consequences of this week's understanding.

Source | Authors | Focus

Basilisk: Fundamentals of Simulations | AVS Lab — bskPrinciples series (What is BSK, Process/Task/Module, Messaging) | Architecture foundation

scenarioOrbitManeuverTH | Schwend, 2024 — Thrusted Hohmann Transfer w/ thrusterStateEffector | Primary thruster example

scenarioPatchedConics | AVS Lab — Earth → Jupiter multi-body transfer with SOI switching | Multi-body / interplanetary

scenarioJupiterArrival | AVS Lab — Hyperbolic arrival orbit with capture burn | Hyperbolic dynamics

scenarioMonteCarloAttRW | AVS Lab — MC Controller, dispersions, retention policy | Validation framework

02 — What Basilisk Actually Is

The AVS Lab at the University of Colorado Boulder describes Basilisk as a software framework, not a program. The distinction matters because Basilisk by itself does not simulate anything. Instead, it provides the building blocks from which spacecraft simulations are assembled. The documentation uses a Lego analogy that I think is exactly right: each module is a self-contained brick, and the user's job is to stack and connect them into the spacecraft behavior they need. The modules encapsulate the physics and algorithms; the user wires them together through a message-passing interface and never needs to understand the math inside each module to use it correctly.

The execution model is built around three concepts: Processes, Tasks, and Modules. A SimBaseClass container holds one or more processes, each of which owns one or more tasks, and each task holds a list of modules that are called in order at a specified time step. Every simulation starts the same way:

# The invariant setup pattern in every Basilisk script
scSim = SimulationBaseClass.SimBaseClass()
dynProcess = scSim.CreateNewProcess("simProcess")
dynProcess.addTask(scSim.CreateNewTask("simTask", macros.sec2nano(2.0)))
# All modules are then added to this task
scSim.AddModelToTask("simTask", scObject) # spacecraft dynamics
scSim.AddModelToTask("simTask", thrusterSet) # thruster effector
scSim.AddModelToTask("simTask", sNavObject) # navigation sensor
scSim.AddModelToTask("simTask", attGuidance) # guidance law
scSim.AddModelToTask("simTask", mrpControl) # feedback controller

The time step, which is 2 seconds in the example above, is specified in nanoseconds using the macros.sec2nano() utility. Everything in Basilisk's internal clock runs in nanoseconds, which initially felt unusual, but it turns out it eliminates floating-point precision issues that would accumulate over multi-day mission simulations.

The modules communicate exclusively through typed messages. A module declares one or more output messages, which other modules subscribe to via their input message connections. There is no shared global state and information only flows through explicitly declared message subscriptions. This matters enormously for my project, because when I wrap Basilisk inside a Gymnasium environment, the "observation" that the RL agent receives will be exactly the output of the simpleNav navigation module, and the "action" it outputs will be written directly into the THRArrayOnTimeCmdMsg thruster command message. The message architecture is, conveniently, already the interface I need.

Architecture Insight: Basilisk's message-passing design has the same form as the MDP formulation I will use. The simpleNav output message is the observation. The THRArrayOnTimeCmdMsg is the action. The simulation stepping from one timestep to the next is the transition function. Building the Gymnasium wrapper will largely consist of wrapping these existing interfaces rather than reimplementing them.

03 — The Thruster Example, Taken Apart

The most important script I worked through this week was scenarioOrbitManeuverTH, which simulates a full Hohmann transfer from Low Earth Orbit to Geostationary Orbit using finite thruster burns rather than the instantaneous impulsive delta-v used in earlier simpler examples. This distinction is critical for my project. Prior Basilisk orbit maneuver examples like scenarioOrbitManeuver just add a velocity vector to the spacecraft state at a specific instant, which is convenient but physically unrealistic. For RL-based control, I need the spacecraft to actually fire a thruster for a computed duration, consuming propellant, building up velocity through a force that acts over time, and with an attitude controller that keeps the spacecraft pointed in the correct burn direction throughout. This example does all of that.

Parameter | Value | Role

Thruster 1 | 3000 N | GTO departure burn (perigee)

Thruster 2 | 420 N | Circularization burn (apogee)

Specific impulse (Isp) | 318 s | Both thrusters

Integrator | RKF45 | Variable-step integrator used for the finite burns

The thrusterStateEffector Module

The thruster is modeled via the thrusterStateEffector module, which is the physically correct approach for finite burns. Unlike a simple force vector added externally, the state effector formulation properly couples the thruster dynamics to the spacecraft's equations of motion. The mass changes as propellant is consumed, and this is tracked self-consistently through the integration. Each thruster is defined via the simIncludeThruster.thrusterFactory() utility, which takes a thruster model type, its location in the body frame, its pointing direction, and parameters including maximum thrust, nozzle area, and specific impulse:

# Create the thruster factory and define two thrusters
thFactory = simIncludeThruster.thrusterFactory()
# Departure burn thruster: 3000 N, high thrust for GTO insertion
thFactory.create('Blank_Thruster', [0,0,1], [0,0,1],
cutoffFrequency=.1,
MaxThrust=3000.0,
areaNozzle=0.046759,
steadyIsp=318.0)
# Circularization thruster: 420 N, lower thrust for GEO insertion
thFactory.create('Blank_Thruster', [0,0,1], [0,0,1],
cutoffFrequency=.1,
MaxThrust=420.0,
areaNozzle=0.046759,
steadyIsp=318.0)
# Attach to spacecraft so the thruster set is now a mass-varying state effector
thFactory.addToSpacecraft("GTOThrusterDynamics", thrusterSet, scObject)

The thruster is commanded through the THRArrayOnTimeCmdMsg message, which specifies how long each thruster in the array should fire. The burn duration is calculated analytically before each phase under a constant-mass assumption, which is standard practice for rough planning when the propellant fraction is small:

# Burn time under constant-mass approximation
# t_burn = (Δv · M_sc) / F_thrust
# From the script: first burn calculation
dv0 = v0p - v0 # required delta-v at perigee
t0Burn = dv0 / (3000.0 / scObject.hub.mHub) # burn duration (seconds)
# Write the on-time message: [thruster_1_time, thruster_2_time]
ThrOnTimeMsgData.OnTimeRequest = [t0Burn, 0]
thrOnTimeMsg.write(ThrOnTimeMsgData, time=simulationTime)

This constant-mass approximation becomes less accurate as the propellant fraction increases. For the deep-space missions I am targeting, where propellant mass can be a large fraction of total spacecraft mass, I will eventually need to use a proper finite-burn propagation that accounts for continuously changing mass through the rocket equation. Basilisk's thrusterStateEffector already does this integration internally and the approximation is only used here to estimate how long to command the burn. The simulation itself correctly propagates the mass change.

Attitude Control During the Burn

A spacecraft cannot just fire a thruster in an arbitrary direction and it has to be pointed correctly first. This example uses a full attitude control loop to keep the spacecraft aligned with the velocity vector throughout the burn, which is the correct direction for a Hohmann transfer. Five modules work together to make this happen, and understanding each of them is essential because I will need an equivalent attitude-keeping loop in my RL environment to make the commanded thrust go in the right direction.

# The full attitude control chain:
# 1. velocityPoint: computes the required attitude (reference frame aligned with velocity)
attGuidance = velocityPoint.velocityPoint()
attGuidance.transNavInMsg.subscribeTo(sNavObject.transOutMsg)
attGuidance.mu = mu
# 2. simpleNav: reads the true spacecraft state and outputs a "sensed" state
# This is the observation function for the RL agent
sNavObject = simpleNav.SimpleNav()
sNavObject.scStateInMsg.subscribeTo(scObject.scStateOutMsg)
# 3. attTrackingError: computes MRP error between commanded and current attitude
attError = attTrackingError.attTrackingError()
attError.attRefInMsg.subscribeTo(attGuidance.attRefOutMsg)
attError.attNavInMsg.subscribeTo(sNavObject.attOutMsg)
# 4. mrpFeedback: PD controller on MRP attitude error gives the torque command
mrpControl = mrpFeedback.mrpFeedback()
mrpControl.guidInMsg.subscribeTo(attError.attGuidOutMsg)
mrpControl.K = 3.5 # proportional gain
mrpControl.P = 30.0 # derivative gain
# 5. extForceTorque: applies the commanded torque to the spacecraft
extFTObject.cmdTorqueInMsg.subscribeTo(mrpControl.cmdTorqueOutMsg)

The velocityPoint module deserves special attention here. It takes the spacecraft's translational navigation state and computes a reference attitude that aligns the spacecraft with the orbit velocity frame. For my RL environment, I will need an equivalent velocity-pointing capability for thrust maneuvers. The attitude control loop runs inside Basilisk; the RL policy only needs to specify when and for how long to fire the thrusters, not how to orient itself during the burn.

Key Design Decision for the RL Wrapper: The thruster command interface in Basilisk is THRArrayOnTimeCmdMsg, which is a vector of on-times for each thruster. The RL agent's action will be exactly this: an array of burn durations, one per thruster, at each decision step. This is a bounded continuous action space, where zero means no burn and the maximum is physically constrained by the available propellant and the orbital maneuver constraints. This maps cleanly onto PPO's Gaussian policy head for continuous control.

The Variable-Timestep Integrator

One subtle but important choice in the example is the use of the svIntegratorRKF45 integrator rather than the default fixed-step RK4. The RKF45 (Runge-Kutta-Fehlberg 4th/5th order) is a variable-step integrator that adapts its step size to maintain a target local truncation error. During the coast phases, where dynamics are smooth, it can take large steps. During and immediately after a burn, when the thrust is changing the trajectory rapidly, it automatically tightens the step size to maintain accuracy.

# Set the RKF45 integrator which is critical for accurate finite burn propagation
integratorObject = svIntegrators.svIntegratorRKF45(scObject)
scObject.setIntegrator(integratorObject)
# The integrator adapts its internal step size to maintain accuracy
# during burn transients, while using large steps during coast arcs.
# This is important when thruster on-time durations are much shorter
# than the orbital period (as they are for finite deep-space burns).

For deep-space missions with long coast arcs interspersed with short burns, this variable-step capability is necessary for getting numerically accurate results without prohibitive computational cost. I had already flagged computational efficiency as a major challenge and the RKF45 integrator is one of Basilisk's tools for managing this tradeoff.

State Recovery: rv2elem and elem2rv

The orbitalMotion utility module provides the classical orbit element conversions I studied from Schaub and Junkins last week. The two most important functions are orbitalMotion.rv2elem(mu, r, v) and its inverse orbitalMotion.elem2rv(mu, oe). In the example, rv2elem is called after the first burn to compute the current orbital elements and then find the time to apogee, where the second burn will occur. This is a direct illustration of the state representation question I raised last week: the simulation propagates in Cartesian coordinates, but physically meaningful decisions (like "when is apogee?") require the orbital element representation.

# After the first burn, recover orbital elements from the Cartesian state
rVt = unitTestSupport.EigenVector3d2np(posRef.getState())
vVt = unitTestSupport.EigenVector3d2np(velRef.getState())
oe1 = orbitalMotion.rv2elem(mu, rVt, vVt)
# Use eccentric anomaly to find the remaining time to apogee
E1 = 2*np.arctan(np.sqrt((1-oe1.e)/(1+oe1.e)) * np.tan(oe1.f/2))
T3Full = np.pi / np.sqrt(mu / oe1.a**3)
T3P = np.sqrt(oe1.a**3 / mu) * (E1 - oe1.e*np.sin(E1))
T3 = macros.sec2nano(T3Full - T3P) # coast time to apogee

The RL agent will need access to this conversion. When constructing the observation vector for the policy network, I plan to augment the raw Cartesian state with the derived orbital elements, as I discussed last week. The orbitalMotion module makes this straightforward because with every step, after calling scSim.ExecuteSimulation(), I can call rv2elem on the current state and append the elements to the observation.

04 — Multi-Body Trajectories and Hyperbolic Arrivals

The thruster example gave me confidence in the core simulation loop, but my project requires more than Earth-orbiting maneuvers. I spent the second half of the week studying two examples that directly address the interplanetary and multi-body dynamics I will need: scenarioPatchedConics and scenarioJupiterArrival.

The scenarioPatchedConics example is particularly instructive because it is the closest existing Basilisk scenario to my eventual problem. It simulates a full Earth-to-Jupiter transfer, beginning with a hyperbolic departure from LEO, proceeding through a heliocentric Hohmann transfer ellipse, and arriving at a Jupiter flyby, all in one continuous simulation. Critically, it does this by running Basilisk's N-body integrator throughout, with Earth, Sun, and Jupiter all contributing to the gravitational acceleration simultaneously. The patched-conics model is used only to design the initial trajectory; the simulation propagates the actual dynamics with all bodies included from the start.

Basilisk always propagates the state in a single inertial frame, but the user changes which body is designated the central body for logging purposes as the spacecraft transitions between gravitational regimes. Early in the simulation, position and velocity are logged relative to Earth. After crossing the Earth's sphere of influence, the central body is switched to the Sun for the heliocentric transfer phase. As the spacecraft approaches Jupiter, the central body switches again. The code manages this by stopping the simulation, adjusting the spacecraft's state vector to account for the offset between reference frames, and restarting with a new logging frame.

# scenarioPatchedConics: frame switching at the Earth SOI boundary
# All three bodies are always contributing gravitational acceleration
# This is not patched conics in the dynamical sense, only in the initial
# design sense. The simulation integrates the full N-body problem.
gravFactory = simIncludeGravBody.gravBodyFactory()
earth = gravFactory.createEarth()
sun = gravFactory.createSun()
earth.isCentralBody = True # initial reference which switches to Sun after SOI exit
# After crossing the Earth SOI, switch the central body for logging
sun.isCentralBody = True
earth.isCentralBody = False
# Adjust spacecraft state to Sun-centered coordinates
# (add Earth's heliocentric velocity to the spacecraft's geocentric velocity)

The standard approach of always using heliocentric inertial coordinates for the state vector, and only switching to planetocentric coordinates for reward computation near flyby bodies, seems to be what Basilisk supports naturally. The message logged by simpleNav always reports state relative to the current central body, so I will need to be careful about which body is designated central at each phase when reading out the observation for the RL agent.

The scenarioJupiterArrival example covers a hyperbolic capture maneuver. The spacecraft arrives on a hyperbolic approach trajectory, propagates to periapsis, and then fires an impulsive delta-v burn to slow down into a circular parking orbit. This example demonstrated concretely how Basilisk handles hyperbolic orbits. The state initialization uses the same orbital element structure (orbitalMotion.ClassicElements()), with a negative semi-major axis indicating the hyperbolic conic. The simulation time for the first phase (initial position to periapsis) is computed from the hyperbolic time equation:

# Hyperbolic time equation:
# N = sqrt(μ / (-a)³) · (t - t_p) = e·tan(ζ) - ln[ tan(ζ/2 + π/4) ]
#
# where ζ = arccos( e·a / (a - r) ) is the hyperbolic anomaly,
# the analog to eccentric anomaly on an ellipse.

The orbitalMotion module handles these conversions, so I do not need to implement the hyperbolic time equation myself, but understanding it is essential for computing meaningful observation features (like "time to periapsis") for the RL agent during a flyby approach.

What the Examples Tell Me About Thruster Use: Neither scenarioPatchedConics nor scenarioJupiterArrival uses thrusters. Instead, they use impulsive delta-v. For interplanetary burns, the delta-v is first computed analytically (using vis-viva or the Lambert solver), and then a thruster burn of the appropriate duration replaces the instantaneous velocity change. The scenarioOrbitManeuverTH example demonstrates this translation pattern on an Earth orbit; I will need to generalize it to heliocentric and planetocentric frames across multiple mission phases.

There is also a dedicated scenarioFlybySpice example that uses real SPICE ephemeris kernels to drive the spacecraft through Venus, Earth, and Mars flybys with accurate planetary positions. This example additionally shows the velocityPoint module being used to align the spacecraft with the velocity frame during flybys and is the exact same attitude guidance module that the thruster example uses during burns. For my project, this SPICE-driven approach becomes relevant once I need to simulate actual historical mission geometries for validation purposes, though I will start with simplified fixed-planet models to keep the first iteration manageable.

05 — The Monte Carlo Machinery

Monte Carlo validation is the primary means by which I will demonstrate that the trained RL policy is actually robust to realistic uncertainties, not just a policy that memorized a specific nominal trajectory. Basilisk has a built-in Monte Carlo framework, and studying it this week showed what kind of validation I should plan for.

The scenarioMonteCarloAttRW example demonstrates the complete MC setup. The entry point is the Controller class from Basilisk.utilities.MonteCarlo. The user specifies the number of runs, the random seeds, and the number of parallel threads. Then dispersions, or the statistical distributions over simulation parameters that will be sampled independently in each run, are defined:

# Basilisk Monte Carlo setup (from scenarioMonteCarloAttRW)
monteCarlo = Controller()
monteCarlo.setShowProgressBar(True)
monteCarlo.setSimulationFunction(myScenario)
monteCarlo.setExecutionCount(500) # 500 independent runs
monteCarlo.setThreadCount(4) # parallel execution
# Define dispersions as each run samples from these independently
# For trajectory work I'll adapt these to spacecraft state/control uncertainties
dispMass = UniformDispersion("scSim.scObject.hub.mHub", [712.5, 787.5])
dispInertia = NormalVectorCartDispersion("scSim.scObject.hub.IHubPntBc_B", ...)
# A retention policy specifies which messages to log across runs
retentionPolicy = RetentionPolicy()
retentionPolicy.addMessageLog(transNavOutMsgName, ["r_BN_N", "v_BN_N"])
monteCarlo.addRetentionPolicy(retentionPolicy)
# Execute and collect. Failed jobs should be zero
failures = monteCarlo.ExecuteSimulation()

For my trajectory optimization application, the dispersions I care about are directly analogous to those used in Zavoli and Federici (2021): initial state noise (position and velocity errors representing navigation uncertainty), control actuation errors (thrust magnitude and direction errors), and missed thrust events (a run where a planned burn partially or fully fails). The Basilisk MC framework can accommodate all of these through its dispersion API.

One important implementation note from the MC documentation: all simulation objects that will be varied between runs must be attached as attributes of the simulation container (scSim). If a module is only stored in a local Python variable, the MC framework cannot find it to apply the dispersion. This is a subtle but critical requirement that I will need to respect when designing my Gymnasium wrapper. The wrapper class will need to store all Basilisk modules as attributes that the MC framework can reference.

Monte Carlo & RL: Two Layers of Validation: There will be two separate uses of Monte Carlo methods in my project. The first is during training: each RL training episode is itself a kind of Monte Carlo trial, where the environment is initialized with random perturbations so the policy learns to be robust. The second is during validation: after training is complete, I run a formal MC campaign with hundreds of fixed scenarios to measure policy performance statistics. Basilisk's built-in MC Controller is designed for the second use; the first will be handled by my Gymnasium wrapper's reset function, which will sample from the dispersion distributions each time a new episode begins.

06 — Mapping Basilisk to a Gymnasium Environment

Having now worked through Basilisk's architecture and the relevant examples in detail, I can state concretely what the Gymnasium wrapper needs to do. This is the central implementation task of next week, and it is worth being explicit about the design before writing any code.

A Gymnasium environment implements four methods: reset(), step(action), render(), and a definition of the observation and action spaces. The mapping to Basilisk is as follows.

# Conceptual mapping: Gymnasium ↔ Basilisk
## reset()
# → create a new SimBaseClass
# → configure spacecraft(), thrusterStateEffector(), simpleNav()
# → initialize orbital elements with randomized dispersion
# → call scSim.InitializeSimulation()
# → return initial observation from sNavObject output messages
## step(action)
# action = array of thruster on-times [t_thr0, t_thr1, ...]
# → write action to THRArrayOnTimeCmdMsg
# → call scSim.ConfigureStopTime(t_current + dt_step)
# → call scSim.ExecuteSimulation()
# → read simpleNav output messages → observation
# → call orbitalMotion.rv2elem() → augment observation with orbital elements
# → compute reward (fuel consumed, trajectory error, constraint violations)
# → check termination condition (mission complete / fuel exhausted / lost)
# → return (observation, reward, done, info)
## observation_space
# → Box(low, high): [r_x, r_y, r_z, v_x, v_y, v_z, mass, a, e, i, Omega, omega, nu]
# positions in m, velocities in m/s, mass in kg, angles in radians
## action_space
# → Box(low=0, high=max_on_time): one on-time per thruster

The biggest open question is the granularity of the RL decision step. Each call to step() advances the simulation by some time interval dt_step. If this is set to 2 seconds (the integration timestep), the RL agent makes thousands of decisions per orbit. If it is set to 1 hour, the agent makes coarser decisions but the episode is tractable in terms of total timesteps. The Zavoli and Federici approach used 40 fixed decision points across the entire Earth-Mars transfer, which amounts to roughly one decision every 9 days. I will begin with a similar coarse discretization and refine based on training stability.

Design Question for Next Week: Should the RL agent command a burn on-time (how long to fire this step) or a delta-v vector (what velocity change to achieve)? The first maps directly to the THRArrayOnTimeCmdMsg interface and is physically grounded. The second is more abstract but may be easier for the policy to learn because it separates the propulsion physics from the trajectory planning. Basilisk's velocityPoint module plus a built-in thrust duration calculator could serve as a lower-level controller, letting the RL agent operate at the delta-v planning level rather than the raw actuator level. I need to decide which abstraction level to start with.

My current leaning is to start at the delta-v abstraction level to simplify early training, then verify that the Basilisk thruster simulation correctly executes those burns. Once I have evidence that the policy is learning something sensible, I can lower the abstraction level toward raw thruster commands as the problem complexity increases.

07 — Looking Ahead

The next two weeks are the most concrete implementation phase of the entire project. Week 4 is for finalizing the MDP formulation of state space, action space, and a preliminary reward function, along with starting the Gymnasium wrapper architecture. Week 5 is for getting the basic wrapper class functional: initialization, state observation, and episode termination logic. The baseline from which I want to depart is a working environment where a random policy can run episodes without crashing, even if it achieves nothing useful. Getting to that baseline, combined with what I learned this week about Basilisk's module structure, is the critical path.

A few specific challenges I am thinking about. The first is the cold start problem of the simulation. Every call to reset() needs to create a fresh simulation object, configure all the modules from scratch, and call InitializeSimulation(). Basilisk does not have a built-in checkpoint-restore mechanism, so there is overhead associated with each episode reset. For RL training with millions of environment steps, this overhead will accumulate. I may need to profile how long InitializeSimulation() takes for different mission configurations and decide whether to accept this cost or find a workaround.

The second challenge is the reference frame consistency issue raised by the patched conics example. During training, episodes will span different mission phases, including Earth departure, heliocentric cruise, flyby approach, and the observation representation needs to be consistent across these phases. I am leaning toward always using heliocentric inertial coordinates for the position and velocity components of the observation, and computing planetocentric elements separately only for reward computation near flyby bodies. This avoids the discontinuities that would arise from switching the primary reference frame mid-episode.

One question I would appreciate feedback on from anyone reading this: for building a Gymnasium-compatible wrapper around a simulation framework like Basilisk, is it better to run the Basilisk simulation in a separate process and communicate with it via inter-process messaging (for robustness and clean separation), or to run it in the same Python process as the RL training loop (for lower latency per step)? My current plan is same-process, but I am wondering whether it may be problematic to reset and re-run many times within the same Python session without memory leaks.

Building the initial codebase and passing the baseline tests. Still a long way from a successfull run, however. Visit mazzola.dev​​​​‌‍​‍​‍‌‍‌​‍‌‍‍‌‌‍‌‌‍‍‌‌‍‍​‍​‍​‍‍​‍​‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍​‍​‍​​‍​‍‌‍‍​‌​‍‌‍‌‌‌‍‌‍​‍​‍​‍‍​‍​‍‌‍‍​‌‌​‌‌​‌​​‌​​‍‍​‍​‍‌‍​‍‌‍​‌‌​‌‍‍‌‌​​‍‍‌‍​‌‌‍‌​‌‍‌‌‍‍‌‌‍‍​‍‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​‌‍‌​‍​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌‌‌‌​‌‍‌‌​​‌‍‌‌‌​​‍‌‌‍​​‍‌‌​‌‌‌‌‌‌‌‌​‌‌‌​​‌​‍‌‌​‌‌​‌​‍‌​‌​‌‌‍‍‌​‌​​‌‌‌‌‌​‌‌‍‍‍‌‍‌‌‌‌‌‌‌​‌‌​‌​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌​​‍‌‍​‍‌‌‌​‌‍‍‌​​‌‌​‌‌‍​‌‌‌​‌​‍‌‌‌​‌‌​​‌‍‌‍‌​‌​​‍‌‍‍‌‌​​​‍‌‌‌​‌‌​‌​‍‌‌​‍‌​‍‍‌‌​‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​​‍​‌‍‌‍‌‍‍‌‌‍‌‌‌‍​‌‍‌​‌‌​​‌‍​‌‌‌​‌‍‍​​‌‌‍​‍‌‍​‌‍‌‍‌‌‌‌​​‌‍‌​‌‌​‌​​‍‍​​​‍‍‌‍​‌‌​‍‌‌​‌‍‍‌‌‍​‌‍​‌‍‌‌​‍​‍‌‌ for nicer formatting

01 — From Design to Code

Welcome back. For the past three weeks I have been reading and planning. However, this week I've begun building. I now have a running codebase with a Gymnasium environment wrapping Basilisk, a PPO training pipeline, a validation script, and ran a training run with 200,000 timesteps on a fixed LEO-to-GEO transfer. However, its trajectories are not quite right. This post covers the most significant design decisions I made building the environment.

Source | Authors | Focus

envs/hohmann_env.py | HohmannTransferEnv, Gymnasium wrapper for Basilisk 6-DOF simulation | Core environment

envs/env_config.py | SpacecraftConfig, ThrusterConfig, OrbitConfig, EnvConfig dataclasses | Centralised configuration

utils/bsk_utils.py | build_simulation(), get_cartesian_state(), safe_rv2elem() | Basilisk integration layer

main_train.py / main_validate.py | PPO training loop (SB3) and post-training evaluation with trajectory plots | Entry points

tests/ (2 modules) | test_env_basic.py, Gymnasium API contract; test_env_oracle.py, physics validation | All passing

Debug training run | 200k timesteps, fixed LEO→GEO orbit, single environment, seed 42 | First results

02 — The MDP: State, Action, and a Design Choice Resolved

Blog 3 laid out the conceptual mapping from Basilisk to a Gymnasium environment. This week I locked in the concrete formulation. The three fundamental choices were what the agent observes, what it controls, and what Basilisk runs at per decision step.

The Observation Space

Last week's open question was whether to observe Cartesian state or Keplerian elements. I chose both. The 16-dimensional observation combines normalised ECI position and velocity with the orbital elements most useful for interpreting orbit shape, specifically semi-major axis a, eccentricity e, and the angular elements encoded as sines and cosines to avoid branch-cut discontinuities. Everything is scaled to roughly [-3, 3] for training stability.

# 16-dimensional normalised observation vector
# o = [ r/r_s, v/v_c, a/r_s, e, sin(nu), cos(nu), sin(omega), cos(omega),
# m/m_0, k/N, r_tgt/r_s, e_tgt ]
#
# r_s = 7 x R_earth ≈ 44 590 km (GEO ≈ 0.95, LEO ≈ 0.16 in normalised units)
# v_c = √(μ / r_target) ≈ 3.07 km/s at GEO
# k / N = step fraction — tells the agent how much time it has left

Notice how I included the target SMA and eccentricity directly in the observation. Without them, the policy network can never know what orbit it is trying to reach, which means it cannot generalize to the randomised-orbit curriculum I plan to enable later. The agent needs to read its goal, not memorise it.

The Action Space: Burn On-Times

I went with burn on-times rather than delta-v vectors. Last week's question was whether to have the agent command an abstract velocity change or a physical firing duration. I chose the latter for two reasons. First, it maps directly to Basilisk's THRArrayOnTimeCmdMsg with no additional controller layer. Second, every burn depletes the tank in a way that shows up in the m/m_0 observation term, and the agent is responsible for that tradeoff.

# Action → burn duration → propellant consumed
# t_burn = burn_frac × 0.9 × Δt_step
# 10% coast margin at end of every step for attitude to settle
# Δt_step = 1.5 x T_Hohmann / 40 (episode time / decision steps)
#
# Δm = m · (1 - exp(-Δv / v_e)), Δv = F · t_burn / m_0, v_e = Isp · g_0
# Tsiolkovsky, not linear dm = F·t/(Isp·g0). The linear approximation
# over-estimates propellant use by ~40-60% across a full Hohmann transfer

Basilisk's thrusterStateEffector in this configuration does not feed propellant depletion back into the orbital integrator. The propagated mass stays fixed at the wet mass. I track fuel consumption analytically in Python and use it for observations and the reward. The error from this constant-mass assumption accumulates to under 5% in delta-v over a full Hohmann, which is acceptable for the current stage.

The Decision Pace

Following Zavoli and Federici's discretisation, I use 40 fixed decision steps per episode. The episode length is set to 1.5 × T_Hohmann, giving the agent 50% more time than the ideal two-burn transfer requires. Each step therefore spans a substantial chunk of orbital time, roughly 400 seconds per step for the LEO-to-GEO problem. The agent is as a result making 40 coarse burn-or-coast decisions spread across the entire transfer window.

03 — Inside the Wrapper

The two essential components were managing Basilisk's simulation lifecycle inside Gymnasium's reset-step loop, and handling the attitude settling period at the start of every episode.

The Reset Problem

Basilisk has no checkpoint-restore mechanism, so every call to reset() constructs a new SimBaseClass from scratch, configures the spacecraft dynamics, thruster effector, navigation sensor, and attitude control modules, initialises the orbital state, and calls InitializeSimulation(). Memory leaks were a potential problem from doing this tens of thousands of times in a single Python session, but in practice the Python garbage collector handles the C++ resource cleanup correctly when the previous sim object is dereferenced.

# reset(), the core Basilisk lifecycle, every episode
self._bsk = build_simulation(sc_cfg, thr_cfg, r_init, v_init, sim_dt_s=2.0)
# Coast 600 s: attitude controller converges to prograde before first burn
scSim.ConfigureStopTime(sec2nano(600))
scSim.ExecuteSimulation()
# step(), one decision step
thrOnTimeMsg.write(OnTimeRequest=[t_burn_s], time=sim_time_ns)
scSim.ConfigureStopTime(sim_time_ns + sec2nano(dt_step))
scSim.ExecuteSimulation()
r_vec, v_vec = get_cartesian_state(bsk)
oe = orbitalMotion.rv2elem(MU_EARTH, r_vec, v_vec)

The Attitude Settling Coast

This was the detail I had not fully appreciated until I ran the first episode without it. The velocityPoint guidance law slews the spacecraft to align the thruster with the prograde direction, and the mrpFeedback controller closes the attitude loop with a settling time on the order of several hundred seconds. On an equatorial orbit, the initial attitude misalignment from the identity MRP can be close to 90 degrees. If the agent fires the thruster immediately at step zero, the burn is nearly perpendicular to the velocity vector and the orbit degrades rather than rises. The reference scenarioOrbitManeuverTH handled this with an explicit pre-burn coast; I replicated it with a fixed 600-second settling period before the first observation is returned. Although the agent never sees this coast, the physics depends on it.

In-Process vs. Out-of-Process Resolved: I went with in-process. The latency per step() is dominated by the Basilisk integrator itself, not inter-process communication overhead, so there was nothing to gain from a separate process. The main risk of memory accumulation across thousands of episode resets did not materialise in practice.

04 — Tests: All Green

Before running any training I wrote two test modules and used Gymnasium's built-in environment checker to catch API violations.

test_env_basic.py verifies the Gymnasium contract: observations returned by reset() and step() must lie within the declared observation_space bounds at every step, including degenerate states near the termination conditions. I added clipping as the final operation in _get_obs() to handle edge cases where orbital elements go singular near perfectly circular orbits.

test_env_oracle.py runs a hardcoded burn sequence designed to approximate the optimal Hohmann transfer with full burn fraction at perigee, coast through the transfer ellipse, and a second burn at apogee to circularise. It checks that the oracle agent reaches a final SMA within a reasonable fraction of the target without triggering any catastrophic termination. This is the physics validation test: if a hand-coded near-optimal agent cannot complete the transfer, the simulation is wrong.

pytest tests/ -v
======== test session results =======
tests/test_env_basic.py::TestReset::test_reset_returns_valid_obs PASSED [ 5%]
tests/test_env_basic.py::TestReset::test_reset_twice_does_not_crash PASSED [ 10%]
tests/test_env_basic.py::TestReset::test_reset_deterministic_with_seed PASSED [ 15%]
tests/test_env_basic.py::TestReset::test_reset_info_contains_timing PASSED [ 20%]
tests/test_env_basic.py::TestStep::test_step_zero_burn_does_not_crash PASSED [ 25%]
tests/test_env_basic.py::TestStep::test_step_full_burn_does_not_crash PASSED [ 30%]
tests/test_env_basic.py::TestStep::test_step_observation_in_bounds PASSED [ 35%]
tests/test_env_basic.py::TestStep::test_step_reward_is_nonpositive_for_coast PASSED [ 40%]
tests/test_env_basic.py::TestStep::test_step_reward_finite_for_burn PASSED [ 45%]
tests/test_env_basic.py::TestStep::test_mass_decreases_after_burn PASSED [ 50%]
tests/test_env_basic.py::TestStep::test_mass_unchanged_during_coast PASSED [ 55%]
tests/test_env_basic.py::TestEpisodeFlow::test_episode_terminates_or_truncates PASSED [ 60%]
tests/test_env_basic.py::TestEpisodeFlow::test_random_policy_survives_multiple_episodes PASSED [ 65%]
tests/test_env_basic.py::TestGymnasiumCompat::test_check_env PASSED [ 70%]
tests/test_env_basic.py::TestClose::test_close_is_idempotent PASSED [ 75%]
tests/test_env_basic.py::TestClose::test_reset_after_close_works PASSED [ 80%]
tests/test_env_oracle.py::TestOraclePhysics::test_oracle_achieves_near_geo_orbit PASSED [ 85%]
tests/test_env_oracle.py::TestOraclePhysics::test_no_burn_orbit_unchanged PASSED [ 90%]
tests/test_env_oracle.py::TestOraclePhysics::test_fuel_consumption_physically_reasonable PASSED [ 95%]
tests/test_env_oracle.py::TestOraclePrinting::test_print_oracle_episode_trace PASSED [100%]
======= 20 passed, 332 warnings in 3.60s =======

05 — The Debug Run: Erratic Rewards, Interesting Trajectories

With the environment passing all tests I ran the first training run: 200,000 timesteps, PPO with a 64x64 MLP, fixed orbit, single environment, seed 42. With a deliberately minimal configuration, I wanted a clean diagnostic signal before scaling up.

Parameter | Value

Training timesteps | 200k

Decision steps per episode | 40

Shared actor-critic MLP | 64x64

Turns out, the reward curve ended up being erratic, bouncing up and down with no clear trend across the full run. The training did not diverge catastrophically, but it definitely did not converge. The evaluation results were more interesting. Running 20 episodes with the best checkpoint and plotting orbital SMA over each decision step gives this:

Evaluation trajectories over 20 episodes. The red dashed line is the target orbit. The bottom panel shows the thruster command at each step. Most episodes overshoot the target significantly; Ep8 misses by under 1%.

Every episode starts at LEO and every trajectory rises. The spacecraft is consistently firing prograde and climbing out of LEO. A random policy would produce flat or chaotic SMA traces, so something is being learned.

However, SMA errors range from under 10% to over 300%, with most episodes overshooting the target altitude by a wide margin. Episode 8 rises, levels off near the target, and circularises correctly. It is surrounded by episodes that sail past the target entirely, reaching two or three times the target SMA before the step budget expires.

In the bottom panel we see the thruster command sequence follows the same coarse pattern across every episode, where burn fraction starts near 1.0 and tapers gradually down to near zero by step 40. The agent has converged on a monotonically decreasing burn schedule rather than the two distinct on-off pulses that the Hohmann transfer actually requires. It is learning a heuristic that sometimes produces a good answer and usually overshoots.

What the Optimal Policy Actually Looks Like: A Hohmann transfer uses two short impulsive burns: a prograde burn at perigee to raise apogee to the target altitude, a coast through the transfer ellipse, and then a second prograde burn at apogee to circularise. The optimal command sequence is pulsed: burn, stop, coast for ~20 steps, burn, stop. The debug policy is a continuous taper: burn hard early, bleed off gradually. It has learned to go up, but not when to stop.

06 — Diagnosing the Reward Function

The trajectories show that the reward function never penalizes overshooting. To see why, we look at the SMA progress term:

# Current per-step SMA reward (potential-based shaping)
# r_SMA = (a_k - a_{k-1}) / (a_target - a_init) × 50
#
# This gives positive reward whenever a increases,
# regardless of whether a is already above a_target.

What Needs to Change: The SMA progress term needs to reward proximity to the target, not just upward movement. A potential based on signed distance in (a, e) space, something like Φ(s) ∝ -|a - a_target| / a_target, would penalize overshooting and undershooting symmetrically. The terminal success bonus needs to come down to the same order of magnitude as the accumulated per-step signal, so the agent is not navigating in the dark between episode boundaries. And the eccentricity shaping should only be valuable near the right altitude, not universally.

Getting reward shaping right on a sparse, delayed-outcome problem is notoriously difficult, and a first-pass reward function almost never works. The good news is that the environment is correct, the physics is correct, and the policy is learning. Because of Episode 8's near-perfect trajectory, I now know the task is solvable with this setup.

The Erratic Reward Curve Explained: Episodes that happen to time out near the target altitude receive the massive +200 terminal bonus, spiking the curve upward, episodes that miss receive only per-step shaping, landing near zero, and the occasional crashes add a -50 dip. The critic cannot learn accurate value estimates under this level of variance, so I lose the advantage estimates that PPO's clipped objective trains on. Even with ε = 0.2 limiting per-update policy change, a noisy advantage signal means noisy gradient updates, which is visible in the policy's instability.

07 — Looking Ahead

Next week will start off on the reward function. I am not going to scale up to more timesteps or parallel environments until the reward curve shows a cleanly improving trend on this fixed, deterministic problem.

The specific changes I plan to try include replacing the one-sided SMA progress term with a distance-to-target potential that penalizes overshoot symmetrically with undershoot, rebalancing the terminal bonus to be proportional to the per-step reward scale, and possibly tying the eccentricity shaping weight to proximity to the target altitude so the agent is only rewarded for circularising in the right place. I also want to look at the value function loss in TensorBoard more carefully. If the critic is not learning an accurate value estimate, which seems likely given the variance, the advantage estimates will be garbage regardless of how well the actor is searching.

One thing that keeps sticking with me is that despite the chaos in the reward curve, the EvalCallback managed to find a checkpoint that produced Episode 8, which was a near-perfect Hohmann transfer with under 10% SMA error. That means the task is tractable with this exact environment and this exact policy architecture. Something in the reward signal, however noisy, was enough to push the policy toward the right answer at least once.

The agent has finally learned what it means to burn efficiently. Visit mazzola.dev​​​​‌‍​‍​‍‌‍‌​‍‌‍‍‌‌‍‌‌‍‍‌‌‍‍​‍​‍​‍‍​‍​‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍​‍​‍​​‍​‍‌‍‍​‌​‍‌‍‌‌‌‍‌‍​‍​‍​‍‍​‍​‍‌‍‍​‌‌​‌‌​‌​​‌​​‍‍​‍​‍‌‍​‍‌‍​‌‌​‌‍‍‌‌​​‍‍‌‍​‌‌‍‌​‌‍‌‌‍‍‌‌‍‍​‍‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​‌‍‌​‍​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌‌‌‌​‌‍‌‌​​‌‍‌‌‌​​‍‌‌‍​​‍‌‌​‌‌‌‌‌‌‌‌​‌‌‌​​‌​‍‌‌​‌‌​‌​‍‌​‌​‌‌‍‍‌​‌​​‌‌‌‌‌​‌‌‍‍‍‌‍‌‌‌‌‌‌‌​‌‌​‌​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌​​‍‌‍​‍‌‌‌​‌‍‍‌​​‌‌​‌‌‍​‌‌‌​‌​‍‌‌‌​‌‌​​‌‍‌‍‌​‌​​‍‌‍‍‌‌​​​‍‌‌‌​‌‌​‌​‍‌‌​‍‌​‍‍‌‌​‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​​‍​‌‍‌‍‌‍‍‌‌‍‌‌‌‍​‌‍‌​‌‌​​‌‍​‌‌‌​‌‍‍​​‌‌‍​‍‌‍​‌‍‌‍‌‌‌‌​​‌‍‌​‌‌​‌​​‍‍​‌​​‍‍‌‍​‌‌​‍‌‌​‌‍‍‌‌‍​‌‍​‌‍‌‌​‍​‍‌‌ for nicer formatting

01 — The Problem Left on the Table

Welcome back. Last week, the trained agent was reaching roughly the right altitude but was doing it wrong. Instead of two short impulsive burns (a prograde kick at perigee to elongate the orbit, a coast for roughly 20-30 steps, then a second kick at apogee to circularise), it was bleeding thrust continuously across the entire episode, ending with roughly the right orbit by sheer persistence.

The reward function was the culprit and two specific failures were embedded in its design. First, the SMA progress term rewarded the agent any time its orbit size increased, even past the target, so overshooting was costless. Second, the fuel penalty was calibrated as a "tiebreaker" rather than a real constraint, so the agent had almost no incentive to concentrate burns instead of spreading them out. This week I redesigned both.

Source | Description | Focus

envs/env_config.py | Reward parameters: fuel_reward_scale, decisive_penalty, ta_bonus_scale, success_bonus, SMA potential redesign | Reward redesign

envs/hohmann_env.py | _compute_reward(): symmetric potential, doubly-quadratic ta_bonus, altitude-weighted eccentricity | Implementation

main_train.py | clip_range 0.2→0.1, ent_coef 0.005→0.015, net_arch [128,128], log_std_init=-1.0 | PPO tuning

tests/test_reward_shaping.py | 381 lines of unit tests verifying every reward component in isolation | New test module

Training runs v4, v6, v7 | 10M timesteps each across 12–14 parallel environments, randomised orbit distributions | Results

02 — Redesigning the Reward

Reward shaping in reinforcement learning is a precise tool. Every term you add changes what the agent is maximising, and the changes interact with each other in non-obvious ways.

Change 1: Symmetric SMA Potential

The original progress term rewarded any increase in semi-major axis, regardless of whether the orbit was already larger than the target meaning that burning past the target was free. The fix is to replace the one-sided progress signal with a potential function that peaks at the target and falls off symmetrically in both directions:

# Before: one-sided, only rewards climbing
r_sma = (a_k - a_{k-1}) / (a_target - a_init) × 50
# After: symmetric, penalises both overshoot and undershoot
Φ(a) = -|a - a_target| / (a_target - a_init)
r_sma = [Φ(a_k) - Φ(a_{k-1})] × 50
# Moving toward the target is always rewarded, moving away is always penalised

This uses the potential-based reward shaping framework from Ng et al. (1999), which proves that any shaping reward of the form γΦ(s') − Φ(s) leaves the optimal policy of the original task unchanged.

Change 2: Fuel Penalty Made Real

I raised the fuel penalty scale from 10 to 75. At the old scale, burning inefficiently cost almost nothing relative to the SMA progress reward. A wasteful agent that drained its tank over 40 steps accumulated a total fuel penalty of about -8.7 across the episode. A successful, efficient agent earned around +246. The penalty was barely visible. At 75x, the same inefficient tank-drain costs about -65, while an efficient Hohmann brings in roughly +172 total.

Change 3: Decisiveness Penalty

Even with a stronger fuel penalty, the agent could still satisfy it by consistently burning at a moderate fraction rather than oscillating between full burns and coasting. To explicitly push toward binary coast-or-burn decisions, I added a per-step penalty on intermediate burn fractions:

r_decisive = -0.5 × burn_frac × (1 - burn_frac)
# burn_frac = 0 (coasting): penalty = 0
# burn_frac = 1 (full burn): penalty = 0
# burn_frac = 0.5 (maximum ambiguity): penalty = -0.125
# Accumulates to at most -5 over a full 40-step episode

The shape of this function (x(1−x)) is zero at both extremes and peaks at the midpoint. Any burn that is neither a full commitment nor a coast gets penalised. The total per-episode penalty is small by design, because the goal is to tip the cost-benefit balance instead of dominating the signal.

Change 4: Location-Aware Burn Bonus

This is the most mechanically interesting change. A Hohmann transfer only works if burns happen at specific points in the orbit: burn 1 at the lowest point (periapsis) to raise the opposite side, burn 2 at the highest point (apogee) to circularise. The agent needed to be told that where you burn matters as much as how long you burn.

I added a per-step bonus that rewards burning near the right orbital location at each phase of the transfer. During the first half of the transfer (while the orbit is still small) the bonus rewards burning near periapsis. After the midpoint, it rewards burning near apogee. The bonus is zero at the wrong location and zero during coasting, so it does not punish coast arcs.

# Phase 1 (orbit smaller than midpoint SMA): reward periapsis burns
r_ta = 4.0 × burn_frac² × max(0, cos ν)²
# Phase 2 (orbit larger than midpoint SMA): reward apogee burns
r_ta = 4.0 × burn_frac² × max(0, -cos ν)²
# ν is the true anomaly (position angle in the orbit). cos ν ≈ +1 at periapsis, ≈ -1 at apogee.

Both the spatial term (cos²ν) and the burn-fraction term (burn_frac²) are squared, making this a doubly-quadratic profile. The spatial squaring concentrates the effective firing window to roughly ±70 degrees around the node: at 60 degrees off-node the bonus is only 25% of peak, and at 80 degrees it is 3%. The burn-fraction squaring is the subtler fix. With a linear-in-burn-fraction formula, there was always some positive marginal reward for starting any burn near the correct location, even a tiny one. The squared version makes the gradient zero at burn_frac = 0: no marginal reward for starting a small burn at all. Only committing to a large burn earns meaningful location bonus. This is what breaks the "small burns everywhere near apogee" local optimum that survived the first version of this term.

Change 5: Altitude-Weighted Eccentricity Shaping

Executing burn 1 sends the eccentricity from near zero (circular LEO) up to roughly 0.5 to 0.7 (the elongated transfer ellipse). Under the old reward, this spike was penalised immediately, even though it is a required step toward the final circular orbit. The fix ties the eccentricity shaping weight to how far along the transfer the spacecraft is. Near LEO the eccentricity term is weighted near zero. Near the target altitude the weight rises to full scale. Circularising in the right place earns a strong reward; the unavoidable eccentricity during the climb costs almost nothing.

Why the Changes Work Together: Each change closes one hole in the previous reward. The symmetric SMA potential stops overshooting from being free. The higher fuel penalty makes efficiency a real constraint. The decisiveness term closes the "moderate burns everywhere" loophole. The doubly-quadratic location bonus closes the "small burns near the right place" loophole. The altitude-weighted eccentricity stops burn 1 from looking worse than coasting. Together they leave exactly one local optimum: two concentrated burns at the right orbital locations.

03 — A Physics Problem: Why Not GEO

There was a second problem with the original setup that is not about reward design at all. The previous training target was GEO (geostationary orbit, roughly 42,000 km altitude). While this is a physically meaningful goal, it turned out to be impossible for the agent to reach within the 5% accuracy tolerance given the current thruster configuration.

The issue is gravity losses. When a thruster fires for an extended period, the spacecraft curves through its orbit rather than pointing straight in one direction. Part of the thrust goes toward fighting gravity rather than building orbital energy. The longer the burn relative to the orbital period, the more energy is lost. For a LEO-to-GEO transfer with a 1000 N thruster on a 750 kg spacecraft, reaching the transfer orbit apogee requires burning for roughly 50% of the orbital period. At that fraction, gravity losses eat about 18% of the theoretical energy gain. The actual apogee ends up at around 34,500 km instead of the required 42,164 km, an 18% error, far outside the 5% tolerance.

Metric | Value

LEO orbital period consumed by burn 1 for GEO transfer | ~50%

Gravity loss on GEO transfer, 3.6× the 5% tolerance | 18%

Burn 1 duration for MEO (11,000 km) target, losses drop to 4.3% | ~16%

The fix is to change the target to Medium Earth Orbit (MEO) in the 8,000 to 11,000 km altitude band. Burn 1 for this target spans only about 16% of the orbital period. Gravity losses drop to roughly 4.3%, which falls within the 5% SMA tolerance. The physics of the transfer are identical: same two-burn structure, same Hohmann geometry, but the shorter burn time makes the problem achievable without changing any equipment. GEO remains the long-term goal; it requires either a more powerful thruster or a finer time resolution (more decision steps per episode) to reduce each individual burn's fraction of the orbital period.

04 — Steadying the Learner: PPO Tuning

Beyond the reward, I also changed several parameters of the PPO learning algorithm itself. These changes were motivated by observing what went wrong in the earlier runs.

Clip Range: 0.2 → 0.1

PPO's defining feature is a clipping constraint that limits how much the policy can change in a single update step. The clip range parameter controls this limit. The standard default of 0.2 allows policy changes of up to 20% per update, which is appropriate for most continuous control problems. For this task, the doubly-quadratic reward landscape is unusually sensitive near the second burn: a small change in burn timing or location produces a large change in reward, which can translate into an oversized gradient update. Setting the clip range to 0.1 forces the policy to take smaller, more conservative steps. This alone had a measurable effect on training stability, particularly in avoiding the early collapse seen in the v6 run.

Entropy Coefficient: 0.005 → 0.015

PPO adds a small bonus for policy entropy, a measure of how "spread out" the agent's action distribution is. High entropy means the agent is trying many different actions; low entropy means it has converged to a nearly deterministic strategy. In the v4 run the entropy eventually collapsed around 7.5 million steps: the policy became too confident in one particular (suboptimal) strategy, stopped exploring, and its performance cratered before slowly recovering. Raising the entropy coefficient from 0.005 to 0.015 keeps the agent exploring longer, preventing that collapse from happening so severely.

Network Width: 64×64 → 128×128

I widened both the actor and critic networks from 64 to 128 units per layer. The Hohmann task with randomised orbits presents 16 inputs and requires the agent to discover a phase-sensitive, two-burn strategy. The additional capacity gives the networks more room to represent the structured timing required.

Initial Action Distribution: log_std_init = -1.0

PPO for continuous actions represents the policy as a Gaussian distribution over each action dimension. The default initial standard deviation is 1.0, which produces a near-uniform spread across the full action range. Initialising with a smaller standard deviation (std ≈ 0.37) biases the agent toward more decisive actions from the start, reducing the amount of time spent in the noisy mid-range of the burn fraction during early training.

05 — Three Training Runs, Three Lessons

I ran three full training campaigns of 10 million timesteps each. The changes were introduced incrementally to isolate which factors were responsible for improvements and which introduced new instabilities.

Run | Configuration | Best Result

v4 | clip=0.2, ent=0.005, ta_scale=2.0 | Best: 184 @ 9.8M steps

v6 | clip=0.2, ent=0.015, ta_scale=4.0 | Best: 191 @ 9.9M steps

v7 | clip=0.1, ent=0.015, ta_scale=4.0 | Best: 192 @ 8.4M steps

v4: A Stable Run with a Bad Strategy

v4 introduced the full reward overhaul (symmetric SMA potential, higher fuel penalty, decisive_penalty, and an initial version of the location bonus at half scale) along with the wider 128x128 network. According to the training curve, the policy discovered a useful strategy around 1 million steps and held onto reward levels in the 175-183 range for the majority of training. The brief collapse around 7.5 million steps, followed by full recovery to 183 by 10 million, showed the policy converged too hard to one specific burn schedule, exploration collapsed, and it took another 1.5 million steps to shake loose.

Looking at the evaluation trajectories, v4 reaches the target orbit in every episode. However, we see the problem in the thruster command panel. Turns out, instead of a single decisive kick at apogee, v4 spreads its burn energy across the latter half of the episode with many small-to-medium burns.

v6: A Worse Collapse, a Better Recovery

v6 doubled the location bonus scale (from 2.0 to 4.0) and raised the entropy coefficient to 0.015 to address v4's collapse, but kept the clip range at 0.2. However, rather than staying stable, the training collapsed even earlier, at around 1.5 million steps, and spent the next 6 million steps stuck near a reward of 33. This is far worse than v4's brief late-run collapse.

A doubled location bonus creates steeper gradients near the optimal firing locations. With clip_range=0.2, a single gradient update could shift the policy enough to fall off the narrow reward peak and land in a flat, low-reward region. Once stuck there, the agent had to rediscover the correct strategy from scratch, which took 6 million steps. When v6 finally recovered around 7.5 million steps, it reached a higher peak (191) than v4 (184), confirming that the stronger location bonus is helpful once the policy finds the right territory.

v7: Controlled Steps, Cleaner Burns

v7 kept all of v6's changes but reduced the clip range from 0.2 to 0.1. The smaller clip range forced the policy to take more conservative steps, preventing the large early collapse that plagued v6. The training curve shows a much healthier progression where the policy reaches useful performance around 1.5 million steps and improves steadily to a peak of 192.4 at 8.4 million steps. There is a brief instability around 6.5 million steps, but the recovery is fast. The run ends with a late-training collapse near 9.5 million steps that the training run did not fully recover from. The best model checkpoint from 8.4 million steps is what matters here.

At 192 reward, v7's best model is the strongest result so far. More importantly, the shape of the trajectory has fundamentally changed.

06 — What the Trajectories Actually Look Like

The reward number is one thing. The trajectory is the real test. Here are the evaluation results from v4 and v7 side by side:

v4: Distributed Burns

v4: all episodes reach the target, but the bottom panel shows burn energy spread across steps 15-40. No clear two-burn structure.

v7: Two-Burn Hohmann

v7: burn 1 at steps 2-5 (periapsis), coast for 20 steps, burn 2 at steps 26-30 (apogee). All episodes land within 1.3% of the target SMA.

The bottom panel of the v7 plot showcases a concentrated burst of thrust at the start of the episode (when the spacecraft is near periapsis and prograde burns are most efficient), a long coast through the transfer ellipse, and then a second concentrated burst near the halfway point of the episode (when the spacecraft has reached apogee and a final prograde burn circularises the orbit). This is a Hohmann transfer, executed autonomously by a neural network.

The Oberth Effect in Practice: The reason burns at periapsis are more efficient than burns anywhere else in the orbit is the Oberth effect. It specifies that kinetic energy gained from a burn is proportional to F · v, where v is the spacecraft's current speed. Orbital speed peaks at periapsis and is lowest at apogee. Burning at periapsis for burn 1, where speed is highest, delivers maximum energy to the orbit.

On Eccentricity Errors: The v7 trajectories show SMA errors around 1.3% and eccentricities of roughly 0.01 to 0.03, well within the tolerances of 5% and 0.10. The small residual eccentricity comes from the fact that the agent is making a finite-duration burn rather than an instantaneous kick. Burning across a range of orbital positions inevitably creates a slightly non-circular final orbit. Even the analytical oracle policy produces similar residuals for the same reason.

07 — What Comes Next

With a working Hohmann policy in hand, the next step is to stress-test this one. I'll want to tackle the second part of my central research question for this project, which is whether a trained RL policy can handle uncertainties in real time. The plan is to make the model robust to perturbations by injecting noise directly into the observation and control channels during evaluation.

The long-term aim of the project is a policy that can autonomously correct for uncertainties in real time without ground intervention, which matters most for deep-space missions where communication delays make ground-loop corrections impossible. Testing the current Hohmann policy under perturbations is the first step toward understanding what that requires. If the policy collapses under moderate noise, that tells me something important about what the architecture or training curriculum needs to change before scaling to multi-gravity-assist scenarios.

There is also the late-training collapse in v7 to account for. The best checkpoint is from step 8.4 million, but the run ended in instability.

Next week: observation and control noise sweeps across the v7 best model and a diagnosis of where and how the policy breaks down under perturbation.

Building robustness through noise injection, then untangling why training kept breaking. Visit mazzola.dev​​​​‌‍​‍​‍‌‍‌​‍‌‍‍‌‌‍‌‌‍‍‌‌‍‍​‍​‍​‍‍​‍​‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍​‍​‍​​‍​‍‌‍‍​‌​‍‌‍‌‌‌‍‌‍​‍​‍​‍‍​‍​‍‌‍‍​‌‌​‌‌​‌​​‌​​‍‍​‍​‍‌‍​‍‌‍​‌‌​‌‍‍‌‌​​‍‍‌‍​‌‌‍‌​‌‍‌‌‍‍‌‌‍‍​‍‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​‌‍‌​‍​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌‌‌‌​‌‍‌‌​​‌‍‌‌‌​​‍‌‌‍​​‍‌‌​‌‌‌‌‌‌‌‌​‌‌‌​​‌​‍‌‌​‌‌​‌​‍‌​‌​‌‌‍‍‌​‌​​‌‌‌‌‌​‌‌‍‍‍‌‍‌‌‌‌‌‌‌​‌‌​‌​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌​​‍‌‍​‍‌‌‌​‌‍‍‌​​‌‌​‌‌‍​‌‌‌​‌​‍‌‌‌​‌‌​​‌‍‌‍‌​‌​​‍‌‍‍‌‌​​​‍‌‌‌​‌‌​‌​‍‌‌​‍‌​‍‍‌‌​‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​​‍​‌‍‌‍‌‍‍‌‌‍‌‌‌‍​‌‍‌​‌‌​​‌‍​‌‌‌​‌‍‍​​‌‌‍​‍‌‍​‌‍‌‍‌‌‌‌​​‌‍‌​‌‌​‌​​‍‍​‌‌​‍‍‌‍​‌‌​‍‌‌​‌‍‍‌‌‍​‌‍​‌‍‌‌​‍​‍‌‌ for nicer formatting

01 — Where We Left Off

Welcome back. Last week, the agent finally learned to execute a clean two-burn Hohmann transfer. The v7 model reached the target orbit with under 1.3% SMA error across all evaluation episodes, with two clearly concentrated burns at periapsis and apogee. That was a preliminary answer to the first half of the research question: can an RL agent learn the right orbital maneuver strategy? The second half is more demanding: does the policy hold up when the environment is no longer clean?

Real spacecraft navigate with imperfect components and unmodeled forces. This week I built a five-channel perturbation system to inject all of that into training, then ran four training campaigns to understand what it takes to produce a policy that is both accurate and robust.

Source | Description | Focus

envs/env_config.py | PerturbationConfig dataclass: 5 channels, default parameters, master enable switch | New config

envs/hohmann_env.py | _apply_actuation_noise(), _set_dynamics_noise(), _clear_dynamics_noise(), obs noise in _get_obs() | Perturbation implementation

main_train.py | --anneal-lr flag (default ON), lr/clip_range linear schedules, --clean flag | Training changes

main_montecarlo.py | Monte Carlo evaluation: N stochastic episodes, Wilson CI success rate, 4-panel output plot | New script

tests/test_perturbations.py | 22 tests covering each channel individually, combined, and interaction checks | New test module

Training runs v8, v9, v10, v11 | 10M timesteps each; 3 failed before v11; progressively tightened the training config | Results

02 — Five Ways to Break a Spacecraft

The perturbation system follows the framework from Federici et al. (2020), adapted for the Basilisk simulation. The idea is to inject realistic uncertainty into every layer of the spacecraft-environment interface so that the trained policy learns to handle noise as a permanent feature of the world rather than an occasional nuisance.

Channel 1: Observation Noise

Gaussian noise is added to the raw position and velocity measurements before any other computation happens. Everything the agent observes, including the derived semi-major axis, eccentricity, and true anomaly, is computed from those noisy measurements. Position noise is σ = 500 m and velocity noise is σ = 0.5 m/s, calibrated to ground-tracking accuracy for low Earth orbit. At a typical orbital radius of 7,000 km, 500 m is about 0.007% of the orbital radius, so it barely affects any single observation but creates enough cumulative uncertainty to matter.

Why Noise Must Propagate to Orbital Elements: An earlier version of the observation noise added Gaussian offsets directly to the normalized position and velocity terms in the observation vector, but computed the semi-major axis, eccentricity, and true anomaly from the true, noiseless state. This is not physically correct. A real navigation system derives orbital elements from sensor measurements, not from ground truth. Fixing this so that all 16 observation components derive from the same noisy position and velocity vectors makes the uncertainty self-consistent: the agent cannot know its true anomaly any more precisely than its position.

Channel 2: Actuation Magnitude Noise

The commanded burn fraction is multiplied by a factor drawn from N(1, σ=0.05), clipped to ±10σ. A 5% thruster magnitude uncertainty means a commanded full burn might deliver anywhere from 95% to 105% of the nominal thrust. Both the simulation and the fuel accounting use the perturbed value, so the agent cannot rely on its commands being executed exactly.

Channel 3: Actuation Pointing Noise

Even if the thruster fires at exactly the commanded duration, the attitude controller is not perfect. The burn fraction is additionally multiplied by cos(θ_err) where θ_err ~ N(0, 1°). At 1° the efficiency loss is under 0.02%, negligible for any single burn. Over 40 steps, and combined with the other channels, the agent experiences a consistent slight underperformance relative to what it commanded.

Channel 4: Dynamics Noise

Before each simulation step, a random inertial force is injected via a dedicated force effector in Basilisk and then zeroed immediately after the step completes. The force magnitude per axis is drawn from N(0, m · σ_dv / Δt) where σ_dv = 1 m/s, giving each step a random velocity perturbation of roughly 1 m/s per axis. This represents the aggregate of unmodeled higher-order gravity, solar radiation pressure, and atmospheric drag residuals. It is calibrated to be much larger than the real noise level so that training is actively challenging.

Channel 5: Missed Thrust Events (MTEs)

At each decision step, with some probability, the burn command is silently zeroed regardless of what the agent sent. This models a thruster command dropout or a discrete failure. The event counter is capped per episode. I started with a 5% probability (about 2 missed burns per episode expected) and later reduced it to 2% (about 0.8 missed burns per episode), which matches a realistic 99% thruster reliability rate. This change had a much larger effect on training than expected, and I will come back to it.

Channel 6: J2 Gravity

Optionally, the point-mass Earth model is replaced with a degree-2 spherical harmonics model using the GGM03S gravity coefficients. J2 is the dominant perturbation for low Earth orbit, responsible for nodal precession of roughly 1 km per orbit for the orbit altitudes used here. Basilisk downloads the coefficient file automatically on first use and caches it locally.

Parameter | Value

Position noise σ (ground tracking accuracy) | 500 m

Thruster magnitude uncertainty σ | 5%

Missed thrust probability per step (realistic rate) | 2%

03 — v8: Training Under Perturbations

v8 was the first run with all five perturbation channels active. It used the same PPO hyperparameters as v7 (ent_coef=0.015, clip_range=0.1) but trained in a perturbed environment where the agent never experiences clean physics. After training for 10 million steps across 14 parallel environments, the training curve climbed to a peak reward of roughly 185 and stayed there.

v7, the clean-physics model from last week, peaked at 192 at 8.4 million steps and then collapsed catastrophically at 9.85 million, ending at around 28. The perturbed v8 model not only recovered less than v7 at its peak but held on throughout with no collapse. At the time, the stability looked like a benefit of training under noise. The actual explanation emerged only after looking at v9 and v10.

v8 evaluation trajectories under perturbations. All episodes reach the target orbit, but the burn command panel shows energy spread across multiple steps rather than two sharp pulses.

The trajectories tell a more mixed story. v8 reaches the target orbit reliably, but the burn commands are smeared across many steps rather than concentrated at the two Hohmann nodes. Compared to v7's clean two-pulse structure, v8 looks like it is hedging, distributing thrust across a wider window rather than committing to a single decisive burn at each location.

Why v8 Hedged: With a 5% missed-thrust probability, the agent expected roughly 2 out of 40 burns to silently fail. That led to the agent betting that if any single burn might be zeroed, spreading the total delta-v across multiple steps means one missed burn hurts less. The policy is rational under its training distribution, but that distribution was set too aggressively. If we set a 99% reliability, the thruster should miss roughly 0.8 burns per episode, not 2.

I also built a Monte Carlo evaluation script this week to measure performance statistically rather than just plotting 20 episodes. The script runs N independently-seeded stochastic episodes, reports a success rate with a 95% Wilson confidence interval, and produces histograms of SMA error, fuel use, and reward. Under a 100-run evaluation, v8 shows strong success rates under perturbations, particularly compared to v7 which was never trained to handle noise.

100-run Monte Carlo evaluation comparing v8 (perturbation-trained) against v7 (clean-trained) under live perturbations. Success rate, SMA error distribution, fuel use, and reward distribution.

04 — v9 and v10: Two Failures, Two Diagnoses

First, the observation noise fix described above: noise now propagates from raw position and velocity into all derived orbital elements. This is strictly more correct than the previous formulation.

Second, the missed-thrust probability was reduced from 5% to 2%. At 5%, the agent expected roughly 2 missed burns across a 40-step episode. At 2%, it expects 0.8. The rate of 2% corresponds to 99% thruster reliability, which is a reasonable real-world figure. The 5% rate was chosen for training convenience rather than physical accuracy.

# Old defaults
mte_prob = 0.05 # 5% per step → ~2 missed burns expected per episode
mte_max = 3 # cap at 3 per episode
# New defaults
mte_prob = 0.02 # 2% per step → ~0.8 missed burns expected per episode
mte_max = 2 # cap at 2 per episode; P(≥3) < 1.5% at 2% per step

v9: The MTE-Entropy Coupling

v9 used the corrected obs noise and 2% MTE, with ent_coef=0.015 carried over from v8. The run peaked at 174 around 2.3 million steps and then collapsed to negative reward in about 50,000 steps, spending the remaining 7.7 million steps producing essentially nothing.

The root cause is a coupling between missed-thrust rate and entropy. In v8, every step had a 5% chance of randomly zeroing the action, which means even a near-deterministic policy experienced substantial reward variance from episode to episode. This variance acts like implicit entropy regularisation, keeping the policy stochastic at the level of training outcomes, supplementing the explicit ent_coef term. With 14 parallel environments and 5% MTE, roughly 87% of episodes experienced at least one action override, making the policy's effective exploration level much higher than ent_coef=0.015 alone would produce.

At 2% MTE, only about 55% of episodes have any missed thrust event. The training environment is roughly 2.5 times less stochastic. The policy converged to a sharp strategy much faster, reaching peak reward at 2.3 million steps instead of v8's 7.2 million, and that fast convergence led immediately to entropy death. Once the policy became too deterministic, small policy updates could not push it out of the local basin, and the basin it settled into turned out not to be the global optimum.

v10: The Constant Learning Rate Problem

The diagnosis of v9 suggested raising ent_coef to compensate for the reduced MTE stochasticity. v10 raised it from 0.015 to 0.025. The training curve rose to around 177 reward, but instead of collapsing, it began oscillating. Over the next 5 million steps, the reward flipped 12 times between the good basin at around 177 and the bad basin at around 25, eventually getting trapped permanently in the bad basin around 7.5 million steps.

This failure mode is different from v9's. v10 found the good basin repeatedly but could not stay there. The reason is that PPO updates the policy by a fixed step size regardless of where in training the policy currently is. At 1 million steps the policy has a lot to learn, and a step size of 3e-4 is appropriate. At 7 million steps the policy is near-optimal, and the same update magnitude is large enough to knock the policy off the narrow reward peak entirely. With ent_coef=0.025 actively trying to spread the action distribution, the combination of a non-decaying learning rate and a high entropy coefficient produced a machine that oscillated the policy in and out of the good basin every 600,000 steps until it eventually lost the path back.

v8 Was Not Stable for the Right Reasons: The retrospective analysis of v9 reveals that v8's stability was a side effect of its unrealistically high missed-thrust rate. With 5% MTE, the training environment provided implicit entropy regularisation through random action overrides. The policy could not converge too fast even if it tried, because 87% of episodes had at least one external disturbance to its action. v8 is a useful baseline model for robustness comparisons, but its spread-burn strategy is an artefact of 5% MTE rather than a proper solution.

05 — The Fix: Learning Rate Annealing

The root cause shared by both v9 and v10 is a learning rate that never decreases. In v9, a constant LR allowed fast convergence to a sharp, brittle policy. In v10, a constant LR applied update steps of equal size at step 7 million (where the policy was near-optimal) as at step 1 million (where the policy had little learned). The fix is linear annealing: both the learning rate and the PPO clip range decay linearly from their initial values to zero over the full training run.

# Annealing schedule (progress_remaining goes from 1.0 → 0.0)
lr_schedule = lambda f: 3e-4 × f # 3e-4 at step 0 → 0 at step 10M
clip_schedule = lambda f: 0.10 × f # 0.10 at step 0 → 0 at step 10M
# Effective LR at 7M steps: 3e-4 × (3/10) = 9e-5
# Same update that destabilised v10 is now 3× smaller

By 7 million steps, the effective learning rate is 9e-4, three times smaller than the value that caused v10's oscillations. Entropy noise of magnitude ent_coef=0.015 can no longer generate gradient steps large enough to cross the basin boundary. The clip range annealing has the same effect on the PPO policy ratio constraint, tightening the allowed policy change as the policy matures.

With LR annealing handling late-training stability, there is no longer a need to push ent_coef up to 0.025 as a substitute. Dropping back to 0.015 produces a sharper action distribution, which means more concentrated burns at the two optimal nodes rather than the diffuse spread-burn strategy that v10's high entropy coefficient enforced.

This is standard practice in deep RL, used in the original OpenAI Atari PPO runs. Without annealing, a near-optimal policy is continuously subjected to destabilising gradient updates, which makes late-training collapse or oscillation nearly inevitable for any environment with a sharp reward landscape.

ent_coef, MTE, and Learning Rate: These three variables are not independent. MTE determines how much implicit stochasticity the environment provides during training; ent_coef adds explicit entropy pressure; LR annealing controls how much any single update can move the policy. v8 worked because 5% MTE provided a lot of implicit stochasticity, making ent_coef=0.015 sufficient. v9 reduced MTE to 2% without raising ent_coef to compensate. v10 raised ent_coef to 0.025 without fixing the underlying LR problem. v11 reduces ent_coef back to 0.015 and lets LR annealing handle the stability that was previously provided by either a high MTE rate or a high entropy coefficient.

06 — v11: Annealed Training

v11 uses ent_coef=0.015, LR annealing on by default, 2% MTE, physically-consistent observation noise propagation, and all five perturbation channels active. As a result, we end up with v8 but with burns that are more concentrated at the two Hohmann nodes rather than diffused across the episode.

v11 evaluation trajectories under perturbations. The burn command panel shows tighter, more concentrated pulses than v8, closer to the two-burn Hohmann structure that v7 demonstrated in clean physics.

The key question v11 answers is whether a policy can simultaneously be robust to perturbations and execute a clean two-burn strategy. v7 achieved clean burns but on clean physics. v8 achieved robustness but at the cost of a hedged, diffuse burn structure. The combination of LR annealing and realistic 2% MTE gives the agent a training environment that is noisy without over-regularising the action distribution, and the trajectory plots show the tighter burns near the correct orbital locations.

Run | Result

v8 | Peak ~185 @ 7.2M — Stable, hedged burns

v9 | Peak 174 @ 2.3M — Collapsed to −6.6

v10 | Peaks ~177 — 12 oscillations, stuck

v11 | LR annealing ON — Concentrated burns

07 — What Comes Next

The Hohmann transfer problem for Earth orbit is essentially solved. The agent discovers the correct two-burn strategy, generalises across randomised orbit sizes within the MEO band, and remains robust under realistic sensor noise, thruster uncertainty, and unmodeled forces. There is still room to characterise v11 further through Monte Carlo runs across different noise levels, but the core problem is understood.

The natural next challenge is a longer-range transfer. An Earth-to-Mars Hohmann transfer shares the same mathematical structure with two impulsive burns, one to exit Earth's heliocentric orbit and one to match Mars's, with a coast arc spanning roughly 260 days in between; however, the simulation timescale changes from hours to months and the target is a moving body rather than a fixed orbital altitude.

The question I am most curious about is whether the same reward design, specifically the doubly-quadratic location bonus that teaches the agent to burn at periapsis and apogee, transfers directly to the interplanetary case. In the Earth-orbit setting, periapsis and apogee have fixed geometric meanings relative to the orbit. In the interplanetary case, the departure burn must happen at the right point in Earth's orbit around the Sun, and the arrival burn must happen at the right point in the transfer ellipse, which is determined by the relative positions of Earth and Mars at launch. Whether the reward shaping can guide the agent to discover the correct departure and arrival windows autonomously is the central experiment for next week.

Rebuilding the simulator around a moving planet, then running into a ceiling that tuning cannot fix. Visit mazzola.dev​​​​‌‍​‍​‍‌‍‌​‍‌‍‍‌‌‍‌‌‍‍‌‌‍‍​‍​‍​‍‍​‍​‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍​‍​‍​​‍​‍‌‍‍​‌​‍‌‍‌‌‌‍‌‍​‍​‍​‍‍​‍​‍‌‍‍​‌‌​‌‌​‌​​‌​​‍‍​‍​‍‌‍​‍‌‍​‌‌​‌‍‍‌‌​​‍‍‌‍​‌‌‍‌​‌‍‌‌‍‍‌‌‍‍​‍‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​‌‍‌​‍​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌‌‌‌​‌‍‌‌​​‌‍‌‌‌​​‍‌‌‍​​‍‌‌​‌‌‌‌‌‌‌‌​‌‌‌​​‌​‍‌‌​‌‌​‌​‍‌​‌​‌‌‍‍‌​‌​​‌‌‌‌‌​‌‌‍‍‍‌‍‌‌‌‌‌‌‌​‌‌​‌​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌​​‍‌‍​‍‌‌‌​‌‍‍‌​​‌‌​‌‌‍​‌‌‌​‌​‍‌‌‌​‌‌​​‌‍‌‍‌​‌​​‍‌‍‍‌‌​​​‍‌‌‌​‌‌​‌​‍‌‌​‍‌​‍‍‌‌​‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​​‍​‌‍‌‍‌‍‍‌‌‍‌‌‌‍​‌‍‌​‌‌​​‌‍​‌‌‌​‌‍‍​​‌‌‍​‍‌‍​‌‍‌‍‌‌‌‌​​‌‍‌​‌‌​‌​​‍‍​‌‍​‍‍‌‍​‌‌​‍‌‌​‌‍‍‌‌‍​‌‍​‌‍‌‌​‍​‍‌‌ for nicer formatting

01 — Where We Left Off

Welcome back. Last week closed out the Earth-orbit problem by producing a robust two-burn Hohmann policy that held up under observation noise and thruster uncertainty. The research question for the rest of the project is whether the same can scale to a longer and more demanding transfer, so this week I rebuilt the simulator around an Earth-to-Mars trajectory and ran the first training campaign on the new environment.

The target policy is structurally similar with a Trans-Mars Injection burn at the right moment in Earth's orbit, a coast arc of roughly 259 days, then a correction phase that steers the spacecraft toward Mars's actual position at arrival. The implementation is almost entirely new as almost nothing from the Hohmann environment survived unchanged.

Source | Description | Focus

envs/mars_transfer_env.py | New 952-line Gymnasium env: heliocentric frame, moving planets, two-phase episode, 18D observation | New environment

envs/env_config.py | InterplanetaryConfig dataclass: launch window geometry, departure/arrival planet parameters, phase step budgets | New config

utils/bsk_utils.py | build_interplanetary_simulation(): Sun-centred sim with Earth + Mars as gravity bodies driven by writable state messages | BSK builder

main_train.py, main_validate.py, main_montecarlo.py | Unified behind an --env flag so the Earth and Mars envs share the same training, evaluation, and MC pipeline | Shared scripts

tests/test_mars_env.py | 35 new tests covering the Gymnasium contract, heliocentric physics, reward components, and full-episode survival | New test module

runs/mars_v2, mars_v3, mars_v5, mars_v11 | Four training campaigns diagnosing departure-burn incentive design; v11 is the first converged model (~335 reward, 100% arrival) | Results

02 — A Different Gravity Problem

In the Hohmann env, Earth was front and center. The spacecraft orbited a single point mass and every position was measured from Earth's center. For an interplanetary transfer that frame is wrong. The spacecraft is in a parking orbit around Earth, but Earth itself is orbiting the Sun at 29.78 km/s, which is roughly four times the spacecraft's orbital speed around Earth. Ignoring that motion and placing Earth at a fixed point produces a simulation that is inconsistent with itself the moment the spacecraft leaves low orbit.

The new simulator uses a Sun-centred inertial frame. The Sun sits at the origin with the Earth and Mars modelled as gravity bodies orbiting it on circular paths. Basilisk exposes each planet's state through a SpicePlanetStateMsg, which is normally populated by a SPICE ephemeris file but can also be written directly. I take the second approach: at every simulation sub-step the environment computes each planet's position and velocity from its circular-orbit parameters and writes them into the message.

# Circular-orbit planet state at simulation time t
θ(t) = θ₀ + ωt, ω = v_helio / r_helio
r_planet(t) = r_helio [cos θ, sin θ, 0]
v_planet(t) = v_helio [−sin θ, cos θ, 0]

Refreshing planet positions at every sub-step rather than once per decision step matters during the departure phase. Earth moves about 8.9 million metres in a single 300-second decision step while the parking-orbit radius is 7 million metres, so a stale Earth position would warp the LEO dynamics by more than the radius of the orbit itself, so the planet-state write happens inside every 10-second burn sub-step and at every coast-chunk boundary.

The 104-Day Clock

Basilisk tracks simulation time as an int64 count of nanoseconds, which overflows at roughly 2⁶³ / 10⁹ ≈ 9.2 × 10⁹ seconds. That is more than two centuries in real time, but the internal API for time conversion uses a smaller accumulator that overflows at around 104 days. A full Earth-to-Mars episode needs roughly 337 days of simulated time (the 259-day Hohmann coast plus margin for the agent to arrive late), which is more than three times that limit.

The fix is to chunk the coast phases into segments no longer than 85 days each. At the end of each chunk the environment pulls the current spacecraft position and velocity out of the simulator, tears the Basilisk simulation down, and rebuilds it from the saved state with the time accumulator reset to zero. This is invisible to the RL agent as the observation it receives the moment before the rebuild and the one it receives the moment after are identical.

The Launch Window

In a real Hohmann-type transfer to Mars, the departure moment is chosen so that Mars will be at the spacecraft's aphelion when the spacecraft arrives there. The geometric condition is:

θ_Mars(0) = θ_Earth(0) + π − ω_Mars × T_Hohmann
# Mars has to "lag" by exactly the angle it will sweep during transit

With T_Hohmann ≈ 259 days and Mars's orbital rate of about 0.524°/day, Mars starts the episode roughly 44° ahead of Earth in its orbit, and by the time the spacecraft has traced half an ellipse around the Sun, Mars will have swept forward to the meeting point. Every episode resets to this ideal geometry, which keeps the problem clean for first-run training.

03 — Two Phases at Two Scales

The biggest challenge is the scale mismatch between the departure burn and the coast arc. The Trans-Mars Injection burn is a few minutes long and has to hit a precise point in a 90-minute LEO orbit. The coast lasts nine months. Using a single time resolution for both wastes computation if it is fine enough for the departure, and makes the departure unsolvable if it is coarse enough for the coast.

The episode is split into two phases with independent step durations. Phase 0 covers Earth departure at Δt = 300 s, giving 40 decision steps across roughly 3.3 hours. Phase 1 covers heliocentric transit and Mars approach at Δt ≈ 3.4 days, giving 100 decision steps across the full 336-day transfer. The agent sees a single continuous episode of 140 steps and a single 18-dimensional observation vector; only the internal integrator timescale changes.

Parameter | Value

Phase 0 decision step (LEO departure) | 300 s

Phase 1 decision step (heliocentric transit) | 3.4 d

Total decision steps per episode | 140

Departure: Geocentric Prograde

Phase 0 burns apply thrust along the spacecraft's velocity relative to Earth, using the geocentric frame instead of the heliocentric one. This is the direction that makes the Oberth effect work during escape from Earth's gravity well. The thruster force is raised to 4 kN (the Earth-orbit env used 1 kN) so that a full-duration Phase 0 burn delivers up to 1440 m/s of Δv. The theoretical TMI requirement is roughly 3526 m/s, so 2 to 3 aligned full burns suffice, which fits comfortably inside a single LEO orbit pass.

A lower thrust setting is actually worse for this problem. At 2 kN, the same burn magnitude needs 5 or more steps, and each burn stretches the orbit slightly so the next periapsis arrives later in the episode than the previous one. The agent never sees the same periapsis twice, making it hard to learn the timing pattern. At 4 kN all the necessary Δv can be delivered inside one periapsis window, so the agent only has to learn the single correct moment to fire.

The Alignment Observation: The observation vector includes a scalar tmi_align term, defined as max(0, v̂_geo · v̂_Earth). It is 1 when the spacecraft's geocentric velocity points directly along Earth's heliocentric velocity (the Oberth-optimal moment) and 0 when they are perpendicular or anti-parallel. This tells the agent directly how aligned the current moment is for an escape burn, without having to infer it from the raw position and velocity vectors. Every tested training run that omitted this term failed to discover the correct timing within 10 million steps.

Transit: Cross-Track Corrections

Once the spacecraft has escaped Earth's sphere of influence, the dynamics flip. The Sun becomes the primary gravitational body, Earth becomes a point moving away in the background, and the agent's remaining job is to steer toward Mars's arrival position. Burns here are small (capped at 14.5 s per step for roughly 89 m/s) and are applied in the cross-track direction: the component of the spacecraft-to-Mars vector perpendicular to the current heliocentric velocity.

Cross-track was chosen over a "burn toward Mars" heuristic because Mars is roughly 9° ahead of the spacecraft's velocity vector at arrival. A burn pointed at Mars would put 99% of its Δv along the velocity direction, adding orbital energy and pushing the aphelion past Mars. The cross-track projection rotates the orbit's orientation while leaving its energy untouched, which in principle lets the agent correct for departure direction errors without overshooting. I will come back to this choice in the last section.

04 — Reward Design for a 140-Step Episode

The Hohmann reward recipe from last month does not transfer directly. The agent has to solve two different sub-problems, and the signals that worked for a 40-step episode fire too sparsely to drive learning across a 140-step one. The Mars reward has four shaping terms plus a terminal component.

Alignment-Weighted Departure Bonus

Phase 0 needs its own positive signal because the distance-to-Mars shaping is essentially flat while the spacecraft is still inside Earth's sphere of influence. The bonus rewards mass burned during departure, weighted by alignment:

r_depart = (Δm / m₀) × (base + scale × tmi_align²) + r_fuel
# With base = 0, scale = 200, fuel_scale = 50:
# Net per burn = Δm/m₀ × (200 × align² − 50)
# Breakeven at align ≈ sqrt(50/200) = 0.50

The first three Mars training runs (v2, v3, v5) all failed at this term. v2 had no departure bonus and the agent simply coasted for every episode. v3 added a flat per-burn bonus and the agent burned at every step regardless of alignment, wasting fuel on retrograde kicks. v5 used a bonus that scaled linearly with alignment, which was still net-positive at alignments below 0.5 and produced the same unselective behaviour. v11 uses the quadratic alignment weight with a base of zero, which gives the agent the choice where burns at alignment below 0.50 cost more fuel than they earn, burns above 0.50 are profitable, and the optimal moment (alignment near 1.0) dominates everything else.

Why Random Initial Exploration Still Works: A breakeven at alignment 0.50 raises the question of whether random Phase 0 burns produce enough positive samples to bootstrap learning. Roughly half of Phase 0 timesteps have the spacecraft in the prograde half of its LEO orbit, so about 50% of random burns land above breakeven. That is enough for the value function to notice a correlation between high alignment and positive reward within the first 100,000 steps, which is what v11's training curve shows.

Distance-to-Mars Potential

The heliocentric shaping term is a potential-based function on the distance to Mars, normalised by the initial distance. It follows the Ng et al. (1999) framework from Week 5, which keeps the optimal policy unchanged while supplying a dense gradient. Closing the gap earns positive reward at every step, widening it costs the same reward, and the per-step magnitude is small enough not to overwhelm the departure bonus during Phase 0.

Graded Terminal Bonus

A 20 million kilometre success threshold (0.134 AU, chosen from the oracle's empirical 7-10 Mkm achievable closest approach) is reached in most episodes once the agent learns departure. For the episodes that fall outside the sphere, a graded bonus fills in the missing gradient:

r_prox = 100 × max(0, 1 − d_min / (2 × d_success))
# d_min: closest approach to Mars achieved during episode
# Active inside 40 Mkm, zero beyond

An episode that closes to 25 Mkm (just outside the success sphere) earns a proximity bonus of roughly 38, while one that closes to 50 Mkm earns zero. The success bonus itself is 200, so an actual arrival is still preferred over any proximity outcome by a wide margin. The purpose is to provide a distinguishable reward signal between "close miss" and "way off", which turns out to be the difference between the agent learning anything at all in Phase 1 and flatlining.

No Fuel-Depletion Termination

One of the subtler changes is that running out of fuel does not end the episode. With 4 kN thrust and the action range [0, 1], a random policy empties the tank within 7 steps, and earlier versions of the env then truncated the episode with a failure penalty. The problem is that this cuts off the agent's exposure to Phase 1 entirely before it has learned to be selective in Phase 0, so the distance-to-Mars gradient never reaches the policy. Now an empty tank simply means every subsequent burn command is silently ignored while the episode continues to run. The agent still coasts for all 336 days and receives the distance-based rewards, including the proximity grade, based on whatever trajectory it ended up on.

05 — mars_v11: First Converged Model

v11 trained for 10 million steps on 8 parallel environments with the PPO configuration that emerged from the Hohmann runs: batch_size=512, net_arch=[256,256], gamma=0.999, ent_coef=0.02, clip_range=0.2, linear LR annealing ON. The Mars env uses a larger network and a higher gamma than the Hohmann env because the episode is longer (140 steps vs 40) and the reward arrives further in the future.

mars_v11 evaluation across 12 episodes. Heliocentric trajectory (left), distance-to-Mars curve (top right), and burn commands by phase (bottom right, blue = Phase 0 departure, green = Phase 1 corrections). All episodes arrive within the 20 Mkm sphere.

The left panel shows the heliocentric plane. Earth starts at the bottom of the frame, the spacecraft traces the lower half of a transfer ellipse, and every episode closes on Mars on the right. The distance-to-Mars curve decays monotonically from 240 Mkm at departure to between 14 and 20 Mkm at arrival. The burn command panel shows the expected two-region structure: dense blue bars during the first 20 steps as the agent executes the TMI burn in pieces within a single LEO orbit pass, then 100 steps of near-zero green bars during heliocentric coast, with small corrective pulses near the end of the episode as the spacecraft approaches Mars.

The final reward across evaluation episodes averages around 335. For reference, the v2 run that simply coasted scored near zero, v3's unselective burning scored around 95, and v5 scored around 102 before plateauing. v11 is the first run whose policy is qualitatively doing the right thing.

06 — The 20 Mkm Wall

Every v11 evaluation episode arrives within the 20 Mkm sphere around Mars, but essentially none close to within 10 Mkm. The scatter is tight between 14 and 20. At first I tried to tune it by adjusting the Phase 1 correction cap, the proximity bonus scale, the entropy coefficient, and the network width. None of it moved the miss distance distribution by more than 1 Mkm.

The actual explanation is that the miss distance is set almost entirely by Phase 0, and Phase 1 corrections are provably incapable of removing the residual.

Cross-Track is Energy-Neutral

A cross-track burn changes the direction of the velocity vector without changing its magnitude, so the orbit's specific energy is preserved. Specific energy determines the semi-major axis, which determines the aphelion, which determines where the transfer ellipse meets Mars's orbit. Any miss distance that originates from an energy error in Phase 0 (overshooting or undershooting Mars's orbital radius) is outside the set of things cross-track burns can touch. In a 2D simulation this leaves only the in-plane orientation of the orbit as a controllable degree of freedom, and the 300-second departure step quantisation sets roughly 18 to 19 Mkm of unavoidable orientation error before any cross-track burn fires.

Why Not Burn Toward Mars Instead: Earlier versions of Phase 1 tried "burn along the spacecraft-to-Mars direction" rather than cross-track. However, Mars sits roughly 9° ahead of the spacecraft's velocity vector at arrival, so 99% of any Mars-pointed burn ends up along-track, inflates the orbital energy, and pushes the aphelion past Mars.

One-Dimensional Action Space

The agent commands a single scalar in [0, 1] and the environment picks the direction. For Phase 0 that is fine because geocentric prograde is the right direction for any departure burn. For Phase 1 it is the root of the problem. A real correction maneuver decomposes into three orthogonal components (radial, along-track, cross-track) and the optimal direction depends on where in the orbit the burn is fired. The agent has no way to express any of this. Every Phase 1 burn is forced into the single pre-selected cross-track direction regardless of what the orbit geometry actually calls for at that moment.

No Mars Orbit Insertion

The episode ends the moment the spacecraft crosses into the 20 Mkm sphere. In a real mission this is the start of the hard part: a capture burn at Mars periapsis that drops the spacecraft into a Mars-centred orbit. Without that terminal phase, the agent has no incentive to arrive at Mars with a specific velocity relative to Mars, so closing further than 20 Mkm earns nothing while spending fuel for the capture burn earns nothing either. The reward landscape has no gradient beyond the success sphere.

Hardcoded Circular Orbits

Earth and Mars are treated as idealised circular co-planar orbits at their mean heliocentric radii. Real orbits have non-zero eccentricity (Mars is at 0.0934, which moves its heliocentric distance by roughly 20 Mkm across a Mars year) and an inclination between them of roughly 1.85°. The 1.85° inclination alone contributes about 4.8 Mkm of out-of-plane distance at Mars's orbit radius, which is a substantial fraction of the 20 Mkm success threshold. The current env simulates neither.

Where the Error Actually Comes From: Integrating the oracle's trajectory from the moment of TMI completion shows that the residual miss distance is set to within 2 Mkm by the time the spacecraft exits Earth's sphere of influence. The 300-second Phase 0 decision step quantises the TMI burn timing by about half an LEO orbital period divided by 40, which is a few seconds of firing error, which propagates to millions of kilometres at Mars. Closing this gap requires continuous control over the burn timing; no amount of Phase 1 correction can recover it.

07 — What Comes Next

Next week's work is an environment rewrite rather than another hyperparameter sweep.

Phase 1 gets a multi-dimensional action space. Adding radial and along-track components alongside cross-track gives the agent actual control over the arrival geometry, and with that control the closest-approach distribution should narrow below 10 Mkm without changing anything else about the reward.

A Mars orbit insertion phase gets appended to the episode. The arrival condition becomes "captured into a Mars-centred orbit below some apoapsis" rather than "passed within 20 Mkm". This gives the agent a reason to arrive at Mars with a controlled relative velocity, which closes the loop on the transfer.

Planet positions get pulled from an ephemeris rather than from the circular-orbit approximation. Pointing the planet state messages at real SPICE positions fixes the orbit eccentricities and the mutual inclination.

The simulator itself has always been 3D; only the initial conditions and the target geometry were restricted to the x-y plane to simplify observation design. Once the planets are on real orbits, the out-of-plane component becomes meaningful and the observation vector has to carry it.

The interesting experiment for Week 8 is whether a policy trained on the ephemeris env can still produce the same clean TMI-plus-corrections structure that v11 found, or whether real orbit geometry produces training pathologies that the idealised version hid.

A new environment built around SPICE ephemerides, a 4D action space, and a Mars capture criterion, plus the two debugging weeks it triggered. Visit mazzola.dev​​​​‌‍​‍​‍‌‍‌​‍‌‍‍‌‌‍‌‌‍‍‌‌‍‍​‍​‍​‍‍​‍​‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍​‍​‍​​‍​‍‌‍‍​‌​‍‌‍‌‌‌‍‌‍​‍​‍​‍‍​‍​‍‌‍‍​‌‌​‌‌​‌​​‌​​‍‍​‍​‍‌‍​‍‌‍​‌‌​‌‍‍‌‌​​‍‍‌‍​‌‌‍‌​‌‍‌‌‍‍‌‌‍‍​‍‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​‌‍‌​‍​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌‌‌‌​‌‍‌‌​​‌‍‌‌‌​​‍‌‌‍​​‍‌‌​‌‌‌‌‌‌‌‌​‌‌‌​​‌​‍‌‌​‌‌​‌​‍‌​‌​‌‌‍‍‌​‌​​‌‌‌‌‌​‌‌‍‍‍‌‍‌‌‌‌‌‌‌​‌‌​‌​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌​​‍‌‍​‍‌‌‌​‌‍‍‌​​‌‌​‌‌‍​‌‌‌​‌​‍‌‌‌​‌‌​​‌‍‌‍‌​‌​​‍‌‍‍‌‌​​​‍‌‌‌​‌‌​‌​‍‌‌​‍‌​‍‍‌‌​‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​​‍​‌‍‌‍‌‍‍‌‌‍‌‌‌‍​‌‍‌​‌‌​​‌‍​‌‌‌​‌‍‍​​‌‌‍​‍‌‍​‌‍‌‍‌‌‌‌​​‌‍‌​‌‌​‌​​‍‍​‌​‍‍‌‍​‌‌​‍‌‌​‌‍‍‌‌‍​‌‍​‌‍‌‌​‍​‍‌‌ for nicer formatting

01 — Where Last Week Ended

Last week's mars_v11 model arrived inside a 20 Mkm sphere around Mars in every evaluation episode and could not get any closer. The audit traced the ceiling to the facts that the simulator was 2D, the action was a single scalar, Phase 1 corrections were hardcoded to a cross-track direction that cannot change orbital energy, the episode ended at flyby with no Mars Orbit Insertion, and the planets followed circular co-planar orbits that ignored real ephemeris geometry. None of these are things tuning can fix.

This week's work is the rewrite that addresses all five at once. The new environment lives at envs/mars_nasa_env.py with ENV_NAME = "mars_nasa". The old environment is untouched, and the registry picks up the new one automatically. The new env builds on real SPICE ephemerides, exposes a 4D continuous action, runs a 3-phase finite-state machine that ends in a Mars-centred capture orbit, and is supported by a hand-rolled Lambert solver, a B-plane reference trajectory, an analytic oracle, a behavior-cloning pretrainer, and a 3-stage curriculum.

Two training campaigns came out of the new environment this week. Both failed in interesting ways, and the diagnoses are the most useful thing I have to write down.

Source | Description | Focus

envs/mars_nasa_env.py | New 1100-line Gymnasium env: 3-phase FSM, 4D inertial-frame action, 36D observation, capture-orbit terminal success | New environment

envs/lambert_solver.py, envs/bplane.py, envs/reference_trajectory.py | Izzo Lambert solver, B-plane frame math, cached Lambert-based reference trajectory used by the observation | Trajectory tooling

utils/ephemeris.py + utils/bsk_utils.py::build_nasa_interplanetary_simulation | SPICE wrapper around the BSK-bundled DE430 kernels; new BSK builder using createSpiceInterface for live planet states | SPICE wiring

utils/oracles/mars_nasa_oracle.py | Lambert + scheduled TCMs + analytic MOI; 78/80 captures across 4 launch windows × 20 seeds | Oracle

scripts/generate_bc_dataset.py + scripts/train_bc.py | 1000 oracle episodes → MLP behavior-cloning checkpoint loaded into PPO actor before training | BC pipeline

utils/curriculum.py | 3-stage curriculum on parking radius and launch-window randomisation, gated by rolling capture rate via mp.Value | Curriculum

runs/nasa_v1, nasa_v2 | Two failed training campaigns; root causes diagnosed and documented (sections 06 and 07) | Results

02 — A 4D Action and a 36D Observation

The single-scalar action of the previous environment was the deepest architectural limitation of the old setup. The new action is four-dimensional and lives in the J2000 inertial frame:

a = (u_x, u_y, u_z, m), u ∈ [-1, 1]³, m ∈ [0, 1]
# Direction is normalised; magnitude scales the burn duration
û = u / max(‖u‖, ε), t_burn = m × t_max,phase

The agent picks any direction in 3D, and a separate scalar sets how long to fire. A spherical (α, β, m) form was rejected because it has a coordinate singularity at the poles where the gradients of a Gaussian PPO policy explode. The unit-vector form lets the policy learn ‖u‖ → 1 on its own without any constraint penalty. Pointing noise is layered on top by rotating û around a random axis using the existing perturbation infrastructure from Week 6.

The maximum burn duration is phase-dependent. Phase 0 (pre-TMI) keeps the 0.9 fraction from the old env, Phase 1 (cruise) uses 3.5e-5 (about 10 s/step) for trajectory correction maneuvers, and Phase 2 (MOI) uses 0.5 to allow the multi-hundred-m/s capture burn at Mars periapsis. The same action vector applies in all three phases; only the integrator window changes.

Observation

The observation grew from 18 floats to 36. The additions are heliocentric Keplerian elements (sin/cos of inclination, RAAN, argument of perigee), Mars-relative position and velocity (replaces the single distance scalar from before), B-plane targeting error against a Lambert reference (sections 04), reference-trajectory state deviation, a phase code, and the existing TMI alignment scalar gated to Phase 0 only.

Parameter | Value

4D | Continuous action (u_x, u_y, u_z, m)

36D | Observation, all components in [-3, 3]

3 | Latching phases: pre-TMI, cruise, MOI

Mars Orbit Insertion

The success condition is finally the right one. Once the spacecraft enters Mars's sphere of influence (about 577 Mkm) and is approaching, the env transitions to Phase 2. From then on, the Mars-relative orbital state is converted to Keplerian elements at every step. Capture is declared when:

a > 0, 0.10 < e < 0.98
r_p = a(1 - e) - R_Mars ∈ [250 km, 10,000 km]
r_a = a(1 + e) < 5 × R_Mars,SOI
# All four conditions held for 3 consecutive steps

The 3-step latch prevents a glancing pass from being mistaken for capture. A graded capture-quality score q_moi ∈ [0, 1] is computed alongside the boolean criterion as a product of Gaussian kernels on (e, r_p) relative to a target of e = 0.85, r_p = 500 km. The terminal reward includes both: +500 for the boolean success plus +200 × q_moi for the graded part, so the agent has a smooth gradient to chase even before it lands inside the strict box.

03 — Real Orbits via SPICE

Hardcoded circular planet orbits go away. Earth and Mars are now driven by JPL's DE430 ephemeris through SPICE, which Basilisk happens to ship as a support file at ~/.cache/bsk_support_data/. utils/ephemeris.py is a thin spiceypy wrapper that furnishes the kernels once at process start and exposes get_body_state(body, et) and utc_to_et. The new BSK builder, build_nasa_interplanetary_simulation, swaps the hand-written planet-state messages for gravFactory.createSpiceInterface(...) with epochInMsg=True, which lets BSK push live planet states into the simulation as time advances.

The episode now picks a random launch window from a list of four real Mars opportunities:

Launch Window | Hyperbolic Excess Velocity

2024-10-15 | v∞ ≈ 4.32 km/s

2026-11-11 | v∞ ≈ 3.49 km/s

2028-12-15 | v∞ ≈ 3.61 km/s

2033-05-05 | v∞ ≈ 3.15 km/s

An earlier candidate of January 2031 was dropped after the oracle could not retarget through it: the Lambert solution diverged from a 374 m/s correction to 18.2 km/s during Phase 1 because the Earth-Mars geometry that month produced a transfer arc too sensitive to the departure timing. Replacing it with May 2033 gave a clean Lambert solve and an oracle that captures successfully.

Because of real ephemerides, none of the simplifying tricks from the old env survive. Mars's heliocentric distance varies by ~20 Mkm across its year because of its 0.0934 eccentricity, the orbit plane is inclined 1.85° relative to Earth's, and the launch geometry shifts substantially between 2024 and 2033. The agent has to learn to handle that variability rather than memorize one fixed scenario.

04 — Lambert, B-plane, and the Oracle

To give the agent a stable target, the observation includes the deviation from a Lambert reference trajectory, which is the unique two-body solution that connects the Earth's position at departure to the Mars's position at arrival in exactly the chosen flight time.

The Lambert solver lives in envs/lambert_solver.py, a hand-rolled implementation of the Izzo 2015 algorithm. The solver takes the two heliocentric position vectors, the time of flight, and μ_Sun, and returns the velocity at each end of the arc. envs/reference_trajectory.py wraps it in a 200-waypoint sampled trajectory plus a B-plane target, with a disk cache keyed on the rounded epoch and parking altitude.

B-plane Targeting

The B-plane is a 2D plane perpendicular to the asymptotic incoming velocity vector v∞. NASA mission designers use it because hyperbolic flybys map cleanly onto two scalar coordinates (B · T̂ and B · R̂), and the periapsis altitude and inclination of the resulting capture orbit depend almost linearly on those two numbers. envs/bplane.py computes them from a Mars-relative state vector during cruise, and the observation carries the difference between the current and reference B-plane coordinates as two normalised scalars.

Why a Reference at All: Without a reference, the agent has to learn arrival geometry from raw position and velocity, which scrambles every time the launch window changes. With the Lambert reference in the observation, the agent's job becomes "follow the reference, and do better at the end", which is well-defined across all four launch windows. The reference weight is annealed to zero by 5M training steps so the policy can outperform Lambert during the final stretch of training.

The Oracle

The oracle (utils/oracles/mars_nasa_oracle.py) is a non-RL expert that captures into Mars orbit using analytic methods only. It exists for two reasons: to validate that the env is solvable in principle, and to provide demonstrations for behavior cloning. Phase 0 fires the Lambert Δv in up to 3 burns at peak TMI alignment. Phase 1 schedules trajectory correction maneuvers at t = 5d, 60d, and 180d using a primer-vector first approximation δv ≈ -r_err / t_remaining. Phase 2 computes the analytic MOI Δv from the hyperbolic state at Mars periapsis and splits it across up to 5 sub-burns.

The oracle gate is the most important test in the project. If the oracle cannot capture reliably, no amount of RL is going to. After two iterations on the Phase 1 cap and the B-plane target offset, the gate test now achieves 78/80 captures across the 4 launch windows × 20 random seeds:

# Oracle gate test — 4 windows × 20 seeds
# 2024-10-15: 20/20
# 2026-11-11: 20/20
# 2028-12-15: 18/20
# 2033-05-05: 20/20
capture rate = 78 / 80 = 97.5%

The Phase 1 dv cap had to grow from 1500 m/s to 5000 m/s before the 2028 window stopped failing, and the B-plane target offset had to drop to zero (the original geometry-dependent offset of 5 × 10⁷ m was actively pushing the spacecraft away from Mars in some configurations).

05 — Bootstrapping with BC + Curriculum

A 36D observation, a 4D action, and four launch windows is enough state space that PPO from random initialisation does not converge in any reasonable time on this problem. The training stack adds two pieces around the standard PPO loop.

Behavior Cloning Warm Start

scripts/generate_bc_dataset.py runs the oracle through 1000 randomised episodes and saves (s_t, a_t) pairs to a single .npz file. scripts/train_bc.py trains a [256, 256] MLP with MSE loss on the actions for 50 epochs at lr = 10⁻⁴. The output checkpoint is loaded into the SB3 PPO actor before the training loop starts:

# In main_train.py, when --bc-init is given:
model = PPO("MlpPolicy", env, policy_kwargs=dict(net_arch=[256,256]), ...)
bc_state = torch.load("bc_checkpoints/bc_policy_v1.pt")
model.policy.load_state_dict(bc_state, strict=False)
# strict=False keeps the value head + log_std random

The strict=False matters because only the actor MLP weights match the BC checkpoint, so the value function and the action distribution scale start from PPO's defaults. Loading the value head would corrupt the on-policy advantage estimates from the very first rollout.

Three-Stage Curriculum

utils/curriculum.py manages a 3-stage progression keyed on rolling capture rate via a multiprocessing.Value shared across vectorized environments:

# Stage 1 — fixed 2026-11-11 epoch, fixed parking radius, no perturbations
# Stage 2 — random parking radius (6800–8000 km), still fixed epoch
# Stage 3 — random epoch from all 4 windows, perturbations on

However, the mp.Value object wraps an mmap region under the hood, and the default SubprocVecEnv start method (forkserver on Linux) cannot pickle mmap across the process boundary. Switching to start_method="fork" when the curriculum is enabled fixes it because forked children inherit the memory directly without serialising. main_train.py sets the start method automatically when --curriculum is passed.

06 — nasa_v1: A Reward Landscape Disaster

The first training run was nasa_v1: 10M steps, BC initialisation, curriculum off (to isolate baseline behavior). Eval reward stayed between -605 and -80,000 across the entire run.

1. Unbounded Reference Penalty

The Lambert reference deviation was added as a potential-based shaping term:

r_ref = -w_ref × (Δ‖r_sc - r_ref‖) / (10⁸ m)

Combining w_ref = 20 and any divergent trajectory, the cumulative penalty over a 440-step episode reached -50,000 long before the spacecraft did anything else interesting. The "potential-based" framing makes the term telescope to a constant in the limit, but in practice the spacecraft never returns to the reference, so the telescoping never closes. Setting w_ref = 0 for nasa_v2 disabled the term entirely. Lambert deviation is still in the observation; however, it's not directly penalized by the agent.

2. B-plane PBRS Init Bug

The B-plane error term carried a similar PBRS structure with a previous-step buffer initialised to zero. The first time the spacecraft entered the B-plane evaluation region (Phase 1 onward), the previous error was zero and the current error was several million metres, so the per-step delta was a one-shot reward of around -1000. With w_bplane = 10⁻⁶ in the original config the magnitude was small enough not to matter, but combined with r_ref it pushed several diagnostic episodes into the truly catastrophic range. Setting w_bplane = 0 removed this entirely.

3. Earth Re-impact False Positives

The episode terminated with a -500 penalty if dist_earth < 6528 km after step 40. The intent was to catch a spacecraft that crashed back into Earth after a botched escape. However, random Phase 0 burns from any starting point regularly nudged the LEO perigee below 6528 km altitude (LEO is at 7000 km radius, parking altitude near 600 km, so any retrograde burn fragment crosses the threshold). Even a "do nothing" policy that just orbited Earth got a -500 penalty around step 41 in a meaningful fraction of episodes.

The fix gates the impact check on a latched _has_left_earth flag, set to true the first time the spacecraft passes 50,000 km from Earth. Re-impact is only flagged for spacecraft that actually escaped first.

Reward Distribution Before vs After: Same env, same launch window, same seed, three policies for comparison:

# Before fix
# Oracle: +146 (should be high, but was being beaten by penalties)
# Coast: -598 (should be near zero, but penalties dominated)
# 5 random: ≈ -50,000 each (catastrophic)
# After fix
# Oracle: +621
# Coast: -6
# 5 random: -22 to -85

The signal-to-noise gap between success and exploration went from less than 1 reward unit to about 700, which is the gradient PPO needs to climb.

07 — nasa_v2: The log_std Runaway

With the reward fixes in place, nasa_v2 ran for 10M steps with BC initialisation and the 3-stage curriculum enabled. The policy improved cleanly during the first few million steps, peaked at an evaluation reward of about 154, then plateaued at around 75 for the remaining run. The agent didn't achieve any captures and every evaluation episode terminated by step-budget timeout at the 440-step limit.

nasa_v2 evaluation in the heliocentric (J2000) plane. The policy completes a credible transfer arc and the spacecraft passes within 50 Mkm of Mars (well inside the 577 Mkm SOI), but the burn-command panel shows scattered Phase 1 firing and no concentrated MOI burn at periapsis. No capture in any evaluation episode.

The deterministic mean policy (with action noise stripped out) was actually fine and re-running evaluation in deterministic mode gave +154 reward and brought the spacecraft to within 600 Mkm of Mars at closest approach. The policy knew what to do; the stochastic version it actually trained with was burying every decision under noise.

The Diagnosis

Tensorboard's train/std trace showed the action distribution width grew monotonically from 0.37 to 2.04 over 10M steps. Inspecting the saved model directly, the per-dimension log_std was [1.16, 2.98, 1.99, 1.13]. The action box is [-1, 1] on each dimension, so a Gaussian with σ = 2 effectively samples uniformly from the box and the deterministic mean carries no information.

The cause is the entropy coefficient. The Hohmann env uses ent_coef=0.01, and that was inherited unchanged into the Mars NASA env's training defaults. Entropy gradient scales with the dimensionality of the action though, and a 4D Gaussian's entropy gradient is roughly four times larger per step than a 1D Gaussian's. The same coefficient that kept Hohmann's policy healthy was actively inflating Mars's log_std at every update. By the time training was halfway through, KL ≈ 0 and no further policy learning was happening.

Curriculum Never Promoted Either: A second issue compounded the first. The curriculum's promotion threshold was 60% capture rate to advance from Stage 1, and zero captures means the agent stayed in Stage 1 for all 10M steps. Stages 2 and 3 (parking radius randomisation, all four launch windows, perturbations) never trained. Lowering the thresholds to 5% and 20% along with a hard 3M-step timestep fallback fixes the curriculum, but the log_std runaway has to be fixed first or even Stage 1 will not get anywhere.

The Fix

# 1. Per-env entropy coefficient default (envs/mars_nasa_env.py)
TRAIN_DEFAULTS["ent_coef"] = 0.002
# 2. New CLI flag pinning the initial action distribution width
python main_train.py --env mars_nasa --log-std-init -2.0 ...
# 3. Reset log_std at the start of a resume run, undoing the runaway
python main_train.py --resume runs/nasa_v2/best_model/best_model
--reset-log-std -2.0 ...

Resuming nasa_v2's best deterministic policy with log σ reset to -2.0 (giving σ ≈ 0.135, well inside the action box) and entropy coefficient dropped to 0.002 should let the policy continue improving from where the deterministic mean already is, rather than restarting from scratch.

08 — What Comes Next

The plan for next week is the nasa_v3 result and the tuning that follows it, which will be based on these factors in the eval logs.

Curriculum advancement. With promotion thresholds at 5% capture, Stage 1 should hand off within the first few million steps, and the agent should start seeing parking-radius randomisation. If the agent never gets above 5% even on Stage 1, the entropy fix did not address the root cause and the log_std is still leaking through somewhere else.

Capture rate by launch window. The oracle scores 100% on three of the four windows and 90% on 2028. If the trained agent's per-window split mirrors that, the policy is generalising; if it collapses to one window only, BC is being undone too quickly.

B-plane error distribution. The whole point of the reference trajectory is that the agent learns to track and improve on it. A trained policy with median B-plane error below 5000 km is the release criterion for this env. Anything above 50,000 km means the BC weights have been lost.

Beyond nasa_v3, the next environmental improvement is enabling Mars's GGM2B spherical-harmonic gravity model in Phase 2 only. Point-mass Mars is fine for cruise and approach, but the capture orbit elements are sensitive to oblateness during the periapsis pass. Enabling it everywhere doubles the integrator cost, so gating it to Phase 2 is a cheaper compromise. The perturbation pipeline from Week 6 also needs to be plumbed through the new env, which currently ignores its own cfg.perturb in Stage 3 because build_nasa_interplanetary_simulation does not accept a perturbation config.

The interesting question for next week is whether 4D action plus capture criterion plus real ephemerides actually delivers what was missing from mars_v11, or whether the next ceiling is lurking one layer deeper.

A reverse curriculum, a Phase 2 reward overhaul, the first 5/5 deterministic Mars captures, multi-offset training that scaled to 20/20, and the wall that stopped extension cold. Visit mazzola.dev​​​​‌‍​‍​‍‌‍‌​‍‌‍‍‌‌‍‌‌‍‍‌‌‍‍​‍​‍​‍‍​‍​‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍​‍​‍​​‍​‍‌‍‍​‌​‍‌‍‌‌‌‍‌‍​‍​‍​‍‍​‍​‍‌‍‍​‌‌​‌‌​‌​​‌​​‍‍​‍​‍‌‍​‍‌‍​‌‌​‌‍‍‌‌​​‍‍‌‍​‌‌‍‌​‌‍‌‌‍‍‌‌‍‍​‍‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​‌‍‌​‍​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌‌‌‌​‌‍‌‌​​‌‍‌‌‌​​‍‌‌‍​​‍‌‌​‌‌‌‌‌‌‌‌​‌‌‌​​‌​‍‌‌​‌‌​‌​‍‌​‌​‌‌‍‍‌​‌​​‌‌‌‌‌​‌‌‍‍‍‌‍‌‌‌‌‌‌‌​‌‌​‌​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌​​‍‌‍​‍‌‌‌​‌‍‍‌​​‌‌​‌‌‍​‌‌‌​‌​‍‌‌‌​‌‌​​‌‍‌‍‌​‌​​‍‌‍‍‌‌​​​‍‌‌‌​‌‌​‌​‍‌‌​‍‌​‍‍‌‌​‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​​‍​‌‍‌‍‌‍‍‌‌‍‌‌‌‍​‌‍‌​‌‌​​‌‍​‌‌‌​‌‍‍​​‌‌‍​‍‌‍​‌‍‌‍‌‌‌‌​​‌‍‌​‌‌​‌​​‍‍​‍​​‍‍‌‍​‌‌​‍‌‌​‌‍‍‌‌‍​‌‍​‌‍‌‌​‍​‍‌‌ for nicer formatting

01 — Where Last Week Ended

Last week's nasa_v2 plateaued at evaluation reward 75 after 10M steps with zero captures, traced to a 4D log_std runaway plus a curriculum that never advanced past Stage 1. The fix list was a per-environment ent_coef = 0.002, a CLI flag for initial log_std, and a --reset-log-std escape hatch for resumes.

nasa_v3 trained for another 5M steps from nasa_v2's best checkpoint with ent_coef=0.002 and log_std reset to -2.0. Eval reward climbed cleanly to about 200, then stalled without captures. The mean policy steered the spacecraft to a closest approach of around 250 Mkm, well inside Mars's sphere of influence (577 Mkm), but it never produced a sustained anti-velocity burn at Mars periapsis. Several follow-up runs (nasa_v4 through nasa_v22) cycled through tweaks to alignment bonuses, distance shaping, fuel weights, and BC dataset composition. All of them plateaued in the same place: somewhere between +250 and +450 evaluation reward, with zero captures across roughly 60M cumulative training steps.

That much wasted compute against the same plateau is a sign that the gradient that PPO is climbing does not lead to capture. There is a partial-brake harvesting path worth +400, and there is a full-capture path worth over +1500. The first is reachable from random exploration; the second is not. This week's work tries to close that gap two different ways.

Source | Description | Focus

scripts/precompute_soi_entry.py + bc_checkpoints/soi_entry_states.npz | Cache oracle Phase 0+1 trajectory waypoints (SOI edge, plus 18 cruise offsets back to 120 days) so the env can spawn the spacecraft mid-flight | Reverse curriculum data

envs/mars_nasa_env.py + envs/env_config.py (PHASE_MOI rewards) | w_moi_burn_align, w_mars_energy linear, w_q_moi_dense, w_bplane re-enabled with tanh + lazy init; terminal capture reward doubled | Phase 2 surgery

main_train.py (--p2-only, --cruise-offset-*, --epoch-weights, --clamp-log-std) | Reverse-curriculum CLI, weighted multi-offset sampling, hard log_std clamp box | Training infrastructure

runs/nasa_v23, runs/nasa_v24 | First deterministic capture at Mars: 5/5 success, +1762 reward, window 0 only, SOI-edge spawn | First capture

runs/nasa_v31, runs/nasa_v32 | Multi-offset training across {T=0, T-7d, T-30d, T-40d}: 20/20 captures from a single policy | Multi-offset breakthrough

runs/nasa_v33 ... v37 | Five attempts to extend: reverse curriculum past T-40d, multi-window from v32, window-1 specialist, multi-window from BC. All failed. | The wall

utils/controllers/mars_hierarchical.py + main_hierarchical_validate.py | Phase-dispatching controller: oracle for Phase 0 + early cruise, v32 for the last 40 d + capture. Ships next week. | Next week's plan

02 — Skip the Hard Part

A reverse curriculum trains the last phase of a problem first. For Mars, the last phase is the brake at periapsis: a multi-hundred-m/s anti-velocity burn that converts the hyperbolic approach into a bound orbit. If the agent could be dropped at Mars's sphere of influence with a credible incoming hyperbola already in hand, capture would be the only thing it ever had to learn.

The oracle already constructs that approach state every time it runs, because Phase 0 (TMI) and Phase 1 (cruise + TCMs) come before Phase 2 in the same episode. So scripts/precompute_soi_entry.py runs the oracle once per launch window, snapshots the heliocentric (r, v) at the moment the spacecraft crosses Mars's SOI, and dumps the result to bc_checkpoints/soi_entry_states.npz. The script also samples 18 cruise waypoints earlier in the trajectory (from 1 hour back through 120 days), so the env can spawn at any chosen time-before-SOI.

InterplanetaryConfigV2.start_at_soi_edge flips on the new spawn path. When set, reset() looks up the cached waypoint matching the chosen launch window and cruise offset, writes those numbers into BSK as the initial state, and sets self._phase = PHASE_MOI (or PHASE_CRUISE if the offset is non-zero). The episode budget is shared with the full-trajectory case (440 steps), but with no Phase 0 or early Phase 1 to spend it on, the agent gets ~200 Phase-2 steps inside the SOI to commit to a brake.

Cache Beats Synthesis: An earlier attempt placed the spacecraft on an analytic two-body Mars hyperbola with a tunable B-plane offset. Sun-perturbed dynamics during real cruise put the spacecraft on a different geometry from any clean two-body flyby, so the analytic spawn was off-distribution from the moment Phase 2 began. Caching the oracle's actual SOI state instead of synthesising one fixed it.

03 — Phase 2 Reward Surgery

Reverse curriculum changes the spawn distribution but says nothing about whether the reward surface inside Phase 2 has a path to capture. Inspecting the v17–v22 plateau showed it does not. The agent's locally optimal behaviour at SOI edge is a brief partial brake that nudges the trajectory closer to Mars, harvests the distance-shaping bonus, and lets the periapsis pass without committing the rest of the dv. Several reward terms had to be added or rebalanced to make a full capture strictly better than the harvest path.

MOI Burn Alignment

The Phase 0 alignment bonus pays the agent for firing in the prograde direction relative to Earth's velocity. The Phase 2 analogue rewards firing anti-velocity relative to Mars:

moi_align = max(0, -û · v̂_rel)
r_moi_align = (dm / m_0) × w_moi_align × moi_align²
# w_moi_align = 400; with w_fuel = 100, perfect alignment is +300·(dm/m0)

This is the gesture an MOI burn requires. With the quadratic profile, mid-alignment burns earn very little, so an agent fishing for partial credit gets nothing; a committed periapsis brake earns hundreds.

Mars-Relative Energy

The boolean capture criterion (bound, latched for 3 steps) is a step function on a continuous quantity. Adding a per-step shaping term on the underlying energy turns it into a smooth gradient:

ε = ½|v_rel|² - μ_Mars / r_rel
r_mars_energy = -w_ε × Δε / E_scale
# w_eps = 40, E_scale = 1e6; total reward over hyperbolic→bound brake ≈ +500

A brake reduces ε, so Δε < 0 and the agent gets positive reward proportional to the energy reduction. An earlier saturating-tanh form gave near-zero gradient at the hyperbolic baseline because ε / E_scale ≈ 2.5 saturates the tanh derivative to about 0.07. The linear form gives proportional signal at any ε. Telescoping is preserved (sum equals -w_ε × (ε_final - ε_init) / E_scale), so an oscillating policy that pumps and brakes the same energy gets nothing.

Dense Capture Quality

The boolean criterion only fires once per episode at the moment of capture; until that moment, the Gaussian-kernel proximity score q_moi ∈ [0, 1] (last week's terminal grading) carries no per-step weight. A tiny per-step contribution turns it into a continuous shaping signal:

r_q,dense = w_q × q_moi, w_q = 5
# q_moi = 0 when e ≥ 1, so this only pays once the SC is on a bound orbit

The gating on e < 1 is what makes this safe to add. The agent earns no continuous reward by getting "almost captured"; it only earns once it has actually bent the hyperbola into an ellipse. Across ~50 to 100 in-SOI steps, this sums to +250 to +500 of reward for the post-capture coast.

B-plane Done Right

Last week the B-plane shaping term was off entirely because of two interacting bugs. A saturating tanh potential plus a lazy first-entry initialisation lets it come back online:

φ(b) = -tanh(|b - b_ref| / B_scale), B_scale = 5×10⁷ m
r_bplane = w_b × (φ_curr - φ_prev), w_b = 200
# _bplane_err_prev set lazily on first entry into the 5·R_SOI band

The total contribution is bounded by 2 w_b = 400, so a divergent trajectory cannot accumulate the -50000 penalty that wrecked nasa_v1. The lazy-init removes the spurious one-shot delta the first time the spacecraft enters the evaluation region.

Terminal Rebalancing

Value | Description

+1000 | Capture success bonus (was +500)

+500 | Capture quality scale (× q_moi, was +200)

−100 | Mars impact penalty (was -500)

−100 | Heliocentric escape penalty (was -300)

A committed MOI brake risks both crashing into Mars (too much dv) and slipping back to heliocentric (too little). The old penalty magnitudes dwarfed the energy reward, so any committed brake was net-negative in expectation under exploration noise. The success bonus and quality scale go up by 2.5× in aggregate; the failure penalties drop to roughly the same magnitude as a single bad partial brake. Earth re-impact stays at -500 because that's a Phase 0 anomaly and has nothing to do with the brake-attempt risk profile.

04 — v24: First Capture

With the new reward surface and SOI-edge spawning, training proceeded in two stages. v23 ran with a relatively wide exploration setting (ent_coef = 0.005, log_std_init = -1.0, no clamp) for 4M steps from BC initialisation. Eval reward bounced between 600 and 900 with deterministic captures still at zero, but rollout reward (which includes stochastic action noise) showed individual episodes occasionally crossing the 1500 mark. Some action samples were threading the capture; the deterministic mean at the centre of the action distribution was not.

v24 forces the issue. It resumes from v23's best checkpoint and tightens the action distribution: ent_coef = 0.0005, --reset-log-std -1.5, and a hard log_std clamp of [-2.5, -0.5]. The clamp is the new piece. Without it, gradient updates that nudge log_std outside the box just keep going; with it, every update gets the freedom to move the mean inside a sigma window of 0.08 to 0.61. The constrained policy can no longer hide capture in stochastic outliers because the outliers don't exist anymore.

The first eval after the resume showed a single capture (eval reward 1763 at step 400k). Subsequent evals slipped back to the 700–900 range, but the best_model checkpoint preserved the 400k peak. Re-running deterministic evaluation on that checkpoint gave 5/5 captures across 5 random seeds, all on launch window 0 (2024 OCT 15), with mean reward +1762 and mean closest-approach to Mars under 1 Mkm (versus a 577 Mkm SOI). The model is saved at successful_models/v24_window0_capture.zip.

Where the Reward Came From: Decomposing a capture episode by reward source clarifies which terms actually drove the policy:

# v24 capture episode (window 0, seed 0)
# Per-step alignment (sum over MOI burns): +280
# Per-step energy (Δε ≈ -1.2e7 / 1e6): +480
# Per-step q_moi dense (post-capture coast): +320
# Phase 2 entry bonus: +200
# Terminal capture success: +1000
# Terminal q_moi quality (≈ 0.7): +350
# Fuel + survival + distance: -870
episode reward ≈ +1762

The energy term and the dense quality term together contribute about as much as the terminal success bonus does, which is the right ratio. If terminal reward dominated, exploration would never reach it; if shaping dominated, the policy would farm shaping without capturing.

05 — One Policy, Many Spawn Points

v24 captures only when spawned exactly at the SOI edge. Spawning even 7 days earlier and giving the agent the same Phase 2 budget produces zero captures. Several rounds of sequential per-offset training (v25 through v30) crystallised separate policies for different cruise offsets, but each new offset training run forgot the previous one (v28 was 5/5 at T-30d but 0/5 at T=0; v30 was 5/5 at T-40d but 0/5 at T-30d). A different approach was needed.

cruise_offset_random_choices_s and cruise_offset_random_weights let one PPO run train against multiple spawn points simultaneously. Every reset() picks an offset by weighted random sampling, so a single network sees all four spawn distributions and is required to handle each. v31 used uniform weights ({1, 1, 1, 1}) and plateaued at 10/20 captures across {T = 0, T - 7d, T - 30d, T - 40d}: window 0 and T-30d only. The two failing offsets (T-7d and T-40d) are at adjacent positions in offset space, suggesting the policy got stuck in a local mode that couldn't satisfy them.

v32 keeps the same spawn distribution but oversamples the failing offsets via {1, 3, 1, 3} weights. The extra exposure to T-7d and T-40d perturbs the policy enough to escape the v31 attractor:

Captures | Spawn Point

5/5 | T = 0 (SOI edge)

5/5 | T - 7 d

5/5 | T - 30 d

5/5 | T - 40 d

20 out of 20 deterministic captures from a single network at 1.35M steps, on launch window 0, across all four spawn points. The model lives at successful_models/v32_window0_multioffset_0to40d_5x5_5x5.zip and represents the best Mars policy this project has produced.

06 — The Wall

The plan from there was to extend v32: more cruise offsets, more launch windows, eventually a single policy that handles the full Earth-to-Mars run. Five attempts followed, and all five failed in different ways.

Reverse Curriculum Past T-40d

v33 added T-50d and T-60d to the spawn pool with {1, 1, 1, 1, 2, 2} weights so the new offsets dominated. v34 was a gentler version with only T-45d added at weight 1. Both runs destabilised v32's anchor: eval reward bounced between 700 and -142 with no recovery, and the saved best_model stayed pinned at v32's checkpoint without ever progressing past it. The agent forgets the inner offsets faster than it learns the outer ones.

Multi-Window from v32

v35 and v35b loaded v32 and turned on launch-window randomisation with {1, 2, 2, 2} weights (oversampling the unseen windows 1, 2, 3). The eval reward locked at 125.66, which is exactly the truncate-no-capture-no-crash terminal value. Window 0's capture behaviour was destroyed in the first million steps, while no capture behaviour for windows 1-3 ever appeared. The policy collapsed to a coast.

Window-1 Specialist

v36 isolated window 1 (2026-11-11) with weights {0, 1, 0, 0}. The first eval at 100k steps produced a single capture (reward +2331), suggesting the geometry was reachable. Subsequent evals dropped to negative reward by 1.5M steps. v32's window-0 weights are too crystallised to generalise; the gradient that would push them toward window 1 also pushes them away from window 0, and there is no path that gets close to both.

Multi-Window from BC v5

v37 sidestepped v32 entirely. A new BC checkpoint (bc_policy_v5_p2heavy.pt) was trained on a Phase-2-heavy oracle dataset and used to initialise PPO with all four launch windows enabled. Eval reward started at -907 and degraded to -2589 over 1.15M steps. The BC v5 baseline alone is poor (the dataset is too narrow), and PPO from a weak baseline degrades further rather than improving.

The Pattern: PPO retraining from the v32 anchor breaks v32 faster than it learns anything new. PPO from BC v5 has a baseline reward too low to climb from. Either the policy crystallises and refuses to budge, or it never crystallises in the first place. None of the four extension paths pushed evaluation reward above v32's baseline at any point.

07 — Compose, Don't Retrain

The full Earth-to-Mars problem has been solved twice already this project. The oracle does Phase 0 + 1 + 2 at 97.5% across all four launch windows. v32 does Phase 2 plus the last 40 days of cruise at 100% on window 0. What is missing is one PPO policy that does both. The plan for week 10 sidesteps that requirement: compose the two solvers behind a single predict interface, dispatching by phase.

# MarsHierarchicalController.predict(obs)
if env._phase == PHASE_PRE_TMI: return oracle()
if env._phase == PHASE_CRUISE and t_to_soi > 40d: return oracle()
if env._phase == PHASE_CRUISE: return v32.predict(obs)
# PHASE_MOI
return v32.predict(obs)

The hand-off boundary is exactly T - 40 days, the outermost cruise spawn point v32 was trained against. bc_checkpoints/soi_entry_states.npz was produced by running this same oracle, so the oracle's natural trajectory state at T - 40d matches v32's training distribution by construction. There is no out-of-distribution stretch where the controller has to bridge between two unrelated regimes.

The PBRS state lives entirely on the env (_ref_dev_prev, _bplane_err_prev, _mars_energy_prev) and updates every step regardless of which solver supplied the action. Telescoping is preserved across the hand-off without any controller-side bookkeeping. The interface is the two-line SB3 contract every other validate script uses (predict(obs, deterministic=True) → (action, None)), so the controller drops into main_validate.py and main_montecarlo.py unchanged.

Verification Plan

Both pieces are individually 100% on clean physics, so anything below ~19/20 captures on clean window 0 will be a wiring bug rather than a design flaw. The gating diagnostics: print Δr and Δv at the hand-off step against the cached cruise_r_helio[0, idx_40d], expect sub-km / sub-cm/s under clean physics. Mean Phase 0 dv against the oracle's reported TMI dv (3.5–4.3 km/s for window 0). Mean q_moi against v32's training distribution. Optional perturbed run: 20 seeds, expected 85–95% capture rate, anything below treated as a known limitation rather than a regression.

08 — What Comes Next

Week 10's deliverable is the hierarchical controller running end-to-end on window 0, clean physics, with a full Earth-to-Mars trajectory plot and a Monte Carlo distribution of capture quality. Two failure modes worth watching for: drift in v32's training distribution due to small ephemeris differences between the precomputed cache and the live SPICE state at evaluation time, and any unexpected sensitivity of the PBRS terms to the hand-off timing.

Beyond that, three open paths from this week's failures are worth keeping on the list. A bigger network ([512, 512] instead of [256, 256]) might fit a multi-window policy where [256, 256] over-specialises. An explicit window-id observation component would let a single network condition on which window it's running, replacing implicit feature learning with a direct switch. Per-window v32 specialists plus an ensemble dispatcher would extend the "compose, don't retrain" idea from across-phase to across-window. The hierarchical controller is the fastest route to a working full transfer; whichever of those three the next iteration explores depends on what the validation reveals.

The interesting question for next week is whether the stitched controller is enough to call it solved, and whether anything qualitatively new shows up at the 40-day hand-off boundary that the individual training of either piece would never have surfaced.

Stitching the oracle and v32 into one end-to-end controller, finding the fuel budget that broke the hand-off, and training a perturbation-robust capture specialist that beats the oracle alone. Visit mazzola.dev​​​​‌‍​‍​‍‌‍‌​‍‌‍‍‌‌‍‌‌‍‍‌‌‍‍​‍​‍​‍‍​‍​‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍​‍​‍​​‍​‍‌‍‍​‌​‍‌‍‌‌‌‍‌‍​‍​‍​‍‍​‍​‍‌‍‍​‌‌​‌‌​‌​​‌​​‍‍​‍​‍‌‍​‍‌‍​‌‌​‌‍‍‌‌​​‍‍‌‍​‌‌‍‌​‌‍‌‌‍‍‌‌‍‍​‍‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​‌‍‌​‍​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌‌‌‌​‌‍‌‌​​‌‍‌‌‌​​‍‌‌‍​​‍‌‌​‌‌‌‌‌‌‌‌​‌‌‌​​‌​‍‌‌​‌‌​‌​‍‌​‌​‌‌‍‍‌​‌​​‌‌‌‌‌​‌‌‍‍‍‌‍‌‌‌‌‌‌‌​‌‌​‌​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌​​‍‌‍​‍‌‌‌​‌‍‍‌​​‌‌​‌‌‍​‌‌‌​‌​‍‌‌‌​‌‌​​‌‍‌‍‌​‌​​‍‌‍‍‌‌​​​‍‌‌‌​‌‌​‌​‍‌‌​‍‌​‍‍‌‌​‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​​‍​‌‍‌‍‌‍‍‌‌‍‌‌‌‍​‌‍‌​‌‌​​‌‍​‌‌‌​‌‍‍​​‌‌‍​‍‌‍​‌‍‌‍‌‌‌‌​​‌‍‌​‌‌​‌​​‍‍​‍‌​‍‍‌‍​‌‌​‍‌‌​‌‍‍‌‌‍​‌‍​‌‍‌‌​‍​‍‌‌ for nicer formatting

01 — The Setup

Welcome back. v32 from last week captures Mars 20/20 across four cruise spawn points on window 0, and the oracle solves the full Earth-to-Mars problem at 97.5% across all four launch windows. Both pieces work; what was missing was a way to glue them together. Five attempts to extend v32 by retraining had failed, so the plan for this week was to compose the two solvers behind one predict interface and skip the retraining altogether.

Phase 0 and most of Phase 1 are oracle territory; the last 40 days of cruise plus the MOI brake go to v32. The hand-off boundary is v32's outermost training cruise spawn point, and the waypoint cache that defines that spawn was produced by running this same oracle. Under deterministic clean physics, the oracle's state at T - 40 days matches v32's training distribution by construction.

Source | Description | Focus

utils/controllers/mars_hierarchical.py (hand-off ET fix + fallback gate) | Switch hand-off from a fixed-offset rule to a cached-ET equality so the SC lands exactly on the trained spawn state | Stitching, refined

envs/env_config.py + envs/mars_nasa_env.py (spawn_mass_kg_override) | New config field for realistic depleted-mass spawn so the capture specialist trains under the same fuel budget it gets in the live hand-off | Fuel budget fix

runs/v33_realistic_t40d | PPO resumed from v32 with mass override 252 kg, T-40d spawn pinned, lr=1e-4, ent_coef=0.001; first capture at 200k steps | v33 specialist

utils/controllers/mars_hierarchical.py (obs spoof at hand-off) | Dynamic TIME_TO_GO recomputed from a latched hand-off step so PPO sees ttg≈1 even though env._step≈197 | Distribution match

envs/mars_nasa_env.py (_mars_perturb + 5 channels) | All five Federici-style perturbation channels plumbed through MarsNasaEnv with Mars-realistic defaults gentler than hohmann's | Perturbation framework

runs/v34f_window0_perturbrobust | BC v5 init + multi-offset spawn + 30k km spawn noise + mars perturb defaults; 48/50 = 96% capture under perturbations | Perturb-robust breakthrough

main_hierarchical_validate.py | End-to-end driver with per-episode trajectory plot, Vizard playback, and perturbation toggle | Validation harness

02 — The Hand-Off Misaligned

The first version of the stitched controller dispatched on time-to-SOI: oracle while t_to_SOI > 40 d, PPO once it dropped under. Clean-physics runs showed the moment the controller switched, the spacecraft's heliocentric position differed from the cached T - 40 d waypoint by roughly 2.3 Mkm and the velocity by about 16 m/s. This meant the hand-off was off by exactly one decision step.

The cache is built at decision-step granularity by scripts/precompute_soi_entry.py. The oracle does not stop precisely at T - 40 d; it stops at whichever decision-step boundary happens around that moment, which can be up to a full step on either side. So the cached waypoint labelled "T-40d" actually corresponds to T - 40.32 d, and any fixed-offset rule lands the live spacecraft on the wrong step's worth of cruise.

The fix is to anchor the hand-off on the cached ET itself rather than on a relative time-to-SOI. The npz file records the ephemeris time of every cruise waypoint it stores. The controller now loads cruise_epoch_et[i, idx] alongside the position and velocity, and switches the moment the live ET reaches that exact value. Under deterministic clean physics the result is Δr and Δv at hand-off both round to zero.

# Old (fixed-offset rule):
if t_to_soi(env) > 40*86400: dispatch_to_oracle()
# Live hand-off was off by ~1 decision step → Δr ≈ 2.3 Mkm at switch
# New (cached-ET rule):
if env.et_now() < cruise_epoch_et: dispatch_to_oracle()
# Clean hand-off lands exactly on the training spawn state (Δr ≈ 0)

Fallback Gate

Under perturbations the live state will not match the cached waypoint perfectly even with the corrected dispatch rule. A pos / vel quality gate guards against off-distribution hand-offs: if |Δr| > 150,000 km or |Δv| > 50 m/s, the controller stays on the oracle for the rest of the episode. The gate is a safety net for outliers and is not exercised in v34f's main results, but it gives the system a graceful fallback when the upstream drift is too large for the trained PPO to recover from.

03 — The Fuel Budget

With the hand-off fix in place, the stitched v32 controller still failed at end-to-end. Eval runs all reached Mars's SOI with the right approach geometry, but ~80% of episodes hit fuel depletion before the brake closed and the eccentricity hovered at 0.995, barely hyperbolic. The MOI burn was being initiated at full commitment but cut short. Something about the live state at hand-off was inconsistent with v32's training distribution after all.

Turns out the culprit was the fuel budget. v32 was trained from a SOI-edge spawn with the spacecraft at full tank: m_0 = 1200 kg, all 1100 kg of propellant available for the MOI brake. The policy learned a burn profile that uses roughly 811 kg of that propellant over the 50–80 steps inside Mars's SOI. In the hierarchical run, by contrast, the oracle has already spent fuel on the TMI burn (~888 kg) and a small TCM allowance (~61 kg) during cruise. By the time the spacecraft reaches T - 40 d, total mass is about 251 kg and the available propellant for capture is 171 kg. v32 then over-spreads its trained burn profile and exhausts the tank before periapsis.

A useful aside on Basilisk physics surfaced while debugging this. BSK's spacecraft hub keeps its inertial mass constant, so the dynamics integrator computes Δv = F · t / m_0 regardless of how much propellant has actually been depleted. The propellant bookkeeping (Tsiolkovsky) runs in parallel but doesn't feed back into the integrator. The practical consequence is that a commanded burn magnitude produces the same per-step Δv whether the tracked tank is full or nearly empty; what changes is only the fuel counter and the latched "fuel exhausted" termination flag. The fix therefore lives entirely in the training-time fuel counter that the policy sees, while the simulator dynamics stay unchanged.

04 — v33: Realistic-Mass Specialist

A new cfg.spawn_mass_kg_override field plus a matching --spawn-mass-kg-override CLI flag let the SOI-edge spawn write any chosen total mass into the env at reset, while leaving the BSK hub mass alone (so per-step Δv is unaffected). v33 is a fine-tune of v32 with the override set to 252 kg, the launch window pinned to 0, and the cruise offset fixed at T - 40 d to match the hand-off boundary:

# v33 fuel-budget specialist (resumes from v32 weights)
python main_train.py --env mars_nasa \
--resume successful_models/v32_window0_multioffset_0to40d_5x5_5x5.zip \
--cruise-offset-seconds 3456000 --spawn-mass-kg-override 252 \
--epoch-weights "1,0,0,0" \
--lr 1e-4 --ent-coef 0.001 --reset-log-std -1.5 \
--clamp-log-std --clamp-log-std-min -2.5 --clamp-log-std-max -0.3

First capture event landed at the 200k-step evaluation: +1521 reward, ep_len 186, deterministic. Training plateaued from there but the best checkpoint was preserved. The recipe is the same v24 / v32 recipe (low entropy, clamped log_std) with the fuel budget as the only new variable.

Obs Spoof at the Hand-Off

v33 is trained from env._step = 0 but at the live hand-off env._step is around 197 (Phase 0 plus 196 cruise steps). The observation vector contains a TIME_TO_GO component that v33 saw as roughly 1.0 during training and would see as ~0.55 in deployment if fed raw. To keep the live distribution close to training, the controller latches the step index at the first PPO dispatch and recomputes the observation's TIME_TO_GO dynamically:

TIME_TO_GO_spoofed = max(0, 1 - (env._step - step_at_handoff) / n_decision_steps)
# For v33 (trained at mass_kg=252): MASS_FRAC is already correct, leave it.
# For v32 (trained at full tank): also spoof MASS_FRAC = 1.0.

A --ppo-realistic-mass / --no-ppo-realistic-mass toggle picks the right spoofing behaviour for whichever PPO is being run, since v32 and v33 have different training-mass assumptions and the controller has to match each one.

05 — End-to-End

The stitched controller with the v33 PPO, the cached-ET hand-off rule, and the obs spoof produces a complete Earth-to-Mars run from LEO to a bound Mars orbit. On clean physics, window 0, fixed parking orbit, 5 random seeds: 5/5 captures, mean reward +1565, mean capture eccentricity 0.940, mean periapsis altitude 84,614 km. The original project goal of an autonomous Earth-to-Mars transfer is met.

Five-episode evaluation of the hierarchical controller on window 0, clean physics. Heliocentric plane (left) with Earth's and Mars's full-year orbits as dashed circles, the spacecraft traces (coloured by episode), and burn events marked with cyan dots (oracle) or red stars (PPO). Distance to Mars (top right) decays from 170 Mkm at TMI to inside the SOI by step 200. Burn Δv per step (bottom right, symlog) shows the TMI as a ~3.5 km/s spike at step 0, small oracle TCMs through cruise, and a concentrated PPO MOI burn around step 400.

The trajectory plot makes the role of each solver visible. The blue Phase-0 spike is the oracle's Lambert TMI, sized at one decision step. The cyan dots scattered through cruise are oracle TCMs, each a few m/s of trajectory correction. The red star cluster near step 400 is v33's MOI burn, concentrated into roughly six high-magnitude steps right at periapsis. The cyan-to-red transition is the hand-off boundary, T - 40 d before SOI entry.

What's Doing What: Both solvers contribute to the final +1565 reward, but not symmetrically:

# Hierarchical capture (window 0, seed 0)
# Oracle Phase 0 TMI: ~3.5 km/s, ~888 kg fuel
# Oracle Phase 1 cruise TCMs: ~40 m/s, ~61 kg fuel
# PPO Phase 1 endgame (40d): ~3 m/s, ~5 kg fuel
# PPO Phase 2 MOI brake: ~430 m/s, ~165 kg fuel
capture: e = 0.940, r_p_alt = 84,614 km

06 — Federici-Style Perturbations

The clean-physics 100% rate is the easy half of robustness. The harder question is what survives when the simulator stops being a clean ephemeris integrator. Roberto Federici's Deep Learning-Based Spacecraft Optimal Guidance dissertation defines a five-channel noise model for exactly this scenario with navigation sensor noise, thruster magnitude noise, thruster pointing noise, unmodelled-force impulses, and stochastic missed-thrust events. The Hohmann env already had these wired in while the Mars NASA env had been ignoring its cfg.perturb entirely.

All five channels are now plumbed through MarsNasaEnv:

Value | Channel

1 km | Position obs noise, 1σ per axis

0.5 m/s | Velocity obs noise, 1σ per axis

0.5% | Thruster magnitude, 1σ multiplicative

0.3° | Thruster pointing, 1σ axis tilt

0.05 m/s | Dynamics impulse, 1σ per step

0.5% | Missed-thrust prob per step (capped at 1/ep)

The Mars defaults are roughly 10× milder than the Hohmann defaults. A first attempt at "v34_robust" used the Hohmann-scale perturbations on Mars and the result was a 100% failure rate: a 5% thruster magnitude error on the 3.5 km/s TMI burn translates to about 175 m/s of Δv uncertainty, which by the time the spacecraft reaches Mars compounds to roughly 1.2 Mkm of SOI miss. The oracle then cannot recover the trajectory inside its Phase 1 Δv budget. Gentler defaults at the Mars scale let the oracle do its job over the 240-day cruise.

Sweeping the oracle alone across {0×, 1×, 5×, 10×} the Mars defaults gives capture rates of 100%, 90%, 80%, 70% across 20 seeds each. At the realistic 1× scale, the oracle drops to 90%. At Federici's full strength (10×) it drops to 70%. The Mars defaults are the regime our launch windows actually face.

07 — v34f: Perturbation-Robust Capture

With perturbations on, v33 alone in the hierarchy collapses from 100% to 25% capture. The drift accumulated through 240 days of cruise puts the spacecraft well off the cached spawn waypoint at hand-off, and v33's training distribution is too narrow to recover. The obvious move is to widen that training distribution, which means training a new capture specialist under the same noise the deployment will see.

Three earlier attempts to do this failed. v34c fine-tuned v33's anchor with 50 km/axis spawn-state noise and collapsed to a coast policy. v34d turned perturbations on without spawn noise and collapsed to the same coast policy. v34e used a narrower 5 km/axis spawn noise on v32 instead of v33, and the basin barely widened (10% capture, deployment-shifted by an average of 230k km from training spawn). This shows that fine-tuning from v32 or v33 anchors with anything resembling realistic spawn-distribution width destroys the policy faster than the new gradient can rebuild it.

v34f changes the initialisation. Instead of resuming from a previous PPO checkpoint, it warm-starts from a behaviour-cloning baseline trained on a Phase-2-heavy oracle dataset (bc_policy_v5_p2heavy.pt). The BC loss landscape is flatter than v32's or v33's sharp local optimum, so PPO has room to find a broader basin without first having to escape an existing one:

# v34f training recipe
python main_train.py --env mars_nasa \
--bc-init bc_checkpoints/bc_policy_v5_p2heavy.pt \
--cruise-offset-random-choices-s "0,604800,2592000,3456000" \
--cruise-offset-random-weights "1,3,1,3" \
--spawn-pos-noise-km 30000 --spawn-vel-noise-ms 20 \
--spawn-mass-kg-override 252 --epoch-weights "1,0,0,0" \
--lr 1e-4 --ent-coef 0.002 --log-std-init -2.0 \
--clamp-log-std --clamp-log-std-min -2.5 --clamp-log-std-max -0.3 \
--timesteps 1_000_000 --n-envs 8

The 30,000 km / 20 m/s spawn-state noise covers the actual deployment drift distribution (mean ~230k km cumulative miss, with the long tail captured by the noise's 3σ wing). 5 km was too narrow, 50 km from a sharper anchor was too wide, 30 km from a BC baseline lands in the sweet spot.

Save Best Model Early

Evaluation reward peaks around step 200,000 (mean 985, std 524 across episodes, with the large variance reflecting captured-vs-not-captured mixing) and then degrades back toward a coast policy by 1M steps (eval reward 528). The standard SB3 best-checkpoint callback caught the peak and the late degradation didn't overwrite it.

Result

Over 50 perturbed episodes on window 0:

Captures | Configuration

48 / 50 | v34f hier (perturbed, this work)

19 / 20 | Oracle alone (perturbed, 95%)

5 / 20 | v33 hier (perturbed, 25%)

5 / 5 | v33 hier (clean physics)

v34f under perturbations not only matches the clean-physics rate of the oracle alone, it produces a markedly tighter capture orbit: mean eccentricity 0.59 versus the oracle's 0.86, mean periapsis altitude 46,000 km versus the oracle's 102,000 km. The PPO uses the extra information in its observation (Mars-relative state, B-plane error) to commit a sharper MOI brake than the oracle's open-loop Lambert solution.

Why BC Initialisation Worked: v32 and v33 are sharp local optima of the Phase-2 reward surface, fitted to narrow spawn distributions. A gradient step in any direction that helps generalisation also moves the policy away from one of those optima, and the reward signal under exploration noise is dominated by the immediate loss of capture probability. PPO retreats. BC v5 sits in a flatter region of the same loss landscape with smaller second derivatives, so PPO can shift the policy along the generalisation direction without paying as much in the local capture rate.

The model is window-0 only. The same training recipe with --epoch-weights set differently is the natural starting point for the remaining three launch windows, but the proof-of-concept for the project's headline goal is complete with a learned policy delivers an Earth-to-Mars transfer under realistic deployment noise.

Eleven weeks of work compressed into one pass. Visit mazzola.dev​​​​‌‍​‍​‍‌‍‌​‍‌‍‍‌‌‍‌‌‍‍‌‌‍‍​‍​‍​‍‍​‍​‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍​‍​‍​​‍​‍‌‍‍​‌​‍‌‍‌‌‌‍‌‍​‍​‍​‍‍​‍​‍‌‍‍​‌‌​‌‌​‌​​‌​​‍‍​‍​‍‌‍​‍‌‍​‌‌​‌‍‍‌‌​​‍‍‌‍​‌‌‍‌​‌‍‌‌‍‍‌‌‍‍​‍‍‌‍‌​‌‍​‌‌‌​‌‍‌‍​‌‍‌‌​​‍‍‌‍​‌‍‌‍‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​‌‍‌​‍​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌‌‌‌​‌‍‌‌​​‌‍‌‌‌​​‍‌‌‍​​‍‌‌​‌‌‌‌‌‌‌‌​‌‌‌​​‌​‍‌‌​‌‌​‌​‍‌​‌​‌‌‍‍‌​‌​​‌‌‌‌‌​‌‌‍‍‍‌‍‌‌‌‌‌‌‌​‌‌​‌​‍‌‍‍‌‌‌​‌‍‌‌‌‍‌‌​​‍‌‍​‍‌‌‌​‌‍‍‌​​‌‌​‌‌‍​‌‌‌​‌​‍‌‌‌​‌‌​​‌‍‌‍‌​‌​​‍‌‍‍‌‌​​​‍‌‌‌​‌‌​‌​‍‌‌​‍‌​‍‍‌‌​‌​‍‌‍‌‌‌‍‌​‌‍‍‌‌‌​​‍​‌‍‌‍‌‍‍‌‌‍‌‌‌‍​‌‍‌​‌‌​​‌‍​‌‌‌​‌‍‍​​‌‌‍​‍‌‍​‌‍‌‍‌‌‌‌​​‌‍‌​‌‌​‌​​‍‍​​‌​​​​‍‍‌‍​‌‌​‍‌‌​‌‍‍‌‌‍​‌‍​‌‍‌‌​‍​‍‌‌ for nicer formatting

01 — Where We Started

Eleven weeks ago this blog opened with the cost of mass to orbit ($7K per kg to LEO, $1.2M per kg to the Moon's surface) and a question of whether a neural network could learn to fly a deep-space transfer well enough to matter. The setup was PPO against a high-fidelity Basilisk simulator with the LEO-to-MEO Hohmann transfer as the warm-up, and the standing goal was Earth-to-Mars under realistic deployment noise.

The Mars problem asks for three things in series. A 3.5 km/s departure burn aligned with Earth's orbital velocity has to put the spacecraft on a Lambert arc, a 240-day cruise has to hold that arc through dynamics-impulse noise and the occasional missed thrust, and a final ~430 m/s brake at Mars periapsis has to land inside a Gaussian-kernel capture region around a target B-plane.

02 — What Got Built

The final system is four layers. At the bottom sits a 6-DOF Basilisk core with SPICE DE430 ephemerides, n-body Sun-Earth-Mars dynamics, solar radiation pressure, and gravitational harmonics. Above that lives a custom Gymnasium environment with a three-phase finite-state machine (Pre-TMI / Cruise / MOI), a 4-D inertial-thrust action space, and a B-plane targeting reward. The PPO capture specialist is a [256, 256] MLP trained with Stable-Baselines3, warm-started from a behaviour-cloning baseline on an oracle dataset and refined through a reverse curriculum against cached hand-off states. A hierarchical orchestrator sits on top, dispatching an energy-tracking Lambert oracle for Phase 0 and most of Phase 1, then handing the spacecraft off to the PPO inside Mars's sphere of influence.

Count | Description

3 | Custom Gymnasium envs (Hohmann, Mars 2D, Mars NASA)

5 | Perturbation channels plumbed end-to-end

4 | Production launch windows (2024, 2026, 2028, 2033)

170+ | pytest tests across environments, oracles, and perturbations

The decomposition into oracle plus learned controller carries the project's central design call. Phase 0 is well-served by classical methods (Lambert solvers have been in the trajectory-optimisation literature since the 1960s). Phase 2 is where the learned policy actually contributes, because the brake has to commit under uncertain incoming state, finite fuel, and a non-trivial reward landscape that punishes both under-braking (escape) and over-braking (Mars impact). The orchestrator routes each phase to whichever solver is better suited.

03 — The Path There

Most of the project's training compute was spent inside a flat plateau. nasa_v1 collapsed to evaluation reward near -50,000 across all 10M of its training steps, traced to three independent reward-shaping bugs in the env (unbounded reference PBRS, an init-time delta on the B-plane term, and an Earth-impact penalty that triggered on natural LEO perigee drift). nasa_v2 fixed those, and then plateaued at +75 reward for another 10M steps with zero captures. A 4-D log_std runaway turned the deterministic mean into noise; a curriculum threshold pinned at 60% capture meant the agent never advanced past Stage 1. v3 through v22 cycled through alignment bonuses, distance shaping, fuel weights, and BC dataset composition with all the obvious variants tried. None of them captured.

Schematic of the project's training arc across four representative runs. Roughly 60M cumulative PPO steps of plateau preceded the first capture at v24. Each subsequent breakthrough was a recipe change rather than a scale-up.

The first capture (v24) came after three structural changes that arrived together. scripts/precompute_soi_entry.py cached the oracle's heliocentric state at Mars SOI edge so the env could spawn the spacecraft mid-flight; a Phase-2 reward overhaul added per-step alignment, Mars-relative energy, and a dense capture-quality term gated on e < 1; and a hard log_std clamp of [-2.5, -0.5] forced the deterministic mean to commit a brake instead of hiding capture probability in stochastic outliers. v24 produced 5/5 deterministic captures on launch window 0 at the SOI-edge spawn.

v32 generalised that to 20/20 captures across four cruise spawn points {T = 0, T - 7d, T - 30d, T - 40d} from a single network, by oversampling the failing offsets with {1, 3, 1, 3} weights to break a multi-modal plateau. Five attempts to extend v32 further (reverse curriculum past T - 40d, multi-window from v32, window-1 specialist, multi-window from a fresh BC) failed for the same reason in different shapes: the v32 anchor was too crystallised to absorb new spawn distributions without losing the original ones. v34f sidestepped the anchor entirely by initialising from a behaviour-cloning baseline (bc_policy_v5_p2heavy.pt) with a flatter loss landscape, and reached 96% capture under realistic perturbations.

04 — The Headline

The hierarchical controller running v34f under Mars-realistic perturbations is the project's terminal result. Across 50 perturbed Monte Carlo runs on launch window 0, the system delivers a bound Mars orbit on 48 of them. The two failures both correspond to large upstream cruise drift that the fallback state-quality gate intercepted, sending the episode back to the oracle and producing a partial (heliocentric, no impact) result instead of a capture.

Metric | Description

100% | Capture under clean physics (5/5, v33 hier)

96% | Capture under perturbations (48/50, v34f hier)

e = 0.59 | Mean capture eccentricity (oracle alone: 0.86)

46k km | Mean periapsis altitude (oracle alone: 102k km)

Capture quality across 50 perturbed Monte Carlo runs per controller. Periapsis altitude on the horizontal axis, capture eccentricity on the vertical. The TOAST cluster ends up lower-left of the oracle, with tighter and more circular orbits.

The TOAST cluster also beats a Lambert-only baseline on the orbital quality of its captures. The mean eccentricity drops from 0.86 to 0.59 and the mean periapsis altitude drops from 102,000 km to 46,000 km. A bound orbit at e = 0.6 / r_p = 46,000 km asks for considerably less circularisation dv than e = 0.86 / r_p = 102,000 km, so the learned policy's contribution shows up as both more frequent captures and lower fuel cost for the orbit it leaves the spacecraft in. The PPO uses the extra information in its observation (Mars-relative state, B-plane error, time-to-go, mass fraction) to commit a sharper brake at periapsis than the oracle's open-loop Lambert solution.

05 — Where It Sits

Roberto Federici's 2020 dissertation Deep Learning-Based Spacecraft Optimal Guidance is the closest prior art for this style of work. The Federici thesis trains a small policy network against simplified two-body dynamics with deliberately strong noise channels (5% thrust magnitude, 1° pointing, ~1 m/s impulse per step, 2% missed-thrust probability), and demonstrates a robust brake under those conditions. The five-channel perturbation framework in this project is a direct port of the Federici noise model, and the v34f recipe inherits the spirit of training against the noise the deployment will face.

Schematic regime map for the three approaches. Classical open-loop trajectory optimisation lives at low perturbation and easier dynamics; Federici 2020 explores extreme noise but on a toy 2-body simulator; this project pushes the problem into flight-grade Basilisk physics under mission-realistic perturbation magnitudes.

The substantive difference is the physics layer underneath. TOAST runs on 6-DOF Basilisk with full SPICE DE430 ephemerides, n-body perturbations from the Sun, gravitational harmonics, and solar radiation pressure, while the Federici simulator is a two-body integrator with parameterised noise. The Federici noise levels are also two orders of magnitude stronger than what a modern flight system actually sees in cruise; a 5% magnitude error on the 3.5 km/s TMI burn translates to about 175 m/s of Δv uncertainty, which compounds to 1.2 Mkm of SOI miss by the time the spacecraft reaches Mars, and no recoverable trajectory exists inside a typical Phase-1 Δv budget. Scaling the noise back to mission-real magnitudes (0.5% magnitude, 0.3° pointing, 0.05 m/s impulses, 0.5% missed-thrust) puts the project inside the regime that an actual planning team would care about.

What This Adds: Two things this project contributes on top of prior art: a flight-grade physics stack for an RL-driven Mars capture (Basilisk + SPICE + n-body + harmonics + SRP, all running in the training loop), and a hierarchical decomposition that lets a learned policy improve on a Lambert baseline in capture-orbit quality without having to relearn the well-understood Phase 0 from scratch. Neither piece is novel on its own; their combination on a single end-to-end trained controller is, as far as I can tell, new.

06 — Carry-Outs

The single-window v34f result is a proof of concept rather than a deployed-grade policy. Three obvious extensions follow from the infrastructure that's already in place:

# Path 1: All four launch windows
python scripts/precompute_soi_entry.py --window 1 # ditto for 2, 3
python main_train.py --env mars_nasa --bc-init bc_checkpoints/bc_policy_v5_p2heavy.pt \
--epoch-weights "1,1,1,1" --cruise-offset-random-choices-s "..." \
--spawn-pos-noise-km 30000 --spawn-vel-noise-ms 20 --spawn-mass-kg-override 252 \
--lr 1e-4 --ent-coef 0.002 --log-std-init -2.0 --timesteps 2_000_000
# Path 2: Variable parking orbit
# Regenerate the SOI-edge cache per parking-radius bin, train an ensemble,
# or add parking-radius as an explicit observation component.
# Path 3: Bigger network + window-id conditioning
# Replace [256, 256] with [512, 512] and feed window-id as a one-hot.
# Lets one network condition on which window it's running.

Path 1 is the natural follow-on for next semester. The recipe that worked on window 0 should transfer with minimal recipe changes once the cache is regenerated per window; the failure modes from week 9's multi-window attempts came from fine-tuning the v32 anchor rather than from the BC-from-scratch recipe v34f uses. Path 2 trades flexibility for cache cost and is probably needed before any operational deployment. Path 3 is the most speculative; the empirical evidence so far is that the [256, 256] MLP fits a single-window manifold cleanly, but whether it generalises across four launch windows is open.

The deeper open question is whether the hierarchical decomposition itself is necessary. End-to-end PPO from a fresh BC baseline never matched v24's anchor across 60M cumulative steps; the cause was reward-surface curvature rather than network capacity. A reward function that explicitly credits long-cruise stationarity (instead of telescoping all PBRS terms to zero across coasting steps) might let a single policy do Phase 0 + Phase 1 + Phase 2 without the Lambert oracle. That route would replace the orchestrator with a single network and would let the agent discover non-Lambert TMI strategies (low-thrust spirals, multi-revolution transfers, gravity-assist routes), but the reward design alone is a research project's worth of work.

Beyond Mars, the same stack should extend straight to any high-thrust capture problem: lunar orbit insertion, near-Earth asteroid rendezvous, Jovian moon capture. The only Mars-specific code is the launch-window catalogue and the SPICE body ID; the rest of the pipeline is parameterised on body mass, SOI radius, and reference epoch.

07 — Sign-Off

A learned neural-network policy can fly a spacecraft from low Earth orbit to a bound Mars orbit through 240 days of perturbed cruise, 96% of the time, on a single launch window. The orbital state at capture is tighter than what the classical Lambert-based open-loop baseline produces on the same problem. The training recipe runs on a desktop with a single GPU in under a week of wall-clock time. That answers the question this blog opened with.

Most of those eleven weeks were spent on failures. nasa_v1 through nasa_v22 plateaued or collapsed across roughly 220M cumulative training steps. v33 through v37 produced four different ways to break v32. v34c through v34e produced three more ways to break v33. The single breakthroughs (v24, v32, v34f) each took about one week of focused recipe iteration on top of accumulated debugging. None of them were predicted by classical RL theory; they came from staring at training curves, decomposing reward into per-source contributions, writing diagnostic scripts that compared deterministic mean against rollout sample, and refusing to trust any single number.

Thanks for reading along. End of transmission.