Kevin E. Witzberger

Glenn Research Center, Cleveland, Ohio

Tom Zeiler

The University of Alabama, Tuscaloosa, Alabama

Direct Method Transcription For a Human-Class 

Translunar Injection Trajectory Optimization

NASA/TM—2012-217215

January 2012

AAS 11–124
NASA STI Program . . . in Profile

Since its founding, NASA has been dedicated to the 

advancement of aeronautics and space science. The 

NASA Scientific and Technical Information (STI) 

program plays a key part in helping NASA maintain 

this important role.

The NASA STI Program operates under the auspices 

of the Agency Chief Information Officer. It collects, 

organizes, provides for archiving, and disseminates 

NASA’s STI. The NASA STI program provides access 

to the NASA Aeronautics and Space Database and 

its public interface, the NASA Technical Reports 

Server, thus providing one of the largest collections 

of aeronautical and space science STI in the world. 

Results are published in both non-NASA channels 

and by NASA in the NASA STI Report Series, which 

includes the following report types:

• TECHNICAL PUBLICATION. Reports of 

completed research or a major significant phase 

of research that present the results of NASA 

programs and include extensive data or theoretical 

analysis. Includes compilations of significant 

scientific and technical data and information 

deemed to be of continuing reference value. 

NASA counterpart of peer-reviewed formal 

professional papers but has less stringent 

limitations on manuscript length and extent of 

graphic presentations.

• TECHNICAL MEMORANDUM. Scientific 

and technical findings that are preliminary or 

of specialized interest, e.g., quick release 

reports, working papers, and bibliographies that 

contain minimal annotation. Does not contain 

extensive analysis.

• CONTRACTOR REPORT. Scientific and 

technical findings by NASA-sponsored 

contractors and grantees.

• CONFERENCE PUBLICATION. Collected 

papers from scientific and technical 

conferences, symposia, seminars, or other 

meetings sponsored or cosponsored by NASA.

• SPECIAL PUBLICATION. Scientific, 

technical, or historical information from 

NASA programs, projects, and missions, often 

concerned with subjects having substantial 

public interest.

• TECHNICAL TRANSLATION. Englishlanguage translations of foreign scientific and 

technical material pertinent to NASA’s mission.

Specialized services also include creating custom 

thesauri, building customized databases, organizing 

and publishing research results.

For more information about the NASA STI 

program, see the following:

• Access the NASA STI program home page at 

http://www.sti.nasa.gov

• E-mail your question via the Internet to help@

sti.nasa.gov

• Fax your question to the NASA STI Help Desk 

at 443–757–5803

• Telephone the NASA STI Help Desk at

443–757–5802

• Write to:

 NASA Center for AeroSpace Information (CASI)

 7115 Standard Drive

 Hanover, MD 21076–1320
Kevin E. Witzberger

Glenn Research Center, Cleveland, Ohio

Tom Zeiler

The University of Alabama, Tuscaloosa, Alabama

Direct Method Transcription For a Human-Class 

Translunar Injection Trajectory Optimization

NASA/TM—2012-217215

January 2012

AAS 11–124

National Aeronautics and

Space Administration

Glenn Research Center

Cleveland, Ohio 44135

Prepared for the

2011 Space Flight Mechanics Meeting

sponsored by the American Astronautical Society (AAS)

New Orleans, Louisiana, February 13–17, 2011
Acknowledgments

The authors would like to acknowledge Dr. Rajnish Sharma for his review and comments.

Available from

NASA Center for Aerospace Information

7115 Standard Drive

Hanover, MD 21076–1320

National Technical Information Service

5301 Shawnee Road

Alexandria, VA 22312

Available electronically at http://www.sti.nasa.gov

Level of Review: This material has been technically reviewed by technical management. 
NASA/TM—2012-217215 1

Direct Method Transcription For a Human-Class 

Translunar Injection Trajectory Optimization 

Kevin E. Witzberger 

National Aeronautics and Space Administration 

Glenn Research Center 

Cleveland, Ohio 44135 

Tom Zeiler 

The University of Alabama 

Tuscaloosa, Alabama 35487 

Abstract 

This paper presents a new trajectory optimization software package developed in the framework of a 

low-to-high fidelity 3 degree-of-freedom (DOF)/6-DOF vehicle simulation program named Mission 

Analysis Simulation Tool in Fortran (MASTIF) and its application to a translunar trajectory optimization 

problem. The functionality of the developed optimization package is implemented as a new “mode” in 

generalized settings to make it applicable for a general trajectory optimization problem. In doing so, a 

direct optimization method using collocation is employed for solving the problem. Trajectory 

optimization problems in MASTIF are transcribed to a constrained nonlinear programming (NLP) 

problem and solved with SNOPT, a commercially available NLP solver. A detailed description of the 

optimization software developed is provided as well as the transcription specifics for the translunar 

injection (TLI) problem. The analysis includes a 3-DOF trajectory TLI optimization and a 3-DOF vehicle 

TLI simulation using closed-loop guidance. 

Introduction 

Ever since the National Aeronautics and Space Administration’s (NASA’s) Apollo program ended in 

the early 1970’s, humans have not traveled beyond Earth orbit (BEO). Various issues have hampered 

BEO travel for humans including the fact that human space flight (HSF) is costly and requires continued 

commitment from each new administration. HSF is also risky, as evidenced in 2003 with the Space 

Shuttle’s Columbia disaster. In response to the Columbia Accident Investigation Report (CAIB), in 

January 2004, President George W. Bush announced a new United States space policy—Vision For Space 

Exploration (VSE)(Refs. 1 and 2). The VSE established the goal of returning humans to the Moon by 

2020 as preparation for a future human Mars mission. As a result of the VSE, NASA established the 

