The variable that makes Lambert's problem easy

Lambert’s theorem says that the time to fly between two points on a conic depends only on the semi-major axis, the chord between the points and the sum of their distances from the centre. So the transfer between two positions in a given time is the root of one equation in one unknown. That makes Lambert’s problem sound finished: write the equation, find the root.

It is finished as physics. As computation it took two centuries to become easy, and it became easy not through a better root-finder but through a better unknown. The same equation, written in the semi-major axis, has a feature that makes every root-finder struggle near the most important transfer of all — the minimum-energy one. Written in a variable introduced by Lancaster and Blanchard in 1969, it has no such feature, and a plain Newton iteration converges in a handful of steps from the same starting point for every time of flight. The difficulty was in the label on the orbits, not in the orbits.

Two branches and one join #

For a given chord and pair of radii there is a smallest semi-major axis an ellipse joining the two points can have — the minimum-energy ellipse, whose semi-major axis $a_m$ is half the semi-perimeter of the triangle formed by the two points and the centre. Every larger semi-major axis gives two ellipses through the two points, because an ellipse’s second focus can lie on either side of the chord. One of them gets from the first point to the second quickly, the other slowly. At $a_m$ the two coincide.

Lagrange’s form of the time of flight makes the two branches explicit. With angles defined by $\sin(\alpha/2) = \sqrt{s/2a}$ and $\sin(\beta/2) = \sqrt{(s-c)/2a}$,

$$\sqrt{\mu}\,t = a^{3/2}\left[(\alpha - \sin\alpha) - (\beta - \sin\beta)\right]$$

on the fast branch, with $\alpha$ replaced by $2\pi - \alpha$ on the slow one. As $a$ falls towards $a_m$, $\sin(\alpha/2)$ rises to one, $\alpha$ reaches $\pi$, and the two expressions meet. Near that point $\alpha$ varies as the square root of $a - a_m$, and so does the time of flight: the curve of time against semi-major axis has a vertical tangent at the join, the shape of $\sqrt{a - a_m}$ turned on its side.

The minimum-energy transfer is not an obscure special case. It is the transfer that needs the least energy of all those joining the two points, and for most interplanetary launch windows the transfers actually flown have semi-major axes close to it. A solver that struggles near $a_m$ struggles near the answers that matter most.

The ellipse at the join #

The minimum-energy ellipse has a geometric description that explains why it sits at the join. For any ellipse through the two points, the sum of the distances from each point to the two foci is $2a$. The occupied focus is the centre of attraction, fixed; the empty focus can be anywhere that satisfies both conditions, which places it at an intersection of two circles, one round each point. When $a$ is large the circles cross at two places, one on each side of the chord, and those are the two ellipses — the fast and the slow. As $a$ shrinks the circles cross closer to the chord, and at $a_m$ they just touch, on the chord itself. The minimum-energy ellipse is the one whose empty focus lies on the line between the two points, and below $a_m$ the circles do not meet at all.

That is the vertical tangent seen geometrically. Near $a_m$ a small increase in $a$ moves the two intersection points apart along a direction perpendicular to the chord, and the distance they move goes as the square root of the increase — the same square root that appears in $\alpha$. The time of flight depends on where the empty focus is, so it inherits the square root.

The same ellipse is the familiar one in a special case. When the two points are on opposite sides of the centre — a transfer through 180 degrees between circular orbits — the chord passes through the centre, the empty focus on it lies on the line of apsides, and the minimum-energy ellipse is the Hohmann transfer, tangent to both orbits. Lambert’s minimum-energy ellipse is the Hohmann transfer’s generalisation to any two points, which is why the cheapest windows cluster near it.

What a vertical tangent does to Newton’s method #

Newton’s method finds the root of $f(u) = 0$ by stepping $u \to u - f(u)/f'(u)$ — following the tangent line down to zero. It converges quadratically, doubling the number of correct digits each step, when it starts close enough and when $f'$ is well-behaved. A vertical tangent breaks the second condition in the worst way.

Near the join the slope in $a$ is huge, so a Newton step divides the error in the time by a huge number and moves $a$ hardly at all when it is on the steep side; from further away the tangent line, extended, crosses zero far beyond the join, into values of $a$ smaller than $a_m$ where no ellipse exists and the formula has no real value at all. A solver has to detect that, retreat, and try again, and the steps it takes near the join are small and slow. It also has to know in advance which branch it is on, which means knowing whether the requested time is above or below the minimum-energy time — which in turn needs the minimum-energy time computed first.

A well-determined root and a badly conditioned method #

