I've been playing around linear tangent based guidance laws (specifically the old Atlas-Centaur Surveyor guidance, Apollo IGM and Space Shuttle PEG code) and while PEG has an option for a coast phase before the final burn, I do not there's any modes where the predictor-corrector will optimize that coast. It seems to be an external parameter.

I've found papers like, e.g.:

Ping Lu, Brian J. Griffin, Gregory A. Dukeman, and Frank R. Chavez. "Rapid Optimal Multiburn Ascent Planning and Guidance", Journal of Guidance, Control, and Dynamics, Vol. 31, No. 6 (2008), pp. 1656-1664.

But ripping out the guts of the PEG algorithm and replacing it with an approach like that is going to be fairly aggressive (and I'm not sure I'd be capable of that right now). I also don't necessarily require optimality down to the last dV of fuel and simplicity and convergence would be more important to me -- so the algorithms that NASA currently uses are likely more extensive than what I'm asking for.

I suspect there's some simple late-60s/early-70s era algorithms for producing at least good-enough answers to:

  • coast-phase before final burn to orbital insertion
  • time of ignition before an on-orbit finite-thrust manuever
  • time of ignition before a descent burn to land on a body (airless)

And this is ultimately a Kerbal Space Program question. And unfortunately my last calculus of variations course was about 20 years ago, so I'm a bit rusty on math.

I vaguely understand primer vector theory, and it certainly seems like the bang-bang switching function is the thing I'd like to get to, but I don't know how to get it out of PEG, or how to bolt it onto PEG.

Some possible clarification: I don't think I'm concerned with closed-loop optimization. I suspect what I want is closer to the OPGUID algorithm and the SWITCH algorithm defined in:

Brown, K.R., Harrold, E.F., and Johnson, G.W. 1969. Rapid optimization of multiple-burn rocket flights. Technical Report NASA CR-1430, NASA, Marshall Space Flight Center.

