User:0xcrab

Second stage guidance

This is where the loop is closed. The guidance system has a specific target orbit it is aiming for, which in its final form is represented by a final altitude, angular momentum vector, and vertical speed. The guidance system maintains a closed-form approximate solution to the equations of motion of powered flight, which include a guidance program, which is the guidance system's current best estimate of how the vehicle will be flown in the future. Unlike the first stage pitch program, the second stage guidance program is continually modified by the guidance system. It solves these equations of motion for the guidance program parameters, and hands those parameters to the control system to fly for a short distance. This process is repeated continuously until the end of powered flight. So, the closed loop is like this:

1. Navigate: Calculate the current state vector of the vehicle, and the osculating orbit, since it is useful for estimation:
2. Estimate: Given the current guidance system parameters, estimate the remaining time of flight, as well as certain other terminal parameters.
3. Guide: Solve the equations of motion for the guidance parameters. Start using these parameters instead of the ones previously computed.
4. Control: Steer the engines and other actuators such that the requested guidance program happens.

Back in the day, on Saturn and such, the process was divided cleanly into two loops, which ran effectively as two separate processes. The navigate/estimate/guide sequence forms one loop, which runs at a relatively low frequency, since in the grand scheme of things, the vehicle just doesn't move very far in a millisecond. Typical frequencies are around 2sec, as used on the Saturn launch vehicle. This was run on the IU Digital Computer in the background. This is called the major loop. It took two seconds because that's just how slow computers used to be. It took most of the computer power of the digital computer for the entire two seconds to perform the calculation.

The control loop runs at a much higher frequency, or in an analog system, continuously. It takes commands from the major loop, which are in the form of a desired pitch and pitch rate, and implements them by steering the engines. On Saturn, this was run on the IU Digital Computer as an interrupt task triggered by a timer running at a fairly high rate. This is the minor loop. The minor loop is responsible for balancing the rocket on its end. Most rockets are aerodynamically unstable, or stable flying in reverse. The control loop keeps the rocket flying straight. The minor loop was kept very simple in order to run quickly enough.

In Orbiter, it is relatively easy to separate the major and minor loops. The major loop is tied to the MFD update message, which runs once per second, and the analog loop is attached to the timestep message, which runs once per frame. Alternatively, the minor loop runs every frame, but keeps count of when to run the major loop again. Both loops run so fast on modern computers that it is not necessary to break up the major loop and run it in the background or in pieces.

Navigation

In all math to follow, the arbitrary vector ${\vec {y}}$ is drawn with an arrow over it. This vector's (scalar) length is $y$ and direction is ${\hat {y}}$ so the equation ${\vec {y}}=y{\hat {y}}$ is always true for any vector. Also, functions are written as $f[t]$ instead of $f(t)$ because there are lots of parentheses already.

The job of the navigation system is to calculate the state vector, as well as any auxiliary parameters that the estimator or guidance systems will need. Currently this includes:

Vehicle state

A real navigation system calculates its state vector based on watching the inertial measurement unit (IMU) and integrating acceleration or velocity to get the current position. In Orbiter, we just call the appropriate API function. In any case we end up with the following vectors, in any convenient coordinate system.

${\vec {r}}$ is the position vector relative to the center of the central body.

${\vec {v}}$ is the inertial velocity vector relative to the center of the central body.

Time

The time epoch $t\equiv 0$ is Right Now, the instant that the major loop is being run. Nothing is less relevant in aviation than the amount of runway behind you, the amount of air above you, and the amount of fuel you left at the airport. Similarly, nothing is less relevant in explicit guidance than what has already happened. Only where the vehicle is now, matters. Each iteration considers the current time and vehicle state to be an initial condition.

In the major loop, the only other time which is important in this derivation is $T$, burnout time, the time remaining until the target orbit is reached and the rocket is shut down. This is not necessarily when the rocket runs out of gas. Burnout time is continuously calculated in the estimation phase. $T$ is not constant during the flight. It always counts the time from the present until burnout. In particular, if it is perfectly estimated, it will decrease by the length of the major cycle each major cycle.

