Porkchop plots and Lambert's Problem

[Jovus]

I've taken on the project of writing an interplanetary trajectory optimization tool and a comparison of algorithm efficiency for the problem at arbitrary starting points.

Looking further into the problem, however, I have a question that I can't seem to answer.

When you optimize an interplanetary trajectory in a patched-conics approximation like KSP, how do initial and target orbit influence the problem?

Specifically, I understand how the 'interplanetary' part works. Given the position of two planets, you can calculate the orbit that will intersect one position at one time (the departure date) and the position of the other at another time (the arrival date) easily by cranking through Lambert's Problem for the solar orbit case.

However, how do you account for leaving the orbit of the start planet and arriving in orbit of the destination planet? Put another way, how do you calculate ejection angle or the optimum burn to leave/enter the patched conic gravity well?

[HebaruSan]

I did ejection burns by using an "absolute phase angle" concept to find a time when the vessel would be at a "local noon" position at the time of the burn, then adjusting based on the desired ejection angle, calculated in the usual way.

Code:

• https://github.com/HebaruSan/Astrogator/blob/master/src/TransferModel.cs#L238
• https://github.com/HebaruSan/Astrogator/blob/master/src/PhysicsTools.cs#L140

[Jovus]

Thanks!

If I'm understanding correctly, Lambert's problem applies to the trajectory between the planets, with the planets simply occupying the desired time/space points in the problem, correct? The problem of the trajectories inside the gravity wells of the planets is irrelevant to it.

[HebaruSan]

6 minutes ago, Jovus said:

If I'm understanding correctly, Lambert's problem applies to the trajectory between the planets, with the planets simply occupying the desired time/space points in the problem, correct? The problem of the trajectories inside the gravity wells of the planets is irrelevant to it.

I do not actually know what "Lambert's Problem" means (my code may or may not have anything to do with it), so I won't be a good person to answer this. Is it just how to calculate the timing for a Hohmann transfer?

[Jovus]

Yes. Lambert's problem in differential equation form is how you solve for time-of-flight for a certain orbit (or, conversely, the characteristic energy): https://en.wikipedia.org/wiki/Lambert's_problem

[VikingStormtrooper]

Dear  @Jovus, I did something similar in the past on Matlab for a university project.

The idea is that you can use the ephemerides (there is a Matlab extension incorporating them, or you could visit ssd.jpl.nasa.gov/horizons.cgi ) in order to know where every planet/moon is, in terms of orbital parameters, whenever you like. Then you select specific moments to compute the current position of your bodies, apply a Lambert problem each time, and finally build a vector of delta Vs. By selecting the lowest value in the vector, you have the cheapest transfer.

This is rather simple if you have just one Lambert arc between two bodies rotating around the same attractors; it becomes extremely complex (and a computational nightmare) if you consider also multiple gravity assists, since you have to rely on genetic algorithms and other horrible things (we had to reach Saturn with an incredibly low delta V and GAs were the only way!), but I assume you are not interested in this.

If I remember all of them, you have to worry about:

- Choosing discrete times to check for the current position of your bodies and building a huge table of orbital parameters at any time;

- The implementation of the Lambert problem (the hardest part);

- Managing orbital parameters -> position/velocity vectors conversion, and vice versa (every book of Orbital mechanics tells you);

- Managing the time law, required to relate the length and shape of your Lambert arc to the time needed to cover it (again, books are gold).

If you are stuck at some point, feel free to ask. Good luck!

[HebaruSan]

2 hours ago, VikingStormtrooper said:

Then you select specific moments to compute the current position of your bodies, apply a Lambert problem each time, and finally build a vector of delta Vs. By selecting the lowest value in the vector, you have the cheapest transfer.

So there's nothing better than an exhaustive search of the available times? That seems like the kind of thing mathematicians and physicists would try hard to improve on---is there perhaps a theorem/proof stating that's impossible?

[Steel]

24 minutes ago, HebaruSan said:

So there's nothing better than an exhaustive search of the available times? That seems like the kind of thing mathematicians and physicists would try hard to improve on---is there perhaps a theorem/proof stating that's impossible?

I doubt there's anything as concrete as a proof. The issue I guess is the fact that the system is chaotic, so it's hard to know how a small change will affect the outcome until you actually try it and see. This is a large part of the reason long range space missions use weeks of supercomputer time for trajectory analysis.

[VikingStormtrooper]