Constellation Program1

. In order to return humans to the Moon, a good deal of research and trade studies 

was conducted in support of returning humans to the Moon. The focus of this paper is on the transfer 

trajectory from Earth to the Moon, i.e., the translunar injection (TLI) trajectory. The TLI problem is a 

trajectory optimization problem; deliver the most mass to the Moon while meeting mission constraints. 

This is not a new trajectory problem—it was solved during the Apollo program using indirect 

optimization methods. 

The research presented in this paper leads to a solution of the TLI problem; however, a direct 

optimization method is employed to solve the problem, transcribing it to a nonlinear programming (NLP) 

problem that is solved with a commercial solver, SNOPT (Ref. 3). The robustness and flexibility offered 

by the direct method serve as the primary motivations of this research work as it is amenable to different 

types of trajectory optimizations and mission constraint requirements. Furthermore, the direct method 

 

1

 The Constellation Program effectively ended during the time of this research. 
NASA/TM—2012-217215 2

applied to the problem is implemented using collocation techniques, similar to those described in 

Reference 4, in generalized settings in a 3 degree-of-freedom (DOF)/6-DOF vehicle simulation program 

named Mission Analysis Simulation Tool in Fortran (MASTIF) (Refs. 4 and 5). MASTIF development 

began in 2006 and was successfully used to independently verify certain Ares I trajectories (Ref. 6). 

Typically, vehicle simulation and trajectory optimization (hereafter referred to as just optimization) 

software programs are separate programs. Combining optimization and 3-DOF/6-DOF vehicle simulation 

software into the same tool makes the software tool rare. The authors’ are aware of just one other 

software program with similar capabilities, Program to Optimize Simulated Trajectories (POST) (Ref. 7). 

This combination addresses present issues that arise from having separate programs and offers benefits as 

well. Frequently, the outputs of the optimization program are inputs into the vehicle simulation program. 

Also, the required modeling inputs (initial conditions, propulsion, aerodynamics, etc.) can differ between 

different programs. In usual practice, the best an analyst can do is to ensure equivalent modeling between 

different tools, which is a potential source of errors. Trajectories are sometimes optimized in a piece-wise 

fashion because establishing a first guess for a complete trajectory can be difficult. This issue can 

sometimes be remedied when the same tool has dual capabilities such as vehicle simulation and 

optimization (especially when closed-loop guidance is available for use in the software for vehicle 

simulations). For example, the entire trajectory can be simulated to achieve an initial trajectory guess (as 

input) for the (complete) trajectory optimization and not just a certain portion of the trajectory. 

Generally, trajectory optimization tools are limited to 3-DOF modeling of the vehicle. Extending the 

vehicle model to 6-DOF is not readily attainable because of an additional burden of significant software 

development or redevelopment. However, when implemented in the same tool as a 6-DOF vehicle 

simulation, extending the trajectory optimization software to 6-DOF is expected to be readily attainable 

because the built-in 6-DOF infrastructure (e.g., equations of motion) already exists. Another major benefit 

of combining the vehicle simulation and trajectory optimization software in a single software package is 

that code maintenance is required only on one software package instead of two. The list below 

summarizes the motivations for combining vehicle simulation and trajectory optimization software. 

• Reduce modeling errors 

• Save time 

• Easier/better initial guesses 

• Readily extendable to 6-DOF optimization 

• Fewer programs to maintain 

As stated above, the research presented in this paper has two main parts. The first part is the 

description of a novel application of a direct method within the framework of 3/6-DOF vehicle simulation 

software tool and its application to the TLI trajectory optimization problem. The vehicle configuration 

includes an Earth Transfer Stage (ETS), an orbiter, lander, and a bipropellant chemical engine. The 

modeled trajectory initiates from a circular low Earth orbit (LEO) with main engine ignition; the ETS and 

propellant tanks are jettisoned an hour after ignition. The trajectory ends with the spacecraft flying over 

the north pole of the Moon. The second part of this paper highlights the dual capability of the software 

tool by utilizing the results of the TLI optimization (the open-loop solution) in the form certain orbital 

states at main engine cutoff (MECO) as input data into a vehicle simulation (the other “mode” in 

MASTIF) that assesses a potential candidate for the closed-loop guidance algorithm. The closed-loop 

guidance algorithm assessed is a variant of the Space Shuttle's guidance named Powered Explicit 

Guidance (PEG); the fact that it is flight-proven makes it a good candidate (Ref. 8). Finally, the results 

from the open-loop and closed-loop solution are compared to show approximate agreement and efficacy 

of MASTIF’s capabilities. 
NASA/TM—2012-217215 3

Nomenclature 

Isp specific impulse 

Tvac vacuum thrust 

MR mixture ratio 

m0 initial mass 

m ETS ETS mass 

m LOX LOX oxidizer tank mass 

m LH2 LH2 fuel tank mass 

m FPR flight performance reserve propellant tank mass 

m RCS RCS propellant tank mass 

gEarth,SL magnitude of Earth’s gravity vector at sea level (9.80665 m/s) 

GM Moon’s gravitational constant (4902.801076 km3

/sec2) 

Translunar Injection Optimization Problem 

The Earth to Moon transfer problem to be solved consists of the following sequence of events and 

modeling assumptions. Starting in circular 29° LEO, at the epoch of January 1, 2018, (midnight), perform 

a main engine burn to place the spacecraft on the correct 3-day transfer ellipse while maximizing the mass 

at MECO. The start-up and shut-down transients of the engine are ignored, but can be included in the 

future. The atmospheric drag is also ignored for now. Approximately 1 hr after the start of the burn, the 

spacecraft jettisons the ETS. Flight performance reserve (FPR) propellant left over from the ascent to 