In the minor loop, $t$ measures the time since the current steering constants became valid, and can be greater than zero.

We have to be careful with time, because for any particular major cycle, greater $t$ means towards the future. However, when considering many consecutive major cycles, the current time $t$ is always zero, and $T$ decreases, in effect counting down to cutoff.

Basis vectors

For convenience, several basis vectors are calculated, forming a local horizontal coordinate frame. These basis vectors themselves are expressed in the same frame as ${\vec {r}}$ and ${\vec {v}}$.

${\hat {r}}\equiv {\frac {\vec {r}}{r}}$ (1a) Vertical direction

${\vec {h}}\equiv {\vec {r}}\times {\vec {v}}$ Angular momentum vector

${\hat {h}}\equiv {\frac {\vec {h}}{h}}={\frac {{\vec {r}}\times {\vec {v}}}{\left|{\vec {r}}\times {\vec {v}}\right|}}$ (1b) Crosstrack direction

${\vec {\theta }}\equiv {\vec {h}}\times {\vec {r}}$ (1c) Downtrack direction (parallel to the local horizontal)

${\hat {f}}$ is the direction the thrust vector is pointed

$\omega \equiv {\frac {{\vec {v}}\cdot {\hat {\theta }}}{r}}$ is the angular velocity, defined by the current tangential component of velocity ${\vec {v}}\cdot {\hat {\theta }}$ (Note the hats and vecs!) divided by the current radius

This is a rotating frame, which therefore has centrifugal force.

Vehicle Measurement

The navigation calculations need several quantities related to the vehicle performance, such as vehicle mass $m$, thrust magnitude $f$, and specific impulse $v_{e}$ (calculated with calls to the Orbiter API)

The quantity $I_{sp}$ is the specific impulse, and is the amount of force obtained per unit of propellant flow. It has a dimension of ${\frac {({\mbox{Force}})({\mbox{Time}})}{\mbox{Mass}}}$, ${\frac {{\mbox{N}}\cdot {\mbox{s}}}{\mbox{kg}}}$ in SI units.

The quantity $v_{e}$ is the effective exhaust velocity. This is related to the specific impulse as follows: For a rocket, specific impulse and exhaust velocity are equivalent concepts. The units of specific impulse happen to cancel out into velocity

${\frac {{\mbox{N}}\cdot {\mbox{s}}}{\mbox{kg}}}$	$={\frac {\not {\mbox{kg}}\cdot {\mbox{m}}}{{\mbox{s}}^{2}}}{\frac {\mbox{s}}{\not {\mbox{kg}}}}$
	$={\frac {\mbox{m}}{{\mbox{s}}^{\not {2}}}}\not {\mbox{s}}$
	$={\frac {\mbox{m}}{\mbox{s}}}$

and this velocity is the effective exhaust velocity of the rocket $v_{e}$. However, in the literature, specific impulse is most often expressed in "seconds", by which is meant ${\frac {\mbox{kgf s}}{\mbox{kg}}}$. The force units on top are cancelled with the mass units on the bottom, even though this is not strictly correct. In this form, it says how long you can burn one mass unit of propellant to support one mass unit in standard gravity. Multiply this by the standard acceleration of gravity at Sea Level $g_{0}=9.80665{\mbox{m}}/{\mbox{s}}^{2}$ (Not $g$ at the current location of the rocket!) to get the exhaust velocity.

$v_{e}=I_{sp}g_{0}$

Wikipedia recommends this terminology to avoid confusion, so $v_{e}$ is in meters per second, and $I_{sp}$ is in "seconds". I will use these symbols as recommended, but I will call them both specific impulse, since specific impulse and effective exhaust velocity are both exactly the same.

In a real vehicle, acceleration is measured directly from the IMU, and specific impulse is treated as a constant. Usually no attempt is made to measure thrust and vehicle mass directly, since acceleration is measured so much more accurately.

In Orbiter, the reverse is true, and vehicle thrust, mass, and specific impulse are available and completely accurate. In Orbiter, the definition of specific impulse used is consistently effective exhaust velocity, $v_{e}$, measured in meters per second.