The idea is that you look at the orbital periods of your bodies, and then you select the most reasonable time discretization. Example: you are considering the Jupiter moon system, and you want to move from one of the four moons discovered by Galileo to another one (e.g. from Europa to Io); since these moons experience orbital resonance, it makes sense to consider the time interval between a certain angular configuration and the next time you have (almost) the same configuration. Then you can choose if you want to consider the intermediate positions every day, every six hours, etc.

If you do not have resonances, you can still select significant time intervals in order to have a good statistics of the possible intermediate phase angles, in such a way that your optimization covers a large variety of cases. Take the Earth-Mars transfer: launch windows occur approximately every two years, because every two years the Earth and Mars are approximately in the same relative configuration. It makes sense to subdivide this two years interval into months. If you want to be more precise, you select weeks, then days, and so on. Look for the best compromise between accuracy and computational load: you don't want the optimization to be too shallow but, at the same time, you don't want your pc to take fire after two months of wild computations just to understand that launching your spacecraft one second later grants a gain of 20 cm/s of delta V.

Edited February 19, 2017 by VikingStormtrooper

[Hannu2]

Solving the Lambert equation gives you velocity vector in solar coordinates. Ship should have that velocity at boundary of SOI (SOI can be approximated negligible small). You can then transform that location and velocity to planet's coordinates. It gives inclination and other orbital parameters. It should be quite straightforward to calculate burn which gives correct orbit from your parking orbit. The most easy and energy efficient is to use circular parking orbit with correct inclination and LAN.

[Jovus]

On 2/19/2017 at 4:38 AM, VikingStormtrooper said:

Dear  @Jovus, I did something similar in the past on Matlab for a university project.

This is in fact exactly what I am doing, except I'll probably be using C or Python, because I like punishment.

I'm certainly going to keep you on speed-dial.

Oh, and the purpose of the project is to compare efficiency of different Lambert solvers. I just want to dress it up as an actually useful tool, because I don't like making toys.

On 2/19/2017 at 0:04 PM, VikingStormtrooper said:

 Take the Earth-Mars transfer: launch windows occur approximately every two years, because every two years the Earth and Mars are approximately in the same relative configuration. It makes sense to subdivide this two years interval into months. If you want to be more precise, you select weeks, then days, and so on. Look for the best compromise between accuracy and computational load: you don't want the optimization to be too shallow but, at the same time, you don't want your pc to take fire after two months of wild computations just to understand that launching your spacecraft one second later grants a gain of 20 cm/s of delta V.

Or you can be tricksy and do an adaptive sectioning technique for optimization down to [insert reasonable section here]. For example, first you do a search on the two-year scale; then you throw out the biaunnum that is longer. Then you do a search on the six-month scale, and throw out the worst one. Then you do a search on the month scale inside the best half-year, then a week search inside the best month, etc.

Then you can draw your contours with lower specificity toward the outside, or just end up with a smaller plot and a disclaimer that 'everything outside this window is outside desired parameters'.

@VikingStormtrooper You already know this, of course; I'm just pointing it out for posterity/other readers.

[VikingStormtrooper]

Sure, you can refine the optimization algorithm as much as you like! I cited the manual process I worked on, but I am sure that it would be easy to make it automatic (probably better for your needs!). Just take into account the orbital periods of your bodies, and fix some kind of tolerance.