orbit and reaction control system (RCS) propellant (carried in two separate propellant tanks) are also 

jettisoned with the ETS. The propellant in these two tanks is not available for use in the 3-DOF 

optimization; it reserved for 6-DOF simulations. Propellant remaining in the tanks from the main engine 

burn is carried along to the Moon, but could easily be jettisoned, as modeling is refined for other future 

needs. 

The spacecraft is to arrive at the Moon flying over the North Pole at a 100 km pericynthion altitude. 

The relatively small ΔV required to capture the spacecraft into a lunar orbit is not currently modeled. 

Table 1 lists the LOX/LH2 bipropellant engine’s specific impulse, Isp, vacuum thrust, Tvac, and mixture 

ratio, MR. The spacecraft’s initial mass, m0, includes the mass of the ETS, mETS, and the mass of the filled 

fuel and oxidizer tanks, mLH2 and mLOX, respectively. Two additional tanks store propellant set aside for 

the reaction control system (RCS), mRCS, and flight performance reserve (FPR), mFPR; no propellant is 

drawn from these tanks for the TLI trajectory optimization. Partitioning the propellant this way makes it 

convenient for future 6-DOF Monte-Carlo analyses, if required. 

TABLE 1.—PROPULSION AND MASS PROPERTIESa

Isp 

(sec) 

Tvac 

(kN) 

MR m0 

(kg) 

mETS 

(kg) 

mLOX 

(kg) 

mLH2

(kg) 

mFPR 

(kg) 

mRCS 

(kg) 

448 1,070.6 5.5 187,901 25,463 81,887 14,889 1,449 342 a

See Nomenclature

Gravitational effects and first order harmonics (J2) from the Earth, Sun, and Moon are modeled. The 

values for the physical constants for these three bodies came from a Jet Propulsion Laboratory (JPL) lunar 

modeling document (Ref. 9). Ephemeris data (Earth, Moon, and Sun time-dependent positions) came 

from the JPL’s SPICE software.2

 

2

 http://naif.jpl.nasa.gov/naif/ 
NASA/TM—2012-217215 4

Discretization of Translunar Injection Problem 

Two vehicles and three (total) phases are utilized to model the problem. The first vehicle includes the 

spacecraft and the ETS and is modeled with two phases. The first phase is the TLI burn and the second 

phase is a nearly 1 hr coast phase. The second vehicle does not include the ETS and the FPR and RCS 

propellant tanks. Its primary role is to easily facilitate the jettison of the ETS and two propellant tanks 

within the software design of MASTIF. The second vehicle has one long coast phase. The first two phases 

are collocation phases (described next) whereas the last phase is an explicitly integrated phase using a 

high order Runge-Kutta integrator (Ref. 10). 

The equations of motion in each collocation phase are discretized by segments that contain nodes. 

Dimensional time is transformed by a Legendre-Gauss-Lobatto (LGL) distribution method that results in 

placing the nodes on a non-dimensional time grid spanning from –1 on the left segment boundary to +1 

on the right segment boundary (Ref. 11). Each segment contains at least three nodes to enable interpolated 

state and control values to be obtained at every other node. There is no maximum limit on the number of 

nodes; however, it must be an odd number due to the mathematical formulation adopted. The number of 

nodes (per segment) corresponds directly to the degree of the interpolating polynomial. The expressions 

for the polynomials are provided in the Appendix; an efficient numerical implementation is the ACM 

algorithm 211 (Refs. 12 and 13). 

Figure 1 shows a schematic of the vehicle and phase layout (with minimal number of nodes and a 

single segment per phase for clarity) and Figure 2 depicts the double nodes that are utilized at all phase, k, 

boundaries to enable discontinuities in the vehicle states, s, and/or controls, u , and will be discussed in 

further detail in the next section. 

Figure 1.—Vehicle and phase schematic. 
NASA/TM—2012-217215 5

Figure 2.—Double nodes at phase boundaries. 

The discretization utilizes a global numbering scheme for the nodes, i, that results in odd-numbered 

and even-numbered nodes starting with the first node (i = 1) contained in the first segment within the first 

phase. The global numbering scheme is an implementation design; all interpolations are local (segment 

level). Potential solutions to the states and controls are provided by SNOPT at each segment’s oddnumbered nodes. These states are used to calculate the true state derivatives at these odd-numbered 

nodes. At each segment’s even-numbered nodes, Hermite polynomial interpolations (using values of 

states and derivatives at the segment’s odd-numbered nodes) are used to estimate the states, which can 

then be used to estimate the state derivatives (Ref. 13). The difference between the estimated state 

derivatives and the true state derivates at each segment’s even-numbered nodes are known as defects, Δ. 

These defects are nonlinear constraints. 

In general, there exist different discretization schemes, expecting the same results within the range of 

accuracy that might vary from user to user. Discretization is advantageous because it enables the user to 

control solution accuracy, which is especially important during the initial attempts at a solution when less 

accuracy is required. Permitting user control over the number of nodes per segment and the number of 

segments per phase is analogous to controlling the order of integration method and integration step-size, 

respectively, in an explicit integration. This numerical control essentially allows accuracy to vary from 

phase to phase according to the modeled dynamics of the problem. In this paper, the discretization process 

is iterative and to accurately represent the trajectory, there must be enough grid points (placement is also 

important) considered for the applied algorithm. 

Problem Formulation 

Now, a general optimal control problem with prescribed dynamical system states, x (t), and a control 

input, u (t), minimizes a cost function 

J = J(x (t),u (t),t), x ∈Rn×1, u ∈Rm×1 (1)

subject to the system’s equations of motion, 

x  ̇ = f(x (t),u (t),t) (2)