Rocket equations

The rocket equations define the ideal acceleration of the rocket under current conditions, at any time in the forseeable future, but not across staging events, or similar changes in vehicle performance.

$\tau ={\frac {v_{e}}{a}}$ (3b) A sort of normalized mass - time to burn the vehicle completely as if it were all propellant

$a[t]={\frac {a}{1-{\frac {t}{\tau }}}}$ (3a) Acceleration at any time in the future. This goes to infinity as $t$ approaches $\tau$

You also need certain integrals, here expressed as functions. The first integral is the integrated acceleration from now till burnout, otherwise known as $\Delta v$ or ideal change in speed from now to then. Ideal in the sense that this is the velocity change if the acceleration is continually applied in the same direction, not considering outside forces, relative to an unpowered body sharing the initial state. It is fairly simple to calculate:

$b_{0}[T]=\int \limits _{0}^{T}a(t)dt=-v_{e}\ln \left(1-{\frac {T}{\tau }}\right)=\Delta v$ (7a)

So basically this part is saying that we will use the fraction $A$ of our acceleration to go straight up and hit our vertical velocity target

The second integral isn't quite as intuitive. It is the first moment of $\Delta v$. The clever guys who came up with PEG also came up with a recursive formula for the nth moment of $\Delta v$ because we will end up needing higher moments later:

$b_{n}[T]=\int _{0}^{T}t^{n}a(t)dt=b_{n-1}\tau -{\frac {v_{e}T^{n}}{n}}$ (7b)

The terms $c_{n}$ are the moments of the ideal distance traveled under power with acceleration in a constant direction, neglecting outside forces, relative to an unpowered body with the same initial state. This third important integral is

$c_{0}[T]=\int _{0}^{T}\int _{0}^{t}a(s)dsdt=b_{0}T-b_{1}$ (7c)

and the fourth is (in nth moment form)

$c_{n}[T]=\int _{0}^{T}\int _{0}^{t}s^{n}a(s)dsdt=c_{n-1}\tau -{\frac {v_{e}T^{n+1}}{n(n+1)}}$ (7d)

Cutoff logic

Even more disturbing, the Titan had a huge range of ballistic trajectories programmed into its guidance computer, which aimed the missile below the horizon if it was headed for a military target or above the horizon if it was headed for space. As the rocket rose, the computer would continually hunt for just the right orientation, causing the missile to wiggle its nose up and down and left and right, bloodhound-fashion, sniffing for a target that might be Moscow, might be Minsk, or might be low Earth orbit, depending upon whether it was carrying warheads or spacemen on that particular mission.

--Lost Moon, p41

As stated in the original document:

As cutoff energy is approached, the equations for calculating $A$, $B$, $T$ and $K$ become indeterminate. If these calculations are continued up to cutoff, the values begin to oscillate and eventually diverge.

The closer and closer you get, the more sensitive the guidance equations to calculate $A$ and $B$ become. Imagine you are headed for a 185km target altitude, zero vertical speed, and the vehicle is 184.9km in altitude, zero vertical speed, with 10 seconds left. In order to zero out this error, the guidance equations will attempt to make $A$ and $B$ really big, so that the vehicle will gain an average of 10m/s, but start and end at zero, and do it all in 10 seconds. This will probably cause ${\hat {f}}\cdot {\hat {r}}$ to go out of range, and cause the time estimator to be invalid. Since it is linearized, it will come up with an answer, but the answer will be well off of what is needed, and $T$ will be wildly wrong also.

I am pretty sure that what is going on is this: The system at every major cycle calculates a course which should bring it near the target condition. However, since there are approximations, even if the vehicle flies this course perfectly, it will not arrive at the desired end state. If there were no approximations and the vehicle flew the target course perfectly, each new major cycle would just calculate new steering which perfectly match the updated old steering constants. Since there are approximations, since the vehicle does not fly the plotted course perfectly, and since this is a converging process in the first place, the new steering constants do not match the updated old ones. Apparently, as the system gets closer to cutoff, the errors cause larger and larger corrections to the steering constants. Eventually the corrections get so big that the control system can't keep up and achieve the commanded rotation rates, the vehicle really starts to fly off course, and the problem gets progressively worse and worse until the vehicle spins out of control.