The vertical tangent has an effect that sounds paradoxical until it is separated into two questions. How precisely does a measured time of flight determine the semi-major axis? And how easily can a root-finder find it?

The answer to the first is: extremely precisely near the join. Near $a_m$ the time changes infinitely fast with $a$, so the inverse function — $a$ as a function of the time — has zero slope there: an error in the time produces almost no error in the semi-major axis. The root is exceptionally well conditioned. What is badly conditioned is the other thing, the branch: a time a little above the minimum-energy time and a time a little below it correspond to nearly the same $a$ but to different branches, and at the join itself a tiny error in the time decides which of two geometrically different ellipses is the answer.

A root-finder working in $a$ has to deal with both. It is searching a variable in which the function is steep and the branch is ambiguous, and its trouble near the join is the trouble of approaching a square-root branch point from one side. A root-finder working in $x$ sees neither: the function is smooth, the branch is the sign of $x$, and the join is an ordinary point at $x = 0$ where nothing special happens. Conditioning is a property of a problem together with the variables it is posed in, and changing the variables can change it without changing a single solution.

The substitution that removes it #

Lancaster and Blanchard introduced a variable, now usually written $x$, defined by $a = a_m/(1 - x^2)$, with the sign of $x$ distinguishing the branches: $x$ between 0 and 1 is the fast branch, $x$ between −1 and 0 the slow one. The minimum-energy ellipse is $x = 0$, the parabola is $x = 1$, hyperbolas have $x \gt 1$, and $x \to -1$ is an ellipse of unbounded size taking unbounded time.

The substitution does two things at once. It joins the two branches into one: a single variable runs continuously from the slowest transfer to the fastest, and the question of which branch the answer is on disappears into the sign of the answer. And it removes the vertical tangent, because near $x = 0$ the quantity $a - a_m$ goes as $x^2$, and the square root that made the time of flight vertical in $a$ becomes a linear function of $x$. The singularity at the minimum-energy ellipse was the square root in the parameterisation, and the new parameter is the square root taken in advance.

That is the same move as the change of variable that regularises a collision orbit: a difficulty that appeared to be in the dynamics turned out to be in the coordinates, and the cure is not a better method but a better description. It is also the idea behind universal variables, which remove the distinction between ellipses and hyperbolas from Kepler’s equation in the same way that $x$ removes the distinction between the two branches here — and in fact $x$ runs smoothly through the parabola into hyperbolic transfers too, so one solver handles every conic.

Counting the steps #

The difference shows most directly in the number of Newton steps each form needs to converge.

The averages are not very different, and that is instructive in itself. Away from the minimum-energy time the curve in $a$ is well-behaved enough, and a Newton iteration in $a$ is only a little worse than one in $x$. The whole cost is concentrated near the join, where the semi-major axis form needs three times as many steps — and that is exactly where the cheapest transfers are.

The dashed curve shows the second half of the modern treatment. In $x$ the remaining difficulty is at the slow end, where the time of flight runs to infinity as $x$ approaches −1 and the curve steepens. Taking the logarithm of the time straightens that end, and Dario Izzo’s 2015 formulation, used in the European Space Agency’s trajectory tools, works in exactly those terms: $x$, the logarithm of the time, an initial guess built from the curve’s known values at $x = 0$ and $x = 1$, and a third-order Householder step instead of Newton’s second-order one. The count falls to two or three iterations for almost every case, and the solver needs no branch logic at all.

The long way round #

A transfer angle greater than 180 degrees — going the long way round the central body — changes the sign of $\lambda = \sqrt{1 - c/s}$, the one geometric parameter the time-of-flight curve depends on, and gives the curve a different shape.

That last property is the deepest simplification the variable brings. In the semi-major axis, each geometry has its own curve and its own join. In $x$ and in the non-dimensional time, every single-revolution Lambert problem in the universe lies on a curve labelled by one number, and a solver can be designed once, for all of them, with an initial guess that is good everywhere. The multiple-revolution problem, where one time of flight admits several ways round, adds a whole number of revolutions to the angle $\psi$ in the formula and gives a family of curves each with its own minimum, and the same variable handles it with the same uniformity.

Gauss’s version of the same trick #

The problem of choosing the right unknown is older than spaceflight. In 1801 Carl Friedrich Gauss had to find the orbit of the newly discovered and then lost asteroid Ceres from a few weeks of observations, and at the heart of his method was a Lambert-type problem: given two positions and the time between them, find the orbit. He did not solve for the semi-major axis. He solved for the ratio of the area of the sector swept out between the two positions to the area of the triangle formed by them and the Sun — a ratio close to one for short arcs, varying smoothly, and expressible as a rapidly convergent series. The iteration on that ratio converged quickly for the short arcs Gauss had, and Ceres was found again where he predicted, from a handful of angles and no distance.