and any initial constraints 

χ(x (t0 ),t0 )= 0 , χ ∈R p×1 (3)
NASA/TM—2012-217215 6

terminal constraints 

Ψ(x (tf ),tf ) = 0 , Ψ ∈Rq×1 (4)

and path constraints 

c (x (t),u (t),t)≤ 0 , c ∈Rr×1 (5)

where t0 ≤ t ≤ tf and t0 and tf are initial and final times, respectively. 

In MASTIF, the optimal control problem is transcribed to a constrained NLP optimization problem of 

the form: 

min J(x ) s.t. L ≤

x 

F (x )

G (x )

⎛

⎝

⎜

⎜

⎞

⎠

⎟

⎟

≤ U (6)

and solved by using SNOPT, where L and U are the lower and upper bounds, respectively. Note that in 

this formulation x is an array containing all the independent variables including the controls, u . F (x) is 

an input array into SNOPT containing linear and nonlinear constraints and the objective function, J(x). 

Lastly, G (x) is an (optional) input array into SNOPT containing the first order partial derivatives of 

F (x). Currently, derivatives are not calculated by MASTIF; SNOPT uses finite differencing as an 

estimate of the partial derivatives. The arrays given in Equation (6) are defined below beginning with the 

equations of motion. 

Equations of Motion 

The equations of motion are modeled in EME2000 Cartesian coordinates (denoted with a ‘j2k’ 

subscript). The velocity vector is 

 ̇ r j2k = v j2k . (7)

With thrust being the only (modeled) non-gravitational force exerted on the vehicle, the acceleration term 

becomes 

 ̇ v j2k = Tvac

m

m ˆ j2k + g j2k (8)

where Tvac is the vacuum thrust (atmospheric effects are neglected) given in Table 1, μˆ j2k is the direction 

of thrust (unit magnitude) expressed in EME2000 coordinates. The mass, m, is given by Equation (9) 

below and is a function of a table lookup of total remaining propellant in all the propellant tanks 

m = f( ) table[ ] mLOX + mLH2 + mFPR + mRCS + mETS if phase ≤ 2

f( ) table[ ] mLOX + mLH2 if phase = 3.

⎧

⎨

⎩

 (9)

The last term in Equation (8), g j2k , is the gravitational acceleration due to the Earth, Sun, and Moon:

g j2k = g j2k,Earth + g j2k,Sun + g j2k,Moon . (10)

To complete the equations of motion, the mass flow rates of each propellant tank are 
NASA/TM—2012-217215 7

m  ̇ tnk =

m  ̇ LOX

m  ̇ LH2

m  ̇ FPR

m  ̇ RCS

⎡

⎣

⎢

⎢

⎢

⎢

⎤

⎦

⎥

⎥

⎥

⎥

 if phase ≤ 2

m  ̇ LOX

m  ̇ LH2

⎡

⎣

⎢ ⎤

⎦

⎥ if phase = 3

⎧

⎨

⎪

⎪

⎪

⎩

⎪

⎪

⎪

 (11) 

where the LH2 and LOX tank mass flow rates are calculated as 

m  ̇ LH2 = −m  ̇ eng

MR+1

 (12)

and 

m  ̇ LOX = −m  ̇ engMR

MR+1 (13)

where 

m  ̇ eng = Tvac

IspgEarth,SL

. (14)

No propellant is drawn from the RCS and FPR propellant tanks; therefore, the mass flow rates for these 

tanks are zero. 

Controls, Objective Function, and Constraints 

 To complete the transcription process, the controls, u , the objective function, J, and the constraints 

are defined in this section. The first phase is the only phase with active controls and a Local-Vertical 

Local-Horizontal (LVLH) Euler set is selected for steering during the burn arc 

u =

φ

θ

ψ

⎡

⎣

⎢

⎢

⎢

⎤

⎦

⎥

⎥

⎥

(15)

and is loosely referred to as roll, pitch, and yaw set, respectively. 

To meet the mission objective of maximizing mass, the objective function is the minimum engine burn 

time, tbeng , which is equivalent to maximizing cutoff mass. It is also the same as minimizing the duration 

of phase 1 (Δt1): 

J = tf1 − t01 = tbeng (16)

where t01 and tf1 are the initial and final time for phase 1, respectively. 

The initial LEO is constrained in a way that ensures the orbit is nearly circular, has the correct size and 

is correctly tilted. The orientation of the orbit plane is prescribed by the ascending node angle, which is 

unconstrained. Optimal timing is achieved by allowing the argument of perigee and true anomaly to be 

unconstrained. The three initial constraints are: 
NASA/TM—2012-217215 8

χ =

a(to1 ) − a

inc(t01 ) − inc

e(t01 ) − e

⎡

⎣

⎢

⎢

⎢

⎤

⎦

⎥

⎥

⎥

 (17)

where the values for the semi-major axis, inclination, and eccentricity are given in Table 2. 

TABLE 2.—LEO INITIAL SEMI-MAJOR AXIS, 

INCLINATION AND ECCENTRICITY 

a 

(km) 

inc 

(degree) 

e 

6,563 29 1.0×10–8

To ensure that the spacecraft arrives on the desired trajectory, three final constraints, relative to the 

Moon, are used in the last phase (phase 3) 

Ψ =

alt p (tf3 ) − alt p

inc(tf3 ) − inc

r  ̇(tf3 ) − r  ̇

⎡

⎣

⎢

⎢

⎢

⎤

⎦

⎥

⎥

⎥

 (18)

where the radial velocity component is given by 

r  ̇ = esinh(H) −aGM

a(1− ecosh(H)) . (19)

The radial velocity component constraint ensures that the pericynthion altitude, altp, is attained and is the 