The original document suggests halting calculations of $A$ and $B$ about 10 seconds before the final state is reached, and just fly with the $A$ and $B$ that comes out of that last iteration.

The system should stop calculating $A$ and $B$ about 7 seconds before cutoff. Actual cutoff time is estimated to the end using the $\Delta v$ formula. Every major cycle, the system calculates the time of cutoff, simply by adding $T$ to the guidance start time. Most of the time this is far in the future, but eventually the time comes when $T$ is less than 1 major cycle, and one of the minor cycles will issue the cutoff command before the major cycle comes back around.

Parking orbit targeting

For typical modern launch vehicles, the boosters, first stage, and some of the second stage is used to put the remainder of the second stage and the payload into a low-earth parking orbit. For many spacecraft, this is the final destination, in which case probably all of the second stage is used up. Other spacecraft may be headed to circular orbits at other altitudes, in which case the remainder of the second stage is used after a coast period to put the spacecraft onto a transfer orbit to the target orbit. The Space Shuttle is put into an initial elliptical transfer orbit with its perigee in the atmosphere, such that the external tank is disposed of. It uses its own OMS engines to circularize the orbit at first apogee. Missions to the moon, planets, and beyond, are put into a highly elliptical or hyperbolic departure trajectory by using the rest of the second stage, and possibly upper stages, at the appropriate moment.

For such missions, the only targeting needed is for the initial parking orbit. Launch vehicles typically aim for a burnout altitude of about 185km, zero vertical speed, and sufficient horizontal speed to be in a circular orbit. The Space Shuttle has a burnout altitude much lower, about 110km, with a significant upward vertical speed, and a horizontal speed such that the first apogee is at the target orbit. The orbit plane normal vector may take on any desired value (with some restrictions) by launching at the right time to the right azimuth.

For such missions, the three target parameters are given, not calculated, and what we have above is enough to achieve the target orbit.

Major loop algorithm

• Guide, based on old (or first-guestimate) burn time remaining T and burnout state target radius $r_{T}$, vertical speed ${\dot {r}}_{T}$, and horizontal speed $v_{\theta ,T}$

  ${\dot {r}}={\frac {{\vec {r}}\cdot {\vec {v}}}{r}}$

  $b_{0}=b_{0}(T,v_{e},\tau )$

  $b_{1}=b_{1}(T,v_{e},\tau )$

  $c_{0}=c_{0}(T,v_{e},\tau )$

  $c_{1}=c_{1}(T,v_{e},\tau )$

  ${\vec {M}}_{B}={\begin{bmatrix}{\dot {r}}_{T}-{\dot {r}}\\r_{T}-r-{\dot {r}}T\end{bmatrix}}$

  $[\mathbf {M} _{A}]={\begin{bmatrix}b_{0}&b_{1}\\c_{0}&c_{1}\end{bmatrix}}$

  Use your favorite matrix solver to solve $[\mathbf {M} _{A}]{\vec {M}}_{x}={\vec {M}}_{B}$ (where ${\vec {M}}_{x}={\begin{bmatrix}A\\B\end{bmatrix}}$) for ${\vec {M}}_{x}$