The sector-to-triangle ratio is badly behaved for long arcs, which Gauss did not need; the semi-major axis is badly behaved near the minimum-energy ellipse, which interplanetary transfers do need; Lancaster and Blanchard’s $x$ is well behaved everywhere. Each generation of the problem chose the variable its own cases needed, and the history of Lambert’s problem is largely a history of finding a variable in which the cases of current interest were easy — a better unknown rather than a better equation.

Why speed of solution matters #

A single Lambert solve takes microseconds either way, and the difference between five steps and fifteen seems unimportant. It becomes important at the scale at which trajectories are actually designed.

A porkchop plot is a map of the energy cost of a transfer over a grid of departure and arrival dates; a fine one for a single launch window is a million Lambert solves. The search for a gravity-assist route through several planets evaluates a Lambert arc for every leg of every candidate sequence of encounter dates, and a global optimiser may evaluate billions. A mission-design tool that runs through the whole of a decade’s opportunities for a sample return, or for a tour of the outer planets, is dominated by the Lambert solver’s inner loop. At that scale, a factor of three in iterations is a factor of three in what can be searched, and a solver that occasionally fails to converge — as solvers in the semi-major axis sometimes do near the join — leaves holes in the map in exactly the regions of lowest cost.

Robustness matters more than speed. A solver that needs branch logic, a separate computation of the minimum-energy time and special handling near the join has more places to fail, and a failure in a global search is silent: the optimiser simply never sees the trajectory it could not compute. The variable $x$ removes those places, and that — more than the reduction in average iterations — is why it replaced the older forms.

Past the parabola #

The curves in the figures stop at $x = 1$, the parabola, but the variable does not. For $x \gt 1$ the same formula, with the square roots and the arc-cosine replaced by their hyperbolic counterparts, gives the time of flight on hyperbolic arcs, orbits with no period that pass between the two points faster than any parabola can. The curve continues smoothly through $x = 1$, falling towards zero time as $x$ grows, with no join and no change of branch.

That matters for two kinds of problem. A spacecraft flying past a planet follows a hyperbola relative to it, and targeting a flyby — arriving at a chosen point near the planet at a chosen time — is a Lambert problem on a hyperbolic arc. And very fast transfers between planets, of the kind a high-energy propulsion system might make, are hyperbolic relative to the Sun. A solver in $x$ handles them with the same code and the same starting guess; a solver in the semi-major axis needs a third branch with a negative semi-major axis and its own formulas, joined to the elliptic ones at the parabola where the semi-major axis is infinite. The parabola, which is a singular point of the semi-major axis, is the ordinary point $x = 1$ of the better variable — the third singularity the substitution removes, after the join and the branch.

What was solved and how #

Every solve here is plain Newton with a numerical derivative, starting from a fixed point — $x = 0$, or twice the minimum-energy semi-major axis. Production solvers use analytic derivatives, higher-order steps and initial guesses fitted to the curve, and converge in two or three steps in both forms away from the join; the comparison drawn is of the variables, not of the best achievable solver in each.

The transfer is single-revolution and its time is non-dimensionalised by the semi-perimeter, so the curves describe any central body and any scale. The Earth–Mars-like radii are illustrative; the pattern holds for any geometry except the degenerate one at exactly 180 degrees, where the plane of the transfer is undefined and a different problem arises.

And the semi-major axis solver is told which branch to use. A solver that had to decide the branch itself would need the minimum-energy time computed first, and near the join an error in that decision sends it to the wrong branch — a failure the figures do not show, and the most common one in practice.

Still open: the last singularity #

The variable $x$ removes the minimum-energy singularity and joins the branches, and the logarithm of the time straightens the slow end. One singular case remains in every formulation: the transfer through exactly 180 degrees, where the two position vectors and the centre are collinear and the plane of the transfer is not determined by them. Near it the solutions are well-conditioned in the time of flight but the orbit’s orientation is not: a tiny change in the positions swings the transfer plane through large angles, and the inclination change the transfer implies becomes a wall across the porkchop plot. No change of variable removes that, because it is not in the parameterisation of a family of conics in a plane but in the choice of the plane itself — a genuine property of the problem rather than of its description, and the one place where Lambert’s problem is still hard for a reason that is not a choice of coordinates.