altitude of closed approach. It is equivalent to arriving at an anomaly of zero. The semi-major axis, a, is 

negative because the trajectory is hyperbolic relative to the Moon. H and GM are the hyperbolic anomaly 

and the Moon’s gravitational constant, respectively. 

Table 3 lists the values for the lunar arrival constraints. 

TABLE 3.—DESIRED LUNAR ARRIVAL CONDITIONS 

altp 

(km) 

inc 

(degree) 

r ̇ 

(km/s)

100 90 1.0×10–8

Now, just for organizational convenience, the vehicle states are grouped together in a new state 

array as 

s =

r j2k

v j2k

m tnk

⎡

⎣

⎢

⎢

⎢

⎤

⎦

⎥

⎥

⎥

 (20)

where it is understood that m tnk is the integral of Equation (11) and refers to all four propellant tanks in 

the first and second phase and only the fuel and oxidizer tanks in the third phase. The discretization 

process results in the NLP states becoming discretized as 
NASA/TM—2012-217215 9

x =

s i

u i

t 0k

Δtk

⎡

⎣

⎢

⎢

⎢

⎢

⎤

⎦

⎥

⎥

⎥

⎥

, i =1,3,5,..., k =1,2,33

 (21)

The remaining constraints are linear. These constraints ensure that the start time of one phase is also 

the end time of the previous phase. These constraints are given as 

η = t01 + Δt

1 − t02

t02 + Δt2 − t03

⎡

⎣

⎢ ⎤

⎦

⎥ (22)

where t01 and t02 indicate the start times for phase 1 and 2, respectively. 

Additional linear constraints given below are required to ensure continuity across phases with the 

vehicle states 

ζs = s 2

+ − s 1

-

s 3

+ − s 2

-

⎡

⎣

⎢ ⎤

⎦

⎥ (23)

and the controls 

ζu = u 2

+ − u 1

-

u 3

+ − u 2

-

⎡

⎣

⎢ ⎤

⎦

⎥ (24)

where s 2

+ and s 3+indicate the vehicle state values at the first node of phase 2 and 3, respectively, and s 1

-

and s 2

-

 indicate the values at the last node of phase 1 and 2, respectively (likewise for the controls). Refer 

to Figures 1 and 2 in the previous section for a schematic representation. Because of the mass jettisons 

that occur at the end of phase 2, mass is allowed to be discontinuous between phase 2 and 3. 

Again for convenience, the linear constraints, given by Equations (22) to (24), and nonlinear 

constraints, given by Equations (17), (18) and the defects, Δ, are grouped into a linear constraint array 

ξ =

η

ζs

ζu

⎡

⎣

⎢

⎢

⎢

⎤

⎦

⎥

⎥

⎥

(25)

and a nonlinear constraint array 

δ =

χ

Ψ

Δi+1

⎡

⎣

⎢

⎢

⎢

⎤

⎦

⎥

⎥

⎥

 i =1,3,5,... 4 (26)

where the defects (evaluated at all of the even-numbered nodes) are the differences between the estimated 

state derivatives,  ̇

 ̃

s, and the true state derivatives, s ̇, as described previously, 

 

3

i is global node number and k is phase number 4

i is global node number 
NASA/TM—2012-217215 10

Δ =  ̇  ̃ s − s ̇ . (27)

This permits the dependant variable array, F (x), in Equation (6) to be represented as 

F (x ) =

J

δ

ξ

⎡

⎣

⎢

⎢

⎢

⎤

⎦

⎥

⎥

⎥

. (28)

The TLI problem is now completely transcribed to a constrained NLP problem using Equations (21) and 

(28). 

Closed-Loop Guidance With Powered Explicit Guidance 

PEG is the closed loop guidance algorithm used on the Space Shuttle to handle all phases of its exoatmospheric powered flight. It uses closed-form equations for the propulsive acceleration term (called 

thrust integrals) to accommodate constant thrust or constant acceleration phases of flight. The fueloptimal time histories of the unit thrust vector and its rate of change are derived using Calculus of 

Variations (CoV). PEG can support up to 7 different end conditions. PEG's simplicity and heritage make 

it attractive. Heritage is an important factor in the design selection process because proven designs can 

significantly reduce risk and costs associated with testing. Slightly different versions of PEG appear in the 

literature and have minor differences in the thrust integral terms (Refs. 14 and 15). The version in 

MASTIF is based on Reference 14. 

Numerical Implementation and Results 

The optimization procedure in MASTIF is iterative, primarily for solution accuracy. A procedure 

known as grid or mesh refinement ensures that the quantity and/or placement of nodes represent the 

solution accurately (Ref. 16). An automated grid refinement routine is a future development effort for 

MASTIF. For now, solution accuracy can still be assessed and even improved by simply adding more 

nodes and/or segments. 

Solution accuracy is assessed by taking the optimal control time histories, u *(t), (captured at all the 

odd nodes), and the optimal initial vehicle states, s ∗ (t0 ), (captured at the first node of the first phase), 

from the optimization mode and then explicitly integrating the equations of motion using s ∗ (t0 ) and 

u *(t) from the initial start time, t0 , to the final time, tf . The resulting trajectory is hereafter referred to as 

the explicit solution whereas the solution obtained from the optimization is hereafter referred to as the 

implicit solution because of the two collocation phases that do not require integration of the equations of 

motion. Any differences between the two solutions are errors with the explicitly integrated trajectory 

regarded as the truth model. Because the last phase is a Runge-Kutta phase, the errors are due to the 

implicit solutions in the first two phases. This process of accessing the solution accuracy revealed that the 

lunar arrival constraints had the largest errors. 