What I don't know is if there's simpler approaches that I could use which could leverage more recent off-the-shelf algorithms (more specifically I already link against http://www.alglib.net/optimization/ for the Levenberg-Marquardt optimizer there, and it has a lot of other optimizers that might make solving the trajectory optimization problem much simpler, but I don't have the background to know what might be useful).

And I think my question comes down to if I should dive into the OPGUID and SWITCH algorithms -- and into dozens of pages of Fortran code which are older than I am -- or if there's a simpler more modern way to go about this which I'm missing? Subject to the external constraint that so far I've been more successful at implementing old NASA TDs than algorithms in the literature.

  • $\begingroup$ The clarifications are helpful! +1 $\endgroup$
    – uhoh
    Commented Aug 7, 2017 at 2:34

1 Answer 1


There are block diagrams available for PEG-4 code here which deals with onorbit and deorbit maneuvers:


It is necessary in the TURNR code to have the phi max limiting to < 45 degrees. And it is also necessary to implement the lambdaDot_xz magnitude limiting to 0.35 of the schuler frequency at the burn. There are more details in the 1979 paper on the PEG algorithm:


The Shuttle code used some subroutines called LTVCON and TRANSTIME to handle burn-coast targeting to integrate through the burn and the coast, including J2 effects of Earth's oblateness, in order to hit a r, v target.

What I've found is that it is sufficient (in KSP) to implement in the corrector:

  1. project rp into the target iy plane normally
  2. set rd = rp (no restriction on the burnout position)
  3. compute the signed angle between rd and the periapsis of the target orbit to get the true anomaly of the burnout position (periapsis cross rd should be in the -iy direction; dot product of that with -iy should be positive; or else flip the sign on the angle)
  4. set vd equal to the velocity of the target orbit at that true anomaly
  5. set the gain equal to 1.0 (since there is no integration through any coast phase)

Most of the LTVCON and TRANSTIME complexity seems to be around using the precise predictor in the shuttle code to account for J2 effects during the coast phase (or most of the coast phase -- as the burnout position changes due to PEG the time of the coast varies by delta-coast and I think that is fed through precise predict with J2=0 -- likely for either numerical efficiency reasons or numerical stability reasons?). I'm also still somewhat baffled by the C1 and C2 constants in the LTVCON routine (particularly what C1 is used for).

For the lead time I've found that simply leading the burn by tgo/2 works fine although it might not be perfectly optimal. I suspect that letting PEG run continuously and then executing the burn at J/L before the burn could be slightly more precise (in the absence of having a better calculus-of-variations-based trajectory planner).

UPDATE: I think I understand a bit better now that the lambdaDot_xz constraints in the TURNR routine are used for Abort Once Around and deorbit burns. What I wound up doing was to use the standard TURNR routine that computes rgo from rgo = rd - ( r + v * tgo + rgrav ) + rbias (which uses rd which is where I believe the position constraint feeds back into the routine) and then computes lambdaDot normally and uses the rthrust and vthrust thrust integrals. Where I ran into difficulty is that even 50-60 seconds out from the end of the burn the lambdaDot computation becomes completely unstable. I wound up clamping lambdaDot to the 0.35 * the schuler frequency value -- so the clamp is applied to lambdaDot and it goes to zero as tgo goes to zero. Forcing that clamp throughout the entire burn (as deorbit or AOA does) meant that the burnout location had a large biasing value so it was not accurate. By allowing the first part of the burn to be driven unclamped then rbias was zero and the position constraints were correct up until they were released in the last 60-ish seconds of the burn as the computed lambdaDot became unstable and exceeded the schuler frequency clamp.

  • 1
    $\begingroup$ C1 and C2 are coefficients in a linear equation which establishes the relation between the horizontal and vertical velocity: V(vert) = C2 * V(horiz) + C1. For shuttle OMS1 or OMS2 burn, the target point was always an apogee or perigee point, which implied that the vertical velocity is zero. C1 and C2 should therefore be zero. For a deorbit burn, the vertical velocity must be a negative number at El, which is the target point. For deorbit maneuvers, C1 would be a large positive number, and C2 would be between 0 and -1. A deorbit from a 150 nm orbit might have C1 = 15,434 and C2 = - 0.6200. $\endgroup$ Commented Nov 27, 2017 at 22:00
  • 1
    $\begingroup$ To be clear there's the normal terminal guidance instability that occurs with PEG and requires releasing the position clamp (rd = rp) 40 seconds before the end of the burn, and then all active guidance seems to really need to be terminated 10 seconds before the end of the burn. This was additional instability on top of that. (I also found more references on C1 and C2 and think I'm starting to understand those better) $\endgroup$
    – lamont
    Commented Dec 14, 2017 at 18:19
  • 1
    $\begingroup$ @OrganicMarble on page 40 Jaggers discusses issues in the last 40 seconds of an OMS burn in the J and Q integrals: ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/19760020204.pdf $\endgroup$
    – lamont
    Commented Dec 15, 2017 at 5:38
  • 1
    $\begingroup$ Discovered that I had accidentally removed the stop-updating-lambdaDot-in-the-last-40-seconds clamp from my code, but still doesn't quite explain lambdaDot freaking out before 40-seconds-to-go, although I'm starting to suspect that 40 seconds was tuned to the shuttle and its engines. $\endgroup$
    – lamont
    Commented Dec 15, 2017 at 18:37
  • 1
    $\begingroup$ Also ran across a statement in a paper that PEG relies on correct sensing of acceleration in the terminal phase to avoid instability, and I tracked down a bug of bugs that I had around that. I was not counting vgo down during the real terminal phase using only sensed dV and by implementing that (so i could see it running without PEG correcting it) and fixing bugs until vgo hit 0.1 when tgo hit 0.0 everything got a lot more stable and orbits got more precise. $\endgroup$
    – lamont
    Commented Dec 16, 2017 at 20:45

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.