• Estimate burn time remaining based on guidance ${\vec {M}}_{x}={\begin{bmatrix}A\\B\end{bmatrix}}$, old burn time estimate T and time since last estimate $\Delta t$
  • Time update of old guidance

    $A^{-}=A+\Delta tB$

    $B^{-}=B$

    $T^{-}=T-\Delta t$

  • Navigation, basis vectors, and additional target conditions - presume ${\vec {r}}$ and ${\vec {v}}$ are available

    ${\vec {h}}={\vec {r}}\times {\vec {v}}$

    ${\hat {\theta }}={\hat {h}}\times {\hat {r}}$

    $h_{T}=v_{\theta ,T}r_{T}$ Final angular momentum magnitude

    $\Delta h=h_{T}-h$ Angular momentum to gain

    ${\bar {r}}={\frac {r_{T}+r}{2}}$Mean radius

    ${\dot {r}}={\vec {v}}\cdot {\hat {r}}$

    $v_{\theta }={\vec {v}}\cdot {\hat {\theta }}$

  • Vehicle performance

    $\tau ={\frac {v_{e}}{a_{0}}}$

    $C={\frac {{\frac {\mu }{r^{2}}}-\omega ^{2}r}{a_{0}}}$ Portion of vehicle acceleration used to counteract gravity and centrifugal force

    $f_{r}=A^{-}+C$ sin(pitch) at current time

  • Estimation

    $C_{T}={\frac {{\frac {\mu }{r_{T}^{2}}}-\omega _{T}^{2}r_{T}}{a[T^{-}]}}$ Gravity and centrifugal force term at cutoff

    $f_{r,T}=A^{-}+B^{-}T^{-}+C_{T}$ sin(pitch) at burnout

    ${\dot {f}}_{r}={\frac {f_{r,T}-f_{r}}{T^{-}}}$Approximate rate of sin(pitch)

    $f_{h}=0$ (Skipped) yaw term

    ${\dot {f}}_{h}=0$ (Skipped) yaw rate term

    $f_{\theta }=1-{\frac {f_{r}^{2}}{2}}-{\frac {f_{h}^{2}}{2}}$ Approximate cos(pitch)

    ${\dot {f}}_{\theta }=-(f_{r}{\dot {f}}_{r}+f_{h}{\dot {f}}_{h})$ Approximate cos(pitch) speed

    ${\ddot {f}}_{\theta }=-{\frac {{\dot {f}}_{r}^{2}+{\dot {f}}_{h}^{2}}{2}}$ Approximate cos(pitch) acceleration

    $\Delta v={\frac {{\frac {\Delta h}{\bar {r}}}+v_{e}T({\dot {f}}_{\theta }+{\ddot {f}}_{\theta }\tau )+{\frac {{\ddot {f}}_{\theta }v_{e}T^{2}}{2}}}{f_{\theta }+{\dot {f}}_{\theta }\tau +{\ddot {f}}_{\theta }\tau ^{2}}}$ Ideal velocity-to-gain in order to get target angular momentum, based on current pitch guidance

    $T^{+}=\tau \left(1-exp\left[-{\frac {\Delta v}{v_{e}}}\right]\right)$ Estimate updated burnout time

To converge the guidance originally, run guidance with a guesstimate of stage burnout, then run the estimate with a Δt of zero. Run guidance and estimation repeatedly until things converge. Real launch vehicles spread this convergence over several cycles, but the computer on your desk running Orbiter is more powerful than any launch vehicle guidance computer, so you can probably run all this guidance simultaneously in a single frame without noticeably slowing down Orbiter.

After guidance has converged, you probably want to run estimate first, then immediately followed by guidance, then pause for the duration of the major loop, then estimate, guide, pause, etc. Each time, the length of the pause will be the Δt used in the estimate routine.

Extension to multiple stages

It isn't practical, until we have something like a fusion-powered DeltaGlider, to get into orbit on a single stage. Therefore, at some point in the ascent, the configuration of the launch vehicle will radically change and must be accounted for.

Frequently, however, we can cheat. The first stage usually runs on first stage guidance, and we only need to start explicit guidance once we have ditched the booster. Therefore, two stage rockets like the Saturn IB and Titan II, and stage-and-a-half vehicles like Mercury Atlas and the space shuttle, only need to worry about a single stage.

Other vehicles, like the Delta II, Titan IV, Saturn V, and Atlas-Centaur, however, have multiple stages even after the boosters are dropped. The middle stage runs for too long to disregard guidance. So, somehow, the guidance program needs to know about and account for the performance change.

Usually lower stages are burned to depletion, and only the last stage is burned to a target orbit. Based on the rocket performance, we then can calculate, and don't have to estimate, the burn time of the lower stages.