TABLE 4.—IMPLICIT VERSUS EXPLICIT ABSOLUTE LUNAR ARRIVAL ERRORS 

nNodesodd altp 

(km) 

inc 

(degree) 

r ̇ 

(km/s) 

21 571.179 7.085 1.157 

35 20.707 1.65 0.183 

61 1.521 0.00016 0.00055 
NASA/TM—2012-217215 11

Table 4 shows how these errors were reduced with additional nodes (nNodesodd). There will be a 

tradeoff between accuracy and execution time. For the purpose of this research, the errors from the 

solution using 61 odd nodes were deemed acceptable and will be used for the remainder of the paper. 

Open-Loop Solution 

The solution to the constrained NLP problem, Equation (6), results in the optimal control, u *, for the 

powered portion of flight. These controls ensure that mass is maximized, while meeting all initial and 

terminal constraints. 

The θ time history is shown as a cubic spline curve fit through the data points at the odd nodes for the 

powered portion of flight in Figure 3. Because nearly all the motion is within the Earth-to-Moon transfer 

plane, out-of-plane motion is minimal (< 0.5°) φ and ψ angles are not shown. 

The mass time history is shown in Figure 4. Mass decreases linearly during the main engine burn 

until MECO. It remains constant during the second phase of flight. One hour after the start of the main 

engine burn, two propellant tanks and the ETS are jettisoned which result in the mass discontinuity. The 

double nodes at phase boundaries as discussed previously enable this discontinuity. The rest of the 

transfer is unpowered; therefore, mass remains constant. 

The spacecraft’s north to south lunar approach trajectory is shown in Figure 5 as well as the Moon’s 

orbital path. 

Figure 3.—Optimal θ time history. 

 	  	  	  	  	



	











 
NASA/TM—2012-217215 12

Figure 4.—Vehicle mass time history. 

Figure 5.—Lunar approach trajectory (600 sat is the name of the spacecraft in visualization tool). 

0 0.2 0.4 0.6 0.8 1 1.2

0.4

0.6

0.8

1

1.2

1.4

1.6

1.8

2

x 105

Time (hour)

Mass (kg)
NASA/TM—2012-217215 13

All constraints were satisfied within acceptable tolerances. The optimized values for the three free 

