TSINGHUA SCIENCE AND TECHNOLOGY
ISSNll1007-0214ll02/21llpp7-11
Volume 14, Number S2, December 2009
Trajectory Optimization of the Exploration of Asteroids
Using Swarm Intelligent Algorithms*
ZHU Kaijian (), LI Junfeng (	), BAOYIN Hexi ()**
School of Aerospace, Tsinghua University, Beijing 100084, China
Abstract: As a typical NP-HARD problem, trajectory optimization has become the hot point in aerospace
research. Considering the limitation of classical optimization algorithms, variant global optimization algorithms for trajectory optimization are applied in an extensive literature. Based on the parameter optimization
of low thrust and impulse maneuver, this paper investigates the difference evolution (DE), particle swarm
optimization (PSO), and genetic algorithm (GA) algorithms and generates two hybrid algorithms. DE proves
to be superior to other non-hybrid algorithms in trajectory optimization problems. The hybrid algorithm of
PSO and DE can improve the optimization performance of DE, which is validated by two given benchmarks.
The results in this paper indicate the validity and feasibility of the hybrid algorithms.
Key words: trajectory optimization; low-thrust; swarm intelligent algorithm; exploration
Introduction
With the rapid development of the aerospace technology and the coming of the second boom of deep space
exploration, the spacecraft trajectory design problems
which can be formalized as optimization problems
have been the hot point in aerospace research. Good
solutions of these problems, usually obtained with local optimization methods, are really a good start for
low-thrust trajectory optimization and the combination
of gravity assist and impulse trajectory optimization.
Global optimization for trajectory problems has been
dealt with recently in many papers. The work of Petropoulos et al.[1] identified gravity-assist trajectories to
Jupiter by using patched-conic techniques. Hartmann
et al.[2] as well as Abdelkhalik and Mortari[3] applied
the genetic algorithm to obtain optimal spacecraft trajectories. Cui et al.[4] utilized the concept of the
Received: 2009-05-08; revised: 2009-06-20
* Supported by the National Natural Science Foundation of China
(Nos. 10832004 and 10602027)
** To whom correspondence should be addressed.
E-mail: baoyin@tsinghua.edu.cn; Tel: 86-10-62795211
Lyapunov feedback control law to derive the semianalytical expressions of suboptimal thrust angles, and
got a near-optimal solution which can approach the
optimal solution by using sequential quadratic programming (SQP) to adjust the five weights in the
Lyapunov function. Ulybyshev[5] used the inner-point
algorithm to optimize the Earth elliptic orbit maneuver
through processing the discrete legs. Betts[6] investigated the orbit transfer optimization by utilizing the
SQP method, which reflects the advantages and features of low thrust. Izzo et al.[7] introduced a deterministic search space pruning algorithm in a combination
of some heuristic global solvers and obtained a very
efficient method to solve MGA problems. Vasile and
De Pascale[8] proposed a global search approach which
hybridized an evolutionary-based algorithm with a systematic branching strategy for MGA problems.
The aim of this paper is to compare various realistic
global optimization algorithms, such as the genetic
algorithm, particle swarm optimization, differential
evolution algorithm and two hybrid algorithms. Two
benchmarks for spacecraft trajectory global optimization are proposed with fully bounds and units of the
state variables and a comprehensive description of
Tsinghua Science and Technology, December 2009, 14(S2): 7-11
given constraints. According to the analysis on the
characteristics of the algorithms, the paper shows that
domain knowledge is the same important as heuristic
knowledge to the performance of global optimizer
based on stochastic search.
Problem Descriptions
1.1
Low thrust trajectory optimization
The spacecraft launched from the Earth must rendezvous with the specific asteroid among the selected ones
during its tour. The initial mass of spacecraft is 2000
kg which can be used as that of propulsion fuel. The
hyperbolic velocity excess of the spacecraft relative to
the Earth is 0.5 km/s and the direction can be set at will.
The maximum magnitude of the engine is 0.15 N with
desired direction and the specific impulse is 3000 s.
For the dynamics model of the two-body problem,
we process the low thrust as perturbation on the small
side. In the inertial coordinate, the parameters change
rapidly and the classical orbit elements exhibit singularities for e=0, and i=0, 90 degrees. Kechichian[9] analysed the near Earth orbit transfer of low thrust spacecraft where dynamical model was defined by using a
modified set of equinoctial coordinates to avoid the
singularities in the classical elements. But the integrate
step must be increased in order to maintain the integrate speed and lower the integrate precision. This paper applies the non-dimensional modified equinoctial
elements to avoid singularities during the numerical
integral. The equations of motion are described as[10]
! ! !nft ,
f ! ! sin (L)f r "[! cos (L) " n(cos (L) " f )] f t # nXg f n ,
g ! #! cos(L)f r " [! sin (L) " n(sin (L) " g )] f t " nXff n ,
1
h ! ns 2 cos(L) f n ,
2
1 2
k ! ns sin (L) f n ,
2
1
L ! nXf n " 2 ,
n!
T
(1)
m!#
Ispg
where
!
n!
, X ! h sin (L) # k cos(L), s 2 !1 "
1 " f cos(L) " g sin (L)
h2 " k 2 .
The modified equinoctial elements can be transformed to the Cartesian state r (rx , ry , rz ) and
V (vx , vx , vx ) according to the expressions:
r
$
2
2
%rx ! s 2 [cos(L) " (h # k )cos (L) " 2hk sin (L)],
%
%r ! r [sin (L) # (h 2 # k 2 )sin (L) " 2hk cos (L)],
% y s2
%
%rz ! 2r (h sin (L) # k cos(L)),
%
s2
%
%v ! # 1 [sin (L) " (h 2 # k 2 )sin (L) # 2hk cos (L) "
& x
!s 2
%
g # 2 fhk " (h 2 # k 2 ) g ],
%
%
1
%v y ! # 2 [# cos(L) " (h 2 # k 2 )cos (L) " 2hk sin (L) #
!s
%
%
f " 2 ghk " (h 2 # k 2 ) f ],
%
%
2
%'vz ! !s 2 [h cos (L) " k sin (L) " fh " gk ]
(2)
where
!2
r!
.
1 " f cos(L) " g sin (L)
Taking account of the representation of asteroid orbit element provided by the organizer, we can transform the modified equinoctial elements into the classical elements:
$! ! a(1 # e2 )/ , f ! e cos(( " ) ),
%
& g ! e sin(( " ) ), h ! tan(i / 2)cos () ),
(3)
%
'k ! tan(i / 2)sin () ), L ! ( " ) " *
1.2
Swing by combined with impulse trajectory
optimization
NASA has achieved a multi-purpose exploration of
Saturn and its moon Titan through Cassini spacecraft
in 1997. Cassini mission from the earth to Saturn includes many trajectory legs, utilizes the gravity assist
of celestial body which can increase or decrease the
orbital energy, and changes the direction of a spacecraft in order to reduce fuel consumption. In general,
gravity assist can contribute to accomplish the long
journey. The preliminary design of multiple gravity-assist trajectories is formulated as a global optimization problem when the sequences of flyby planets
are determined. This paper considers the model with
some assumptions as shown below.
(1) Using the Lambert theory based on the solar
square inverse gravity field between the planets.
ZHU Kaijian () et al.Trajectory Optimization of the Exploration of Asteroids 
(2) Considering the time of spacecraft passing
through the influence sphere of the planet as zero.
(3) Propagating the motion of planet analytically.
(4) Considering the spacecraft as mass point.
The paper applies the impulse maneuver in the interplanetary trajectory so as to satisfy the mission requirements. Generally we divide the trajectory into two
parts and compute the impulse maneuver as two Lamberts problems, i.e. the problem from the departure
body to the gravity assist body and the problem from
the gravity assist body to the arrival body as shown in
Fig. 1. V1 is the vector of impulse at the departure
body; V is the vector of impulse at the periapsis of
hyperbola orbit; V2 is the vector of impulse at the
arrival body; ! in is the angle between V"in and the
connecting line of the celestial body and the main body;
and ! out is the angle between V"out and the connecting line of the celestial body and the main body.
Fig. 1
Single impulse at the periapsis
We can get the basic relation of the variants from
Fig. 1.
in
% % ! & 0;
out
1
'
(
,
in & arccos ) %
2 *
rv
) 1 $ p #in *
)
" p *,
+
1
'
(
,
out & arccos ) %
2 *
rv
) 1 $ p #out *
)
" p *,
+
2" p
p
$ v#2 out %
Optimization Algorithms
2.1
Genetic algorithm
Genetic algorithms (GAs)[11] have been extensively
used as search and optimization tools in various problem domains, such as sciences, commerce, and engineering. The primary reasons for their success are their
broad applicability, ease of use, and global perspective.
The concept of a genetic algorithm was first conceived
by John Holland of the University of Michigan, Ann
Arbor. Thereafter, he and his students have contributed
much to the development of the field.
This paper applies the realcoded GAs to the trajectory optimization. We get ten results of every benchmark through GAs where the goodness of an individual in the population is measured by its fitness value.
GAs select individuals to reproduce, perform crossover,
and mutation to make the offspring, evaluate the individual fitnesses of the offspring, and finally replace the
worst ranked part of the population with the offspring.
From Table 1, we can see that the GAs cannot get the
optimal solution of the two benchmarks and are not
panacea. GAs are developed as a framework for global
searches in the designed space and require many function evaluations. They are more than thirty parameters
to be optimized in trajectory optimization, which results in the invalidation of the traditional binary coded
GAs. Thus some novel techniques must be applied in
GAs before they can be used in benchmarks. Matthew
and Kathleen[12] developed a novel technique for
global low-thrust, interplanetary trajectory optimization through a Gas- and a gradient-based direct method
hybridization.
2.2
! & arccos(V#inV#out );
.v &
2" p
rp-
$ v#2 in
v#in & |||V#in ||;
v#out & ||V#out ||;
(4)
Differential evolution (DE)
Differential evolution[13], introduced by Storn and
Price, is a stochastic direct search method based on the
differential of a set of parameter vectors. It resembles
the EA in the structure, but differs from the traditional
EA in generating new candidate solutions and using a
greedy selection scheme. It has shown great prospect
in many numerical benchmark problems and real world
applications. In many test suites, DE algorithm outperforms other methods and has the advantages of the
fast convergence rate and of a low computational consumption of function evaluations.
Tsinghua Science and Technology, December 2009, 14(S2): 7-11
10
Table 1
Algorithm
GAs
DE
PSO
PSO and SA
PSO and DE
Results of the optimization algorithm
Benchmark
Population/Iteration
Optimum
Frequence
Low thrust
Swingby
Low thrust
Swingby
Low thrust
400/3000
1895 kg
2/10
1860 kg
100/5000
5.84 km/s
1/10
6.0-6.2 km/s
400/3000
1971 kg
2/10
1969 kg
100/5000
5.48 km/s
8/10
5.48 km/s
400/3000
1928 kg
2/10
1905 kg
Swingby
100/5000
5.36 km/s
1/10
5.46 km/s
Low thrust
400/3000
1916 kg
3/10
1914 kg
Swingby
100/5000
5.45 km/s
2/10
5.57 km/s
Low thrust
400/3000
1977 kg
2/10
1970 kg
Swingby
100/5000
5.01 km/s
1/10
5.30 km/s
From Table 1, we can see that the DE can get the
sub-optimal solution of the problems. Based on the
greedy search principle, DE always has global convergence for trajectory optimization as well as in the first
IEEE evolution algorithms competition.
2.3
Particle swarm optimization (PSO)
[14]
PSO has undergone many changes since it was introduced in 1995 and has been applied in almost all the
engineering fields. Full of potential and fertile in new
idea and new perspectives, it may be studied more
widely in the future years and more modified PSO algorithms may be utilized.
From Table 1, we can see that the PSO cannot get
the sub-optimal solution of the problems. With the increase of the complexity of solution space, the selection of parameters will be more significant. PSO is
very sensitive to the selection of algorithm parameters.
2.4
Result with a maximum frequency
PSO and simulated annealing (SA)
SA[15] is a stochastic approach for solving combinatorial optimization problems whose basic idea comes
from the annealing process of solids. In this process, a
solid is heated until it melts, and then the temperature
of the solid is slowly decreased (according to an annealing schedule) until the solid reaches the lowest
energy state or the ground state. If the initial temperature is not high enough or if the temperature is decreased rapidly, the solid at the ground state will have
many defects or imperfections. It was Kirkpatrick et
al.[15] who firstly used SA to solve a combinatorial optimization problem.
This paper provides a hybrid algorithm of PSO and
SA that the principle of simulated annealing is applied
in the process of renewing the position and velocity of
the individuals in the population, the results of which
are listed in Table 1. Metropolis criterion is adopted to
select the offspring by the fitness value.
2.5
PSO and DE
The hybrid algorithm applied to solve the trajectory
optimization is described as the following steps.
Step 1 Initializing the particle swarm with the
population 10 times the parameter variables number,
and the iteration number of 5000.
Step 2 Computing the fitness of each particle and
recording the best global and local particle.
Step 3 Updating the particle according to the generation.
Step 4 If the criterion of convergence is fulfilled,
program will output the results. Otherwise go to Step 2
and update.
From Table 1, we can see that the PSODE can get
the optimal solution of the problems. The hybrid algorithm can synthesize each advantage of each algorithm
and is desired to get a better performance. For the two
benchmarks, PSODE obtains the optimum 1977 twice
and 5.01 once in the ten results where the other mean
values are 1970 and 5.30 which are very close to the
optimum.
Conclusions
This paper applies multiply global optimization algorithms for two trajectory optimization problems that
may serve as reference problems to further develop
models and solvers. Test results for these problems
obtained with the five solvers are reported. The paper
shows that the standard global optimization solvers
ZHU Kaijian () et al.Trajectory Optimization of the Exploration of Asteroids 
which fit the simple problem with low dimension are
not enough to find the optimum solution of large scale
problem which has many local traps. According to the
NO FREE LUNCH theory, a single approach can not
fit all the trajectory problems. The results indicate that
the collaborative solvers outstand among similar approaches in finding the optimum solutions.
References
11
[7] Izzo D, Becerra V, Myatt D, et al. Search space pruning
and global optimization of multiple gravity assist spacecraft trajectories. Journal of Global Optimization, 2007, 38:
283-296.
[8] Vasile M, De Pascale P. Preliminary design of multiple
gravity-assist trajectories. Journal of Spacecraft and Rockets, 2006, 42: 794-805.
[9] Kechichian J A. Optimal low-thrust orbit geostationary
earth orbit intermediate acceleration orbit transfer. Journal
[1] Petropoulos A, Longuski J, Bonfiglio E. Trajectories to
of Guidance, Control and Dynamics, 1997, 20(4): 803-811.
jupiter via gravity assists from Venus, Earth, and Mars.
[10] Gao Y. Advances in low-thrust trajectory optimization and
Journal of Spacecraft and Rockets, 2000, 37(6): 776-783.
[2] Hartmann J, Coverstone-Carroll V, Williams S. Generation
flight mechanics [Dissertation]. University of MissouriColumbia, Columbia, MU, USA, 2003.
of optimal spacecraft trajectories via a pareto genetic algo-
[11] Goldberg D E. Genetic Algorithms in Search, Optimization,
rithm. Journal of the Astronautical Sciences, 1998, 46:
and Machine Learning. New York, USA: Addison-Wesley
267-282.
Publishing Company, Inc., 1989.
[3] Abdelkhalik O, Mortari D. N-impulse orbit transfer using
[12] Matthew A V, Kathleen C H. Global low-thrust trajectory
genetic algorithms. Journal of Spacecraft and Rockets,
optimization through hybridization of a genetic algorithm
2007, 44(2): 456-459.
and a direct method. In: Proceedings of AIAA/AAS As-
[4] Cui P Y, Ren Y, Luan E J. Low-thrust, multi-revolution
orbit transfer under the constraint of a switch function
trodynamics Specialist Conference and Exhibit. Hawaii,
USA, 2008: 2008-6614.
without prior information. Transactions of the Japan Soci-
[13] Storn R, Price K. Differential evolution: A simple and effi-
ety for Aeronautical and Space Sciences, 2008, 50:
cient adaptive scheme for global optimization over con-
240-245.
tinuous space. Technique Report International Computer
[5] Ulybyshev Y. Continuous thrust orbit transfer optimization
using large-scale linear programming. Journal of Guidance,
Control and Dynamics, 2007, 30(2): 427-436.
[6] Betts J T. Very low-thrust trajectory optimization using a
direct SQP method. Journal of Computational and Applied
Mathematics, 2000, 120: 27-40.
Science Institute, Berkley, 1995.
[14] Kennedy J, Eberhart R C. Particle swarm optimization. In:
Proceedings of the 4th IEEE International Conference on
Neural Networks. Perth, Australia, 1995: 1942-1948.
[15] Kirkpatrick S, Gelatt C D, Vecchi M P. Optimization by
simulated annealing. Science, 1983, 220(4598): 671-680.