The pitch program was designed with the help of a tool not available to me, variational calculus. Further analysis using that tool shows that during a staging, or any other step change in acceleration or specific impulse, the thrust direction ${\hat {f}}$ and rate ${\dot {\hat {f}}}$ should be continuous. However, the acceleration is discontinuous. Therefore, the component of thrust necessary to counter gravity and centrifugal force will change discontinuously, and therefore the steering constants $A$ and $B$ need to be discontinuous as well.

Stage designation

Stages are numbered such that the first PEG-guided stage is $\ell =1$, the second is $\ell =2$, etc. This stage number doesn't change as stages drop, so stage 2 is always stage 2, even after stage 1 has gone its way. The number $\ell$ then is the name of any particular stage. It's still numeric, but typically used as an array index or for-loop variable. Also, there is nothing sacred about the name $\ell$. It does not represent the currently running stage, or any particular stage. In fact, the number of the currently running stage is frequently $j$.

With this terminology, on an Atlas-Centaur, stage 1 refers to the Atlas sustainer after booster drop while stage 2 refers to the Centaur. On a Saturn V, stage 1 is the S-II stage (nominally second stage), and stage 2 is the S-IVB (nominally third stage). On a Titan-Centaur, which flies its first stage guidance while only the boosters are firing, stage 1 is Titan 1, 2 is Titan 2, and 3 is Centaur.

As you can see with the Titan above, we are very careful in designating stages so that we may write a generic guidance program which can deal with any number of stages, not just one or two.

Stage timing

As before, we have a progress variable $t$ such that $t=0$ is Right Now and $t=T$ is the time of cutoff, at achievement of the target parameters.

With multiple stages, other instants become important. We define $T_{\ell }$ as the remaining burn duration of each stage. On a currently burning stage, $T_{\ell }$ decreases each major cycle, while any higher stage(s) are constant until it is its/their turn to fire. In order to eliminate certain corner cases, we define $T_{0}\equiv 0$. During an upper stage burn, the lower stage(s) have its/their $T_{\ell }$ set to zero. This way, $T=\sum _{\ell =0}^{N}T_{\ell }$, or just $T=T_{1}+T_{2}$ for a two-stage vehicle.

Stages normally have a short coast period between shutdown of the lower stage engine and start of the upper stage. This coast period is disregarded. The staging event is considered instantaneous, with one stage burning up to the instant $T_{\ell }$ is zero, and the next stage starting up instantly at full rated thrust in the next instant. This causes several effect:

• The guidance program must know what to do when the vehicle is in fact coasting, and must not crash itself or fill feedback variables with NaN or other such problem
• $T_{\ell }$ can never be zero or negative for any running stage. The original document makes sure $T_{\ell }$ is always at least 1 second for the running stage.
• It doesn't matter where in the staging event the instant of $T_{\ell }$ is, just that it is while neither stage is firing.

On a two-stage rocket, predicted duration of first-stage firing is a design parameter of the rocket which the guidance program must be told. The duration of the second stage is calculated by guidance itself. Conversely, performance of the first stage is directly measured during first stage flight, while the performance quantities $a(t_{\ell }=0)$ (or $\tau$) and $v_{e}$ for the upper stage must be provided from the vehicle design. On a three-or-more stage rocket, both timing and performance constants must be provided from vehicle design for each intermediate stage. The performance of all stages but the currently firing one must be provided, along with the burn time of all the stages but the last.

It is sometimes convenient to have secondary stage timeline variables $t_{\ell }\equiv t-\sum _{j=0}^{\ell -1}T_{j}$. These are zero at the start of their respective stages, and the $t_{\ell }$ of the currently running stage is always zero at in the major cycle.

The calculus of variations shows that the thrust direction and rate ${\hat {f}}$ and ${\hat {\dot {f}}}$ should be continuous across staging for optimum performance. Since the acceleration will be discontinuous at staging points, A and B must also be discontinuous in order to satisfy the continuity in ${\hat {f}}$ and ${\hat {\dot {f}}}$.

Algorithm

• If first time around, then

  • Guess $A_{1},B_{1},T_{2}$

  else

  • Update $A_{1},B_{1},T_{1}$ given $\Delta T$ since last iteration (can be zero)

    $A_{1}=A_{1}+B_{1}\Delta T$

    $T_{1}=T_{1}-\Delta T$

  endif