orbital angles (right ascension of ascending node, Ω, argument of perigee, ω, and true anomaly, ν, are 

shown in Table 5. The burn starts nearly 18° (~ 1 min) before perigee passage (ν = 360°). 

TABLE 5.—OPTIMAL LEO DEPARTURE ANGLES 

Ω

(degree) 

ω

(degree) 

ν

(degree) 

100.1550 226.9605 341.8067 

Table 6 shows that the burn lasts for nearly 400 sec and lists the propellant remaining at MECO. 

TABLE 6.—BURN TIME AND PROPELLANT REMAINING AT MECO 

tbeng 

(sec) 

mLOX

(kg) 

mLH2

(kg) 

395.260595 387 70 

The actual ΔV required, ΔVactual, to perform the burn is compared to the ideal, ΔVideal, in Table 7 

(Ref. 17). The difference is due to the gravity loss associated with having a non-zero flight path angle for 

the finite burn and a zero flight path angle for the instantaneously modeled burn (the ideal velocity 

change). 

TABLE 7.—VELOCITY CHANGE FROM THE TLI BURN 

ΔVideal 

(m/s) 

ΔVactual 

(m/s) 

ΔVgLoss 

(m/s) 

3,142.843 3,157.415 14.572 

Closed-Loop Solution 

The targeted end conditions for PEG are the MECO values of flight path angle, γ, radius magnitude, 

R, and velocity magnitude, V, inclination, inc, and right ascension, Ω, which are obtained from MASTIF’s 

optimization mode (the open-loop solution) and are listed in Table 8. Perfect navigation state knowledge 

is assumed. Guidance is called at a frequency of 1 Hz. Numerical integration of the equations of motion is 

performed with a second order Runga-Kutta with a fixed step-size of 0.01 sec. 

TABLE 8.—PEG TARGET MECO PARAMETERS 

R 

(km)

V

(km/s)

γ

(degree)

inc

(degree)

Ω

(degree)

6,730.5465 10.7965 8.0187 28.9387 100.0423 

The closed-loop LVLH θ time history during the TLI burn is compared to the open-loop solution in 

Figure 6. After five iterations, PEG converges to the same initial θ angle as the optimized open-loop 

solution. At each guidance cycle, the converged solution from the previous cycle is used to start the 

iteration. The simulation proceeds in this manner the until the guidance commanded engine shutdown, 

which indicates that the vehicle has arrived at the MECO target. Differences of 1.5° are seen at the end of 

the burn. 
NASA/TM—2012-217215 14

Figure 6.—θ angle comparison. 

An open-loop and closed-loop performance comparison, as characterized by the burn time required 

(and the resulting remaining propellant), is shown in Table 9. The open-loop burn time is rounded to two 

decimal places because the finest resolution achievable for the closed-loop burn time is 100 Hz due to the 

selected integration frequency. The performance of PEG compares favorably. 

TABLE 9.—OPEN-LOOP AND CLOSED-LOOP MECO TIME AND LOX/LH2 REMAINING 

Method tbeng 

(sec) 

mLOX

(kg) 

mLH2 

(kg) 

Open-loop 395.26 387 70 

Closed-loop 395.29 381 69 

In order for the closed-loop solution to have the correct lunar arrival conditions, the errors at MECO 

must be zero, which they are not, as shown in Table 10. Errors are present because closed-loop systems 

use approximate formulations to derive closed-form solutions. Different PEG input parameters can be 

adjusted to minimize the errors at MECO, but result in slightly different θ angle profiles. No adjustments 

could be made to null all the MECO errors. 

TABLE 10.—CLOSED-LOOP MECO ERRORS 

R 

(m)

V

(m/s)

γ

(degree)

inc 

(degree)

Ω 

(degree)

1,318 1.32 0.044 0.00015 10.00044 

 
  
  
  
 	 
















 
NASA/TM—2012-217215 15

Table 11 shows that MECO errors translate to lunar arrival errors (for the closed-loop system). The 

error in pericynthion altitude is 225 km and the inclination error is 6°. To ensure that the spacecraft 

arrives at the desired lunar location, mid course corrections (MCC) would be needed during the transfer 

and/or MECO target parameters could be updated (instead of fixed) during the burn. Both a targeting 

scheme and MCC modeling are potential future tasks. 

TABLE 11.—LUNAR ARRIVAL CONDITIONS 

Parameter Open-loop Closed-loop 

altp (km) 100 325 

inc (degree) 90 96 

r ̇ (km/s) 0 –0.4 

Conclusion 

Trajectory optimization software utilizing a direct optimization method using collocation was 

developed and implemented as a new mode in a low-to-high fidelity 3-DOF/6-DOF vehicle simulation 

tool named MASTIF. This newly developed optimization mode is used to solve a human-class TLI 

trajectory problem by transcribing the trajectory optimization problem to a constrained NLP problem. So 

far as known from the available literature, this paper presents new research work using a direct method to 

solve this human-class translunar trajectory problem. 

MASTIF’s software capability is rare because it can be utilized as a trajectory optimization program 

and a vehicle simulation tool to design and test GN&C algorithms. MASTIF’s 6-DOF architecture can 

readily be leveraged to extend the optimization mode to 6-DOF. Currently, the authors’ are aware of only 

one program (POST) used within NASA with that functionality. To demonstrate MASTIF's versatility, a 

closed-loop guidance scheme with heritage, PEG, is assessed using the results, in the form of fixed 

MECO target parameters, from the optimization (open-loop) solution. Results show the efficacy of the 

discretization scheme and, additionally, that using fixed MECO target parameters would be acceptable for 

the first design cycle analysis, but that one or more MCCs would be required to ensure the spacecraft's 

desired lunar arrival conditions are satisfied. 

The synergy from combining trajectory optimization and 3/6-DOF vehicle simulation software into a 

single software tool is expected to be a more efficient means to designing, testing, and evaluating vehicle 

GN&C software because model fidelity remains consistent for the trajectory optimization and vehicle 

simulation. 

NASA/TM—2012-217215 17

Appendix: Hermite Interpolating Polynomials 

The expressions below are nearly exactly as those given in Reference 12. The interpolated state value 

is given by 

s  ̃ (t) = αi

i=1

n

∑ (t)s(ti) + βi(t)s  ̇(ti)

i=1

n

∑ (A.1) 

where n is a function of the number of nodes per segment, nNodesseg, 

n = nNodesseg +1

2 . (A.2)

The order of the polynomial corresponds to nNodesseg, which is a user input. 

The generic functions, αj(t) and βj(t) are 

β j(t) = t − ti

t j − ti

⎛

⎝

⎜

⎞

⎠

⎟

2

i=1

i≠ j

n

∏

⎡

⎣

⎢

⎢

⎢

⎤

⎦

⎥

⎥

⎥

t − t ( )j (A.3)

α j(t) = t − ti

t j − ti

⎛

⎝

⎜

⎞

⎠

⎟

2

i=1

i≠ j

n

∏

⎡

⎣

⎢

⎢

⎢

⎤

⎦

⎥

⎥

⎥

1+ 2 t ( ) ( ) j − t

1

t j − ti i=1

i≠ j

n

∑ (A.4)

and it is understood that the subscript j is a dummy index. Likewise, the interpolated state derivative value 

is 

 ̇  ̃ s (t) = α ̇ i

i=1

n

∑ (t)s(ti) +  ̇ βi(t)s  ̇ (ti)

i=1

n

∑ (A.5)

where α ̇ (t) and  ̇ β(t) are obtained by the product rule of differentiation. 
NASA/TM—2012-217215 18

References 

1. Columbia Accident Investigation Board, “The Columbia Accident Investigation Board Report,” 

NASA, 2003. 

2. NASA, “The Vision for Space Exploration,” PDF file, 2004. 

3. P.E. Gill, W. Murray and M.A. Saunders, “SNOPT: An SQP Algorithm For Large Scale Constrained 

Optimization,” Society for Industrial and Applied Mathematics, Vol. 47, No. 1, 2005, pp. 99–131. 

4. C.R. Hargraves and S.W. Paris, “Direct Trajectory Optimization Using Nonlinear Programming and 

Collocation,” Journal of Guidance, Control, and Dynamics, Vol. 10, No. 4, 1987, pp. 338–342. 

5. K.E. Witzberger, D.A. Smith, M.C. Martini and T.W. Wright, “MASTIF User’s Manual,” to be 

released, NASA Glenn Research Center. 

6. S.A. Cook and T. Vanhooser, “The Next Giant Leap: NASA’s Ares Launch Vehicles Overview,” 

Aerospace Conference, 2008 IEEE, Big Sky, Montana, March 1–8, 2008. 

7. S.A. Striepe, and et al., “Program To Optimize Simulated Trajectories (POST II),” NASA Langley 

Research Center and Lockheed Martin Corporation, Version 1.1.6.G, January 2004. 

8. R.L. McHenry and et al., “Space Shuttle Ascent Guidance, Navigation, and Control,” Journal of the 

Astronautically Sciences, Vol. XXVII, No. 1, 1979, pp. 1–38. 

9. R.B. Roncoli, “Lunar Constants and Models Document,” JPL D-32296, Jet Propulsion Laboratory, 

Pasadena, CA. 2005. 

10. E. Fehlberg, “Classical fifth-, sixth-, seventh-, and eight-order runge-kutta formulas with step-size 

control,” NASA TR R-287, NASA, 1968. 

11. C. Canuto, A.Q. Hussaini and T.A. Tang, Spectral Methods in Fluid Dynamics, Springer-Verlag, 

1987. 

12. S.W. Paris, J.P. Riehl and W.K. Sjauw, “Enhanced Procedures for Direct Trajectory Optimization 

Using Nonlinear Programming and Implicit Integration,” AIAA/AAS Astrodynamics Specialist 

Conference and Exhibit, Keystone, Colorado, August 21-24, 2006. 

13. G.R. Schubert, “Algorithm 211, Hermite Interpolation,” Communications of the ACM, Vol. 6, No. 10, 

1963, p. 617. 

14. R.F. Jaggers, “An Explicit Solution to the Exoatmospheric Powered Flight Guidance and Trajectory 

Optimization Problem for Rocket Propelled Vehicles,” Guidance and Control Conference, 

Hollywood, Florida, August 8-10, 1977. 

15. R.F. Jaggers, “Asymmetrical Booster Ascent Guidance and Control System Design Study, Volume V, 

Space Shuttle Powered Explicit Guidance,” NASA CR 140191, NASA, 1976. 

16. J.T. Betts, “Practical Methods for Optimal Control Using Nonlinear Programming,” Society for 

Industrial and Applied Mathematics, 2001. 

17. Plain text file. “Exlx_Flt3.0_Option_1_North_to_South_RKF78_MODIFIED.txt.” From NASA 

Johnson Space Center, 2010. 

REPORT DOCUMENTATION PAGE Form Approved 

OMB No. 0704-0188 

The public reporting burden for this collection of information is estimated to average 1 hour per response, including the time for reviewing instructions, searching existing data sources, gathering and maintaining the 

data needed, and completing and reviewing the collection of information. Send comments regarding this burden estimate or any other aspect of this collection of information, including suggestions for reducing this 

burden, to Department of Defense, Washington Headquarters Services, Directorate for Information Operations and Reports (0704-0188), 1215 Jefferson Davis Highway, Suite 1204, Arlington, VA 22202-4302. 

Respondents should be aware that notwithstanding any other provision of law, no person shall be subject to any penalty for failing to comply with a collection of information if it does not display a currently valid OMB 

control number. 

PLEASE DO NOT RETURN YOUR FORM TO THE ABOVE ADDRESS. 

1. REPORT DATE (DD-MM-YYYY)

01-01-2012 

2. REPORT TYPE

Technical Memorandum 

3. DATES COVERED (From - To)

4. TITLE AND SUBTITLE

Direct Method Transcription For a Human-Class Translunar Injection Trajectory 

Optimization 

5a. CONTRACT NUMBER

5b. GRANT NUMBER

5c. PROGRAM ELEMENT NUMBER

6. AUTHOR(S)

Witzberger, Kevin, E.; Zeiler, Tom 

5d. PROJECT NUMBER

5e. TASK NUMBER

5f. WORK UNIT NUMBER

WBS 789822.01.02.03 

7. PERFORMING ORGANIZATION NAME(S) AND ADDRESS(ES)

National Aeronautics and Space Administration 

John H. Glenn Research Center at Lewis Field 

Cleveland, Ohio 44135-3191 

8. PERFORMING ORGANIZATION

 REPORT NUMBER

E-17896 

9. SPONSORING/MONITORING AGENCY NAME(S) AND ADDRESS(ES)

National Aeronautics and Space Administration 

Washington, DC 20546-0001 

10. SPONSORING/MONITOR'S

 ACRONYM(S)

NASA 

11. SPONSORING/MONITORING

 REPORT NUMBER

NASA/TM-2012-217215 

12. DISTRIBUTION/AVAILABILITY STATEMENT

Unclassified-Unlimited 

Subject Category: 13 

Available electronically at http://www.sti.nasa.gov 

This publication is available from the NASA Center for AeroSpace Information, 443-757-5802

13. SUPPLEMENTARY NOTES

14. ABSTRACT

This paper presents a new trajectory optimization software package developed in the framework of a low-to-high fidelity 3 degrees-offreedom (DOF)/6-DOF vehicle simulation program named Mission Analysis Simulation Tool in Fortran (MASTIF) and its application to a 

translunar trajectory optimization problem. The functionality of the developed optimization package is implemented as a new “mode” in 

generalized settings to make it applicable for a general trajectory optimization problem. In doing so, a direct optimization method using 

collocation is employed for solving the problem. Trajectory optimization problems in MASTIF are transcribed to a constrained nonlinear 

programming (NLP) problem and solved with SNOPT, a commercially available NLP solver. A detailed description of the optimization 

software developed is provided as well as the transcription specifics for the translunar injection (TLI) problem. The analysis includes a 3-

DOF trajectory TLI optimization and a 3-DOF vehicle TLI simulation using closed-loop guidance.

15. SUBJECT TERMS

Trajectory optimization; Trajectory analysis 

16. SECURITY CLASSIFICATION OF: 17. LIMITATION OF

 ABSTRACT

UU 

18. NUMBER

 OF 

 PAGES

26 

19a. NAME OF RESPONSIBLE PERSON

STI Help Desk (email:help@sti.nasa.gov) 

a. REPORT

U 

b. ABSTRACT

U 

c. THIS 

PAGE

U 

19b. TELEPHONE NUMBER (include area code)

443-757-5802 

Standard Form 298 (Rev. 8-98)

Prescribed by ANSI Std. Z39-18