LAMBERT'S PROBLEM: MATHEMATICAL INSANITY

Special thanks to /u/jofwu for helping me understand these concepts.

Lambert's problem is the solution to the issue of non-circular orbit transfers. What it does is it takes two time points (A start and an end) and calculates the ΔV for that point. This is why I call it mathematical insanity: You do the same thing over and over again and hope for a better result each time.

Starting off

So, to start, you need to know where the target and interceptor are at the start, and where the interceptor is at the finish. So, we'll need these values for both orbits:

Along with a couple of constants derived from them:

insert image of 2π*sqrt(a^3/k)=p insert image of sqrt(k/a^3)=N

N will be useful later, because it's the factor for each orbit for doing what we have to do next: Stepping time forward and seeing into the future.

Beating King Crimson: Forecasting where an orbiter is forward in time

To get this, we need to introduce two new measures in an orbit: Mean anomaly M (or ma on Ti-92) and Eccentric anomaly E (or ea).

Now, you may have heard of True anomaly. In fact, we define true anomaly. But Mean and Eccentric anomaly exist too, and they are different.

Eccentric anomaly is like True anomaly, only instead of being measured from the focus of the orbit (where the gravity well is), eccentric anomaly is measured from the center of the ellipse. 0 is still periapsis, however.

Mean anomaly is more difficult to explain. While Eccentric Anomaly is simply True Anomaly measured at the center of the ellipse rather than the focus, Mean anomaly can be thought of as the answer to the question, "If I were on a circular orbit with the same period and semimajor axis as the ellipse I'm in, where would I be on that orbit?" It is also measured so that 0 is Periapsis.

Why do we need this? Because we need Mean anomaly to shift our orbit forth, and we can't simply exchange θ into M - we need to go through E.

To start, we need to turn θ into E with this equation:

insert image of arccos((e+cos(θ)))/(1+e*cos(θ)))=ea

This is an arccos, and unfortunately, gives answers only between 0 and 180°. Now, we could do a test to see if θ is between 180° and 360°, and then act accordingly, but there's a simpler answer: make θ between -180° and 180°, and if θ is greater than either of them, shift it one complete revolution or 360° to get it within that range.

This is because it greatly simplifies matters: Instead of determining whether something is between 180° and 360°, simply see if θ is less than 0 and negate ea if it is. See? Simpler.

Turning E to M is, again, simple. You convert it with this equation:

insert image of ea-e*sin(ea)=ma

Simple, right? Well, We'll meet it again on the way back to θ. Make your judgements there.

Now we have M, and we can shift it along. So, we shift it along with:

insert image of N*t+ma=sma

Where t is how many seconds you want to shift the orbit forward. I'm going to denote this shifted Mean Anomaly with sma and Ms from now on. This is the value of M…IN THE FUTURE.

Now, we're on the way back down to get θs/sθ. Unfortunately, we have to use the same equation from earlier to get Es/sea. And it is not something you can just solve. We have to solve this using brute force. So let me introduce you to my pal, Isaac Newton.

Newton-Raphson method of iterating a problem

Now, Isaac Newton invented calculus. This is important because in doing so, he also figured out a way to find the roots of a problem. This was refined by a man named Raphson to create a way to easily find any root given a problem, a guess, and time:

insert image of x-y1(x)/d(y1(x),x)→xn

Where y1(x) is the equation you need to solve for, and that d() is the Ti-92 function for differentiation. Basically, when you take the previous iteration's answer, and subtract the function at that answer divided by the derivative of that function at that answer, you get a number that's closer to the answer than before. Repeat that enough times, and you'll arrive at the answer… eventually. This also assumes the thing in question is set equal to 0. So ea-e*sin(ea)=ma will need to be rewritten as ea-e*sin(ea)-ma=0.

Now, if you needed to do the whole 'between -180° and 180°' thing, and start at 0, you could be expected to take a while. You do have the entire angle spectrum to work with. But luckily, Ms is usually a pretty good guess of what Es is, so as your initial guess, plug sma in. And since you don't need all those decimal places to be correct, you can just do it until it's accurate enough (Usually, people cut it off when the difference between the last answer and your new answer is less than 1 x 10-6.)


So, once you get an accurate enough value for Es/sea, you can now go back to the new θs/sθ with this equation:

insert image for arccos((cos(sea)-e)/(1-e*cos(sea)))=sθ

Which, again, is an arccos equation, and again, we can use the fact that the sign for Es and θs is the same when calculating for the same point in the orbit, we negate if sea is negative.

Then we shift into the range of (-180°,180°], and we're done!


Except we're not done. This is simply the method to figure out where a ship in a given orbit will be at a given point in the future. We'll have to run through this twice - once to figure out where the interceptor is at the start of the transfer to, and again to figure out where the target is at the end of the transfer tf.

Then, we want to get where we actually are in the places we just calculated as the target's finishing point and the interceptor's starting point. Their altitudes can simply be found by putting in the polar ellipse equation:

insert 92-polar-ellipse

And their absolute longitudes (Which we need to find how much true anomaly changes between them) can be found simply by adding the longitude of periapsis and the true anomaly together at each point.

And with this, we have just finished calculating prerequisites. Now we get on to the part where things get strange instead of simply complex.

Lambert's Problem proper

So, now we get to the problem itself. Whoo!

There are many ways to solve Lambert's problem, but they all revolve around iterating it until you get a decently accurate value for what we want. in the case of Lambert's problem, this is E2, or ψ as most write it.

The way I see most to solve it (and the way I was recommended to do) is the Universal Variables method. This is useful because it can use both cartesian and polar definitions of the start and end positions. It essentially involves calculating the change in Δθ between the start of the transfer and the end (as measured according to the planet - so basically, the longitude of periapsis combined with the True anomaly), setting low and high bounds for ψ, then running though to get time, seeing if time is too high or too low, setting the current 'got' ψ as the high or low bound if it is too high or too low, respectively, then averaging ψlo and ψhi to get a new ψ, calculating new values c2 and c3, then repeating until it is - you guessed it - accurate enough. I'm not going to even try explaining it in more detail. There's some good ones From the University of Colorado that explain it well enough.

While there are faster methods to solve Lambert's Problem, the bisection method works well enough. I've got it to deliver a variable on a Ti-92, and the only thing keeping it from working on 84+ seems to be errors when I program it.

So, what are you supposed to get out, in the end? Well, either C3, or the characteristic Energy for that particular transfer, or ΔV, the delta-v (in the case of the U-Colorado paper). But since you don't need just one, you're going to need to assemble a lot of these values. This is called a Pork-Chop Plot. insert image of a Pork-Chop Plot.

A Pork-Chop Plot has its two axes labelled with the departure date and the arrival date, and in the middle is a series of contours of specific C3 or ΔV. The lower the numbers are, the better a transfer it is. (My programs don't have the ability to make Pork-Chop Plots, sadly; Contours are reserved for the third dimension on the Ti-92+, a place where you cannot take your plots, and the Ti-84+ simply can't do it.)