• Vertical state at staging, given guesstimate of first stage guidance $A_{1},B_{1}$. Calculating the powered flight integrals requires the vehicle performance and burnout time for stage 1 only.

  ${\dot {r}}[T_{1}]={\dot {r}}[0]+b_{0,1}[T_{1}]A_{1}+b_{1,1}[T_{1}]B_{1}$

  $r[T_{1}]=r[0]+{\dot {r}}[0]T_{1}+c_{0,1}[T_{1}]A_{1}+c_{1,1}[T_{1}]B_{1}$

• Current required pitch, given guidance $A_{1},B_{1}$ and current state

  $f_{r,1}[0]=A_{1}+{\frac {{\frac {\mu }{r_{0}^{2}}}-\omega _{0}^{2}r_{0}}{a_{1}[0]}}$

• Horizontal state at staging, given same course $A_{1},B_{1}$. Need an initial guess on $\omega (T_{1})$, iterate until convergence

  $f_{h,1}={\dot {f}}_{h,1}=0$ until we feel like calculating yaw

  $f_{r,1}[T_{1}]=A_{1}+B_{1}T_{1}+{\frac {{\frac {\mu }{r[T_{1}]^{2}}}-\omega [T_{1}]^{2}r[T_{1}]}{a_{1}[T_{1}]}}$

  ${\dot {f}}_{r,1}={\frac {f_{r,1}[T_{1}]-f_{r,1}[0]}{T_{1}}}$

  $f_{\theta ,1}=1-{\frac {f_{r,1}^{2}}{2}}-{\frac {f_{h,1}^{2}}{2}}$

  ${\dot {f}}_{\theta ,1}=-(f_{r}{\dot {f}}_{r,1}+f_{h}{\dot {f}}_{h,1})$

  ${\ddot {f}}_{\theta ,1}=-{\frac {{\dot {f}}_{r,1}^{2}+{\dot {f}}_{h,1}^{2}}{2}}$

  $h[T_{1}]={\frac {r[0]+r[T_{1}]}{2}}\left(f_{\theta ,1}b_{0,1}[T_{1}]+{\dot {f}}_{\theta ,1}b_{1,1}[T_{1}]+{\ddot {f}}_{\theta ,1}b_{2,1}[T_{1}]\right)$

  $v_{\theta }[T_{1}]={\frac {h[T_{1}]}{r[T_{1}]}}$ (42d)

  $\omega [T_{1}]=v_{\theta }[T_{1}]/r[T_{1}]$

• Guidance staging discontinuities

  $\Delta A=\left({\frac {\mu }{r[T_{1}]^{2}}}-\omega [T_{1}]^{2}r[T_{1}]\right)\left({\frac {1}{a_{1}[T_{1}]}}-{\frac {1}{a_{2}[0]}}\right)$

  $\Delta B=-\left({\frac {\mu }{r[T_{1}]^{2}}}-\omega [T_{1}]^{2}r[T_{1}]\right)\left({\frac {1}{v_{e,1}}}-{\frac {1}{v_{e,2}}}\right)+\left(3\omega [T_{1}]^{2}-{\frac {2\mu }{r[T_{1}]^{3}}}\right){\dot {r}}[T_{1}]\left({\frac {1}{a_{1}[T_{1}]}}-{\frac {1}{a_{2}[0]}}\right)$