It's ok if you want to use C or Python; I don't know Python, but I am pretty sure that a C solver is more computationally efficient than a Matlab one (there was a paper about this, don't remember what). If you like, I can still give you my Matlab functions for reference; I did not write the Lambert solver personally, but I made the other functions. The optimization algorithm that I used was quite slow, but it worked flawlessly (I had also a gravity assist along the way, meaning two Lambert arcs and two nested for cycles, a nightmare for my CPU!). All that said, I'm waiting for your solver!  ...and don't forget to make the porkchop plots for reference, they are extremely valuable tools.

Oh, ask me if you do not know how to deal with the ephemerides from JPL Horizon (it's not always easy to set the options or understand the outputs!)

Edited February 20, 2017 by VikingStormtrooper

[Jovus]

Thanks!

Is your stuff on github? I'd love to see it.

I'm not inventing any new solvers, of course, but just doing a survey of (some of) those already developed. Specifically, my strongest reference is probably https://arxiv.org/pdf/1403.2705.pdf and its references.

Gravity assists and whatnot are out of scope - mostly because that sounds like a pain in the butt. Maybe if I get into the project and find out it's not as time-intensive as I expect, I'll build that in just to show off. It shouldn't be hard after I have the basic problem architecture, but I'm professionally lazy (which I define as doing the job right the first time, and not frittering away time I need to be doing something else on a pet project).

And yes, I'll be making porkchop plots, because I have to present my data somehow, right? Though I'll probably stick to C3 instead of resolving it into delta-V.

Finally, a pseudo-state solver is also out-of-scope, meaning my time-of-flight calculations might be off up to 30 days when we're dealing with the outer planets. Again, out of scope.

I'm not sure what you mean by the ephemerides, or the reference to JPL's Horizon?

 

Also, thanks to everyone who has contributed to this thread. You've all been immensely helpful. I'll try to keep this thread updated as I work on this project, too, since it might be of niche interest in this community.

[VikingStormtrooper]

No, I did not upload my engineering stuff on Github (but there are still my mods! ); in any case, if you think this could be useful I will take this idea into account!

That paper is nice, and features the basic notions as well! Very good reference.

Your strategy makes sense, especially if you are interested only in the basic Lambert transfers; however, if one day you feel bored, keep in mind that adding a single gravity assist is not terribly complex (although pretty time consuming, if you do not set up the algorithm efficiently).

Well, the most interesting thing to get from your computations is the minimum delta V, and the corresponding maneuver; this is just one point on a 3D plot. However, porkchops can also indicate points that are not minima in the strict sense, but still acceptable in terms of delta V or time-of-flight: this allows you to choose launch windows rather than discrete times. It's a lot of flexibility, isn't it?

Ok then!

Ephemerides are the "refined" way to deal with the problem of orbits and times. They are basically a collection of orbital parameters which allows people to know where the celestial body X is at time Y; JPL created the Horizon tool, a huge database of these parameters computed via accurate n-body simulations for a huge number of bodies. As you may imagine, it is extremely valuable for astronomers and space engineers; keep in mind also that, differently than KSP, real life is not based on patched conics, and the ephemerides allow you to track how the orbital parameters of a given body (the Earth, Titan, Mercury, etc.) change in time, depending of the overall gamma of perturbances they may experience. However, if you deal with KSP stock planets, AFAIK there are no ephemerides in-game, since KSP uses a relative time mode ("count time from the moment in which you start the match") instead of an absolute one ("minute-hour-day-month-year etc.). This means that you plan transfers basing on the current phase angle between your bodies, and then you keep the orbits of your bodies fixed in time up to moment in which you forecast a phase angle which grants a low delta V. You know how much time you have to wait, but not the day you start and the day you finish. The basic mechanism is similar, but I am personally less familiar with this one, and I prefer dealing with ephemerides.

[Jovus]

Got it. Thanks! I've seen references to the ephemerides in further papers, but never explained as well as you just did.

 

Once I'm done I'll probably post my work to github, because after all I might as well.

Unfortunately I haven't done much further work on it, because I'm focused for the moment on building a linear algebraic equation solver. (Without using the solver methods that come with most programming languages.)

[Jovus]

Been a while since I posted here, so I don't expect a response, but I thought I'd ask:

Pretty much all interplanetary launch window solvers use time-of-flight as an input, rather than an output. Does anyone know a good source for various bodies what a reasonable time-of-flight would look like? Or how to get a starting estimate?

I can find Earth/Mars data that seems to indicate a good middle value being ~200 days, and of course I would pad that out to around 100-500 to catch unlikely trajectories. I'm also sure someone somewhere has tables that'll give me reasonable values for Earth/planet trajectories. But what if I want to leave from Venus? Or Jupiter?

Maybe what I should do is take the lower bound of ToF to be half the inner body's orbital period, and the outer bound to be the outer body's orbital period, but that seems like an excessively large window.

Related, does anyone have any recommended sources for orbital parameters for the planets? Especially including true anomaly at some specified date. (I can programmatically update from there; long-term accuracy isn't a concern.)

Again, everyone, thanks for all your help.

[HebaruSan]

20 minutes ago, Jovus said:

Pretty much all interplanetary launch window solvers use time-of-flight as an input, rather than an output. Does anyone know a good source for various bodies what a reasonable time-of-flight would look like? Or how to get a starting estimate?

TransferWindowPlanner starts with a Hohmann transfer, time = π * sqrt(a³ / μ), where a = (Ap + Pe) / 2, then (I think) it checks a range of ±10%.

• https://github.com/TriggerAu/TransferWindowPlanner/blob/master/TransferWindowPlanner/SharedStuff/LambertSolver.cs#L281
• https://github.com/TriggerAu/TransferWindowPlanner/blob/master/TransferWindowPlanner/SharedStuff/LambertSolver.cs#L210

[Jovus]

1 minute ago, HebaruSan said:

TransferWindowPlanner starts with a Hohmann transfer, time = π * sqrt(a³ / μ), where a = (Ap + Pe) / 2, then (I think) it checks a range of ±10%.

• https://github.com/TriggerAu/TransferWindowPlanner/blob/master/TransferWindowPlanner/SharedStuff/LambertSolver.cs#L281
• https://github.com/TriggerAu/TransferWindowPlanner/blob/master/TransferWindowPlanner/SharedStuff/LambertSolver.cs#L210

Perfect answer. Thanks!

[HebaruSan]

5 minutes ago, Jovus said:

Perfect answer. Thanks!

Sorry, the ±10% part was wrong, as I suspected. It's actually ... easier to quote the code than summarize it :

            TravelMin = Math.Max(hohmannTransferTime - cbDestination.orbit.period, hohmannTransferTime / 2);
            TravelMax = TravelMin + Math.Min(2 * cbDestination.orbit.period, hohmannTransferTime);

• https://github.com/TriggerAu/TransferWindowPlanner/blob/master/TransferWindowPlanner/TWPWindowWorkers.cs#L62

cbDestination.orbit.period is the destination body's orbital period. So this is kind of sort of searching a full orbit of the destination body, centered around the Hohmann time.

[Jovus]

That's exactly what I was going to do, except I'll choose the farther body, whether destination or departure body, on the principle that the inner body obviously make a full revolution if the outer body does. But centering on the Hohmann time is the bit I was missing - which sounds obvious once you've said it, but that's how ideas work. Thanks again!

Edited April 24, 2017 by Jovus
My grammar's terrible when I don't get sleep

[Jovus]

Running into a weird bug with my Lambert solver.

I've construed the problem as you can see here: http://www.braeunig.us/space/interpl.htm#gauss

I'm also using their method of linear interpolation to find the correct parameter.

As construed, this method should converge for all cases except where R1 and R2 are colinear. However, that's not at all my experience. I'm clearly getting behaviour where, after a certain time-of-flight bound, the algorithm diverges (and quits with an error due to either taking the square root or a trig function of a number that is out of bounds).

If anyone feels like taking a look, you can see the code in question here: https://github.com/Jovus/Transfer-Window-Finder/blob/master/libs/scratch.py

Note I'm using cartesian coordinates, because ellipses aren't spherically symmetric. Also, I'm more than happy to give credit, despite the slapdash state of the license.

[Jovus]

OK, for a large enough time-of-flight, my algorithm transgresses the parabolic minimum or maximum (depending on which is in place). Now to track down why...

[Jovus]

Finally got it to work!

Heavily based on (and many thanks to) AlexMoon's and @TriggerAu's work in the area (https://github.com/alexmoon/ksp/blob/gh-pages/src/lambert.coffee and https://github.com/TriggerAu/TransferWindowPlanner/blob/master/TransferWindowPlanner/SharedStuff/LambertSolver.cs)

Which leaves me with a licensing question. You can see AlexMoon's original license from 2013/14 and its requirements regarding distribution and modification of the software. However, I'm technically not distributing or modifying that software in any way, but instead writing my own, in a different language, heavily based upon that original work but with some subtle-but-important differences.

Right now my work is ARR, mostly because that's easier in my time crunch than looking into what open-source license I actually want to use when I make this available. And I'm not trying to avoid giving credit where credit is due.

At the same time, to be brutally honest, I'm not sure I want to (at the moment) include a license that could be read as saying effectively, "Hey, this is someone else's work and I just cribbed," because I have to present this as an end-of-semester project, and my Professor may have strange ideas about the acceptability of using someone else's solution to the problem to help me write one, especially when they're as similar as they are.

So I'm asking y'all (and tagging @TriggerAu and would be tagging AlexMoon if I knew what his handle on the KSP forums was), because I don't want to be either a) in breach of license or b) ungrateful for the help, and I would like advice on how to proceed.

This is also the first time I've ever done something where licensing might possibly be an issue. Maybe the easiest thing would be to delete my work from GitHub for the nonce and re-upload when I've fixed the licensing? How bound am I by the original license in AlexMoon's work, since this is to my mind a heavily-inspired-but-technically-original piece? (You can see the repository in question here: https://github.com/Jovus/Transfer-Window-Finder)

[Jovus]

...okay, I'll assume silence is consent, unless I hear otherwise!

Some time soon I'll either update my version of the license, to include credit where I realize credit is due, or take my stuff off GitHub. Probably the former.

Thanks for all your help, guys!