• Solve Explicit Guidance Equations. An estimated stage 2 burn time $T_{2}$ is needed

  $[\mathbf {M} _{A}]={\begin{bmatrix}M_{A11}&M_{A12}\\M_{A21}&M_{A22}\end{bmatrix}}={\begin{bmatrix}b_{0,1}[T_{1}]+b_{0,2}[T_{2}]&b_{1,1}[T_{1}]+b_{1,2}[T_{2}]+b_{0,2}[T_{2}]T_{1}\\c_{0,1}[T_{1}]+c_{0,2}[T_{2}]+b_{0,1}[T_{1}]T_{2}&c_{1,1}[T_{1}]+b_{1,1}[T_{1}]T_{2}+c_{0,2}[T_{2}]T_{1}+c_{1,2}[T_{2}]\end{bmatrix}}$

  ${\vec {M}}_{B}={\begin{bmatrix}M_{B1}\\M_{B2}\end{bmatrix}}={\begin{bmatrix}{\dot {r}}[T]-{\dot {r}}[0]-b_{0,2}[T_{2}]\Delta A-b_{1,2}[T_{2}]\Delta B\\r[T]-r[0]-{\dot {r}}[0](T_{1}+T_{2})-c_{0,2}[T_{2}]\Delta A-c_{1,2}[T_{2}]\Delta B\end{bmatrix}}$

  Use your favorite matrix solver to solve $[\mathbf {M} _{A}]{\vec {M}}_{x}={\vec {M}}_{B}$ (where ${\vec {M}}_{x}={\begin{bmatrix}A_{1}\\B_{1}\end{bmatrix}}$) for ${\vec {M}}_{x}$

• Estimate stage 2 burn time $T_{2}$, using the same method as for the single stage, but from staging state to target state

  $A_{2}=\Delta A+A_{1}+B_{1}T_{1}$

  $B_{2}=\Delta B+B_{1}$

  $\tau _{2}={\frac {v_{e,2}}{a_{2}[0]}}$

  $h[T_{1}]=r[T_{1}]v_{\theta }[T_{1}]$ Specific angular momentum at staging

  $h[T]=r[T]v_{\theta }[T]$

  $\Delta h=h[T]-h[T_{1}]$Angular momentum to gain

  ${\bar {r}}={\frac {(r[T]+r[T_{1}])}{2}}$

  $C_{2}[0]={\frac {{\frac {\mu }{r[T_{1}]^{2}}}-\omega [T_{1}]^{2}r[T_{1}])}{a_{2}[0]}}$ Portion of vehicle acceleration used to counter gravity and centrifugal force

  $f_{r,2}[0]=A_{2}+C_{2}$ Staging sin(pitch)

  $\omega [T]={\frac {v_{\theta }[T]}{r[T]}}$ Target angular speed

  $C_{2}[T]={\frac {{\frac {\mu }{r[T]^{2}}}-\omega [T]^{2}r[T])}{a_{2}[T_{2}]}}$ Gravity and centrifugal force term at target state

  $f_{r,2}[T]=A_{2}+B_{2}T_{2}+C_{2}[T]$ sin(pitch) at target state

  ${\dot {f}}_{r,2}={\frac {f_{r,2}[T]-f_{r,2}[0]}{T2}}$ Approximate speed of sin(pitch)

  $f_{h,2}=0$ Still a placeholder for yaw steering

  ${\dot {f}}_{h,2}=0$

  $f_{\theta ,2}=1-{\frac {f_{r,2}^{2}}{2}}-{\frac {f_{h,2}^{2}}{2}}$ Staging cos(pitch)

  ${\dot {f}}_{\theta ,2}=-(f_{r,2}{\dot {f}}_{r,2}+f_{h,2}{\dot {f}}_{h,2})$ Approximate speed of cos(pitch)

  ${\ddot {f}}_{\theta ,2}=-{\frac {({\dot {f}}_{r,2}^{2}+{\dot {f}}_{h,2}^{2})}{2}}$ Approximate acceleration of cos(pitch)

  $\Delta V={\frac {{\frac {\Delta h}{\bar {r}}}+v_{e,2}T_{2}({\dot {f}}_{\theta ,2}+{\ddot {f}}_{\theta ,2}\tau _{2})+{\frac {{\ddot {f}}_{\theta ,2}v_{e,2}T_{2}^{2}}{2}}}{f_{\theta ,2}+{\dot {f}}_{\theta ,2}\tau +{\ddot {f}}_{\theta ,2}\tau _{2}^{2}}}$Required DeltaV to gain target angular momentum

  $T_{2}=\tau _{2}\left(1-\exp \left[-{\frac {\Delta V}{v_{e,2}}}\right]\right)$