• No results found

A Heuristic Search Algorithm for Asteroid Tour Missions

N/A
N/A
Protected

Academic year: 2021

Share "A Heuristic Search Algorithm for Asteroid Tour Missions"

Copied!
118
0
0

Loading.... (view fulltext now)

Full text

(1)

A Heuristic Search Algorithm for Asteroid Tour Missions

Mohd Bilal

Space Engineering, master's level (120 credits) 2018

Luleå University of Technology

Department of Computer Science, Electrical and Space Engineering

(2)

Cranfield University

Mohd Bilal

A Heuristic Search Algorithm for Asteroid Tour Missions

School of Aerospace, Transport and Manufacturing Astronautics and Space Engineering

MSc

Academic Year: 2017-2018

Supervisor: Joan Pau S´ anchez Cuartielles

August 2018

(3)

Cranfield University

School of Aerospace, Transport and Manufacturing Astronautics and Space Engineering

MSc

Academic Year: 2017-2018

Mohd Bilal

A Heuristic Search Algorithm for Asteroid Tour Missions

Supervisor: Joan Pau S´ anchez Cuartielles August 2018

This thesis is submitted in partial fulfilment of the requirements for the degree of MSc Astronautics and Space

Engineering.

©

Cranfield University 2018. All rights reserved. No part of this publication may be reproduced without the written

permission of the copyright owner.

(4)

Disclaimer

This project has been funded with support from the European Commission.

This publication reflects the views only of the author, and the Commission cannot be held responsible for any use which may be made of the information contained therein.

(5)

Abstract

Since the discovery of Ceres, asteroids have been of immense scientific inter- est and intrigue. They hold answers to many of the fundamental questions about the formation and evolution of the Solar System. Therefore, a mission surveying the asteroid belt with close encounter of carefully chosen asteroids would be of immense scientific benefit. The trajectory of such an asteroid tour mission needs to be designed such that asteroids of a wide range of compositions and sizes are encountered; all with an extremely limited ∆V budget.

This thesis presents a novel heuristic algorithm to optimize trajectories for an asteroid tour mission with close range flybys (≤ 1000 km). The core search algorithm efficiently decouples combinatorial (i.e. choosing the aster- oids to flyby)and continuous optimization (i.e. optimizing critical maneuvers and events) of what is essentially a mixed integer programming problem.

Additionally, different methods to generate a healthy initial population for the combinatorial optimization are presented.

The algorithm is used to generate a set of 1800 feasible trajectories within a 2029+ launch frame. A statistical analysis of these set of trajectories is performed and important metrics for the search are set based on the statistics.

Trajectories allowing flybys to prominent families of asteroids like Flora and Nysa with ∆V as low as 4.99 km/s are obtained.

Two modified implementations of the algorithm are presented. In a first iteration, a large sample of trajectories is generated with a limited number of encounters to the most scientifically interesting targets. While, a poste- riori, trajectories are filled in with as many small targets as possible. This is achieved in two different ways, namely single step extension and multiple step extension. The former fills in the trajectories with small targets in one step, while the latter optimizes the trajectory by filling in with one aster- oid per step. The thesis also presents detection of asteroids for successfully performing flybys. A photometric filter is developed which prunes out badly- illuminated asteroids. The best trajectory is found to perform well against this filter such that nine out of the ten planned flybys are feasible.

(6)

Contents

Disclaimer i

Abstract ii

Acknowledgements vi

List of Figures viii

List of Tables x

List of Algorithms xi

1 Introduction 1

1.1 Background . . . 1

1.2 Asteroids . . . 1

1.2.1 Asteroids and Origin of the Solar System . . . 2

1.2.2 Asteroid Classification . . . 3

1.2.3 Asteroid SMASSII Taxonomy . . . 4

1.3 Asteroid Tour Mission . . . 5

1.4 Aim of the Thesis . . . 8

1.4.1 Objectives of Thesis . . . 8

1.5 Interplanetary Trajectory Design . . . 9

1.5.1 Two Body Orbital Boundary Value Problem . . . 9

1.5.2 Porkchop Plot . . . 11

1.5.3 Excess Velocity and Declination . . . 12

1.5.4 Patched Conics and Planetary Flyby . . . 13

1.5.5 Deep Space Maneuver . . . 15

1.5.6 Minimum Orbital Intersection Distance . . . 15

2 Mathematical Optimization 18 2.1 Mixed Integer Non-linear Programming . . . 19

2.2 Global Optimization . . . 20

(7)

Contents iv

2.2.1 Deterministic Algorithms . . . 20

2.2.2 Metaheuristic Algorithms . . . 22

2.3 Genetic Algorithms for Interplanetary Trajectory Design . . . 26

2.3.1 Analysis of GA Parameters . . . 30

3 Mixed Integer Programming Solver 42 3.1 CASTPath MIP Solver . . . 43

3.1.1 Integer Genetic Algorithm Solver . . . 46

3.1.2 Continuous Optimization Solver . . . 47

3.1.3 Solution Scoring . . . 47

3.2 Extended MIP Solver . . . 49

3.2.1 Single Step Extension . . . 50

3.2.2 Multiple Step Extension . . . 50

4 Case Study and Results 54 4.1 Description of Simulation Cases . . . 54

4.1.1 Case Study 1 : CASTPath MIP Solver . . . 54

4.1.2 Case Study 2 : CASTPath SSE Solver . . . 55

4.1.3 Case Study 3 : CASTPath MSE Solver . . . 55

4.2 Simulation Setup . . . 55

4.2.1 Asteroid Belt Data Distribution . . . 56

4.2.2 Simulation Parameters . . . 57

4.2.3 Launcher . . . 59

4.2.4 Launch date . . . 60

4.3 Results and Discussion . . . 63

4.3.1 Case Study 1 : CASTPath MIP Solver . . . 63

4.3.2 Case Study 2 : CASTPath SSE Solver . . . 72

4.3.3 Case Study 3 : CASTPath MSE Solver . . . 74

5 Asteroid Detection and Optical Navigation in the AMB 78 5.1 CASTAWAY Optical Instrument . . . 78

5.2 Asteroid Photometry . . . 79

5.2.1 Phase Angle . . . 79

5.2.2 Absolute Magnitude . . . 79

5.2.3 Opposition Effect . . . 80

5.2.4 Photometric Model of Asteroids . . . 81

5.2.5 Lumell - Bowell Photometric Model . . . 81

5.3 Photometric Filtering in MIP Solver . . . 83

5.4 Optical Navigation . . . 85

(8)

Contents v

6 Performance Improvement of CASTPath MIP Solver 89

6.1 Traveling Salesman Problem . . . 90

6.2 Heuristic Initial Population Generation Algorithms . . . 90

6.2.1 Randomized Nearest Neighbour . . . 91

6.2.2 k-Nearest Neighbours . . . 93

6.2.3 Opposition Based Learning . . . 94

6.3 Performance Analysis . . . 94

6.3.1 Randomized Nearest Neighbour (RNN) . . . 94

6.3.2 Opposition Based Learning . . . 95

7 Summary 99 8 Future Work 101 9 Conclusions 102 Bibliography 103 Appendix A 103 A.1 Derivation of Detection Range of Asteroids . . . 103

A.2 List of Prominent Asteroid Families . . . 104

(9)

Acknowledgements

I am grateful to my supervisor Dr Joan Pau S´anchez Cuartielles for provid- ing me the opportunity and support to pursue the research of my interest.

His supervision was exemplary and full of novel insights and constructive critique. I had many interesting discussions with Dr S´anchez on the topic of research which led to new ideas and my growth in the field of astrodynamics.

I would also like to thank him for his personal advice during the academic year.

I would also like to thank Dr Jennifer Kingston and Dr Victoria Barabash for their advice and support throughout the duration of the thesis. I also acknowledge the IT service and academic support teams at Cranfield Uni- versity who provided High Performance Computing resources and literature resources respectively. I am also grateful to the EACEA (Education, Au- diovisual and Culture Executive Agency) and the Erasmus+ Programme of the European Union for the Erasmus+ scholarship, without which the thesis would not have been possible.

Finally, I would like to thank my family and fellow SpaceMaster student Damiana Catanoso who have been a source of emotional support.

(10)

vii

(11)

List of Figures

1.1 Prominent asteroid classes . . . 4

1.2 Asteroid tour mission trajectory . . . 7

1.3 Geometry of Lambert’s problem . . . 10

1.4 Lambert arcs : Earth to Mars . . . 11

1.5 Porkchop : Earth - Mars transfer . . . 12

1.6 Declination angle, δ . . . 13

1.7 Mars swingby : deceleration . . . 15

1.8 Mars swingby : acceleration . . . 16

1.9 Tisserand graph : Venus, Earth, Mars . . . 16

2.1 Continuous genetic algorithm . . . 25

2.2 Generic two impulse transfer . . . 27

2.3 Effect of population size : trajectory . . . 31

2.4 Effect of population size : trajectory . . . 32

2.5 Effect of population size : trajectory . . . 32

2.6 Effect of population size : ∆V . . . 33

2.7 Effect of population size : ∆V . . . 33

2.8 Effect of population size : ∆V . . . 33

2.9 Effect of bounds : trajectory . . . 35

2.10 Effect of bounds : trajectory . . . 35

2.11 Effect of bounds : trajectory . . . 36

2.12 Effect of bounds : ∆V . . . 36

2.13 Effect of bounds : ∆V . . . 37

2.14 Effect of bounds : ∆V . . . 37

2.15 Effect of mutation rate : trajectory . . . 38

2.16 Effect of mutation rate : trajectory . . . 39

2.17 Effect of mutation rate : trajectory . . . 39

2.18 Effect of mutation rate : ∆V . . . 40

2.19 Effect of mutation rate : ∆V . . . 40

2.20 Effect of mutation rate: ∆V . . . 40

(12)

List of Figures ix

3.1 CASTPath MIP Solver . . . 45

3.2 Extended MIP Solver : SSE . . . 51

3.3 Extended MIP Solver : MSE . . . 53

4.1 Spatial distribution of asteroids in CATABC database . . . 56

4.2 Spatial distribution of asteroids by size . . . 57

4.3 Spatial distribution of small and large asteroids . . . 58

4.4 Porkchop analysis for Asteroid Tour Mission (CASTAWAY) . 62 4.5 Exploration windows for case studies . . . 63

4.6 Launch region selection method . . . 64

4.7 Solution properties : semi-major axis and diameter . . . 65

4.8 Solution properties : semi-major axis and eccentricity . . . 66

4.9 Solution properties : scores and ∆V . . . 67

4.10 Solution properties : scores and ∆V histogram . . . 68

4.11 Solution properties : ∆V of integer and continuous GA solver 68 4.12 Continuous optimization performance : ∆V improvement . . . 69

4.13 Best science score trajectory . . . 71

4.14 Best science score trajectory : spacecraft’s distance from Sun and Earth . . . 72

4.15 Cast study 2 : Single Step Extension of flybys . . . 73

4.16 Case study 3 : MSE stepwise trajectories . . . 76

4.17 Spacecraft’s distance from Sun and Earth . . . 77

5.1 Phase angle : Sun-asteroid-spacecraft . . . 80

5.2 Saturn’s rings in opposition, courtesy NASA [? ] . . . 81

5.3 Relationship between range and phase angle . . . 83

5.4 MIP Solver with photometric filter after MOID . . . 84

5.5 Optical navigation window and asteroid visibility . . . 85

5.6 Asteroid phase angle and optical navigation window . . . 86

5.7 Earth contact angle : definition . . . 87

5.8 Earth contact angle for each flyby . . . 87

5.9 Distance of spacecraft from the Sun and Earth . . . 88

6.1 Nearest neighbour and short-sightedness . . . 91

6.2 Normalized gradient of tmoid . . . 93

6.3 Randomized nearest neighbour solutions . . . 95

6.4 Scores and ∆V : CASTPath population seeding . . . 96

6.5 Scores and ∆V : RNN population seeding . . . 96

6.6 Fitness evaluation by opposition learning . . . 97

6.7 Runtime performance of opposition learning . . . 98

(13)

List of Tables

1.1 Missions to asteroids . . . 6

2.1 GA parameters : effect of population size . . . 31

2.2 GA parameters : effect of bounds . . . 34

3.1 Integer genetic algorithm set up options . . . 46

3.2 Classification of asteroids by semi-major axis . . . 48

3.3 Weights of scoring parameters . . . 49

4.1 Soyuz performance . . . 59

4.2 Ariane 62 performance . . . 60

4.3 ESA medium mission call and launch years . . . 61

4.4 Trajectories with the best science score and least ∆V . . . 70

4.5 SSE Solver performance . . . 74

4.6 MSE Solver performance . . . 75

5.1 OSIRIS Wide Angle Camera (WAC) [? ] . . . 79

9.1 Prominent asteroid families : names and orbital parameters . . 104

(14)

List of Algorithms

1 Branch and Bound . . . 21

2 Simulated Annealing . . . 23

3 Roulette Selection . . . 29

4 Mating and Crossover . . . 30

5 Integer GA Solver . . . 47

6 Continuous GA Solver . . . 48

7 Randomized Nearest Neighbour . . . 92

(15)

Chapter 1 Introduction

1.1 Background

The existence of the Universe has puzzled humans since the beginning of civilization, and perhaps even before. Early philosophers had proposed some remarkable ideas about its origin, from Plato’s Demiurge being responsible for creating order in chaos to Aristotle’s description of the Universe as a fusion of form and matter. Regardless of which philosophical underpinnings are right, the acknowledgement of the Universe’s existence has motivated humans to explore it and gain as much knowledge as possible using scientific methods. Reducing the focus drastically from the Universe to the Solar System, one could ask many scientifically relevant questions like : Is there evidence of origin of life on celestial bodies other than the Earth? How did the Solar System attain its current configuration? In the spirit of scientific reasoning, a more specific question could also be asked, i.e. what can humans learn about the formation of the Solar System by studying the properties of asteroids and how could this be accomplished?

1.2 Asteroids

Asteroids have long been an area of interest to the scientific community, espe- cially astronomers. The first asteroid discovered, Ceres, in 1801 by Giuseppe Piazzi in Sicily, was thought to be a planet. It was reasoned that Ceres could be the missing planet between Mars and Jupiter as predicted by the Titius-Bode law, which provides a formula to approximate the distance of each planet in the Solar System from the Sun as in equation 1.1.

d(AU) = 0.4 + 0.3 × 2m (1.1)

(16)

Asteroids 2

where m = −∞, 0,1,2,3,.. and so on.

It can be seen from the law that at a distance of 2.8 AU there is no planet corresponding to the law’s prediction. Even though the Bode law is not a natural law, scientists of the time hypothesized a missing planet which was then confirmed by the discovery of Ceres. However, the discovery of asteroid Pallas in 1802 led to the possibility that instead of one single missing planet, there might be a number of smaller ones in that region. Consequently, many discoveries were made in the following years, an exercise that continues to date.

1.2.1 Asteroids and Origin of the Solar System

There are many theories and hypotheses on the formation of the Solar Sys- tem, most notably the planetesimal hypothesis by Thomas Chrowder Cham- berlin and Forest Moulton, who proposed the accretion of planetesimals as integral steps in planetary formation [? ]. Planetesimals are solid bodies from 10-1000 km in size, a majority of which were eliminated by collisions or ejection from the Solar System. Most theories in planetary science con- verge on the formation of planetesimals as the most significant stage in the formation of planetary embryos. However, due to still unknown reasons, planetesimals between Mars and Jupiter never underwent accretion to form one single planet, but led to the formation of asteroids. Thereafter, high speed collisions led to the evolution of the asteroid belt with fragmented as- teroids. These collisions deposited a layer of dust on the asteroids, what is known as regolith. The deeper regolith compacted to form breccias, com- monly found in meteorites on Earth [? ]. These processes are a subject of rigorous discussion and must be read with caution as these are not absolute facts. Furthermore, all of the scientific community does not conform to the same opinion. However, such theories and hypotheses serve as important pre- cursors to the role that asteroids play in understanding the deeper workings of the Solar System, and its origin.

In contrast, observable facts about the chemical composition of asteroids lend support to the idea of exploring the asteroid belt. Asteroids share many compositional trends with planets in the Solar System. In general, the density of planets decreases as the distance from the Sun increases, i.e., inner planets are metallic and contain silicates, while outer planets are icy and contain volatile material like methane and ammonia, and the asteroid belt is a transitional region. In the asteroid belt, innermost asteroids are silicaceous (S-type) having a tinge of red, middle asteroids are metallic (M- type) and outer asteroids are carbonaceous (C-type), and also contain water of hydration. Furthermore, the spectra of C-type asteroids resembles that

(17)

Asteroids 3

of carbonaceous chondrite meteorites, and are bluish in color, containing material remains of the primordial Solar System. With this perspective of the role of asteroids in the origin of the Solar System, and the resemblance of their properties to the planets, a survey of the asteroid belt with a provision to visit asteroids of special interest is of high scientific significance. It would be beneficial to uncovering some of the mysteries of the Solar System, and advancing human knowledge. The next section (1.3) deals with defining one such endeavor, an asteroid tour mission, to explore the asteroid belt in a comprehensive manner. This would lead towards defining the aim and objectives of this thesis. However, before delving into the specifics, it is imperative to identify the orbital characteristics and taxonomy of asteroids in the Solar System.

1.2.2 Asteroid Classification

Asteroids can be classified in many different ways, based on their composi- tion, spectral properties, and spatial distribution in the Solar System. The latter helps to gain an understanding of the properties of asteroids in terms of their orbital parameters, and thus this classification is briefly described below. Spatial distribution of prominent asteroid classes is depicted in figure 1.1.

1. Main Belt : Majority of asteroids are predominantly concentrated be- tween 2.0 and 3.3 AU from the Sun. This region is known as the main belt. One of the most notable features of the main belt is the presence of Kirkwood gaps at about 2.5, 2.82, 2.96 and 3.26 AU from the Sun, where only a few asteroids are present. These locations correspond to fractional parts of the period of revolution of Jupiter which destabilizes these orbits.

2. Trojans : Asteroids that are present at L4 and L5 points of a planet are termed Trojans. There are Earth, Mars, Neptune and Jupiter trojans with the latter being largest in number.

3. Near Earth Asteroids : Commonly referred to as NEAs, these asteroids are present in orbits that either approach the Earth closely or cross the Earth’s orbit (Earth-crossers). Their orbits are characterized by a perihelion of less than 1.3 AU. The exact number of NEAs are currently unknown and is an active area of research.

4. In addition to these groups, notable groups are Hildas, Thules, Greeks, and Cybeles.

(18)

Asteroids 4

Figure 1.1: Prominent asteroid classes

1.2.3 Asteroid SMASSII Taxonomy

In 2002, Bus and Binzel presented a new feature-based taxonomy based on the Small Main-belt Asteroid Spectroscopic Survey (SMASSII) [? ]. This taxonomy has the three main groups, C, S, and X, with additional types defined below.

1. C Group : consists of several types of carbonaceous asteroids, which is the majority type in the Solar System. The most common asteroids in the C Group belong to the C-type. C-type asteroids have a low albedo, and are mostly concentrated beyond 3 AU from the Sun. Other types in C group are B-type, Cg Ch, Cgh, and Cb type.

2. S Group : is a wider complex of silicaceous asteroids with absorption feature around 1 µm [? ]. Other types apart from the standard S-type in this group are A-type, Q-type, R-type, K-type, and L-type.

3. X Group : was formed after Bowell et al. introduced E- and M-class asteroids [? ]. These exhibit red colored spectra. E-type asteroids have high albedo compared to M-type and a flat featureless spectrum.

These can be found in the Hungaria family. Steins, an E-type asteroid was visited by the Rosetta mission. Other types in X-group include P, Xe, Xc, and Xk types.

4. Other types in the SMASSII taxonomy are D Type, Ld Type, and O Type.

(19)

Asteroid Tour Mission 5

1.3 Asteroid Tour Mission

Having established the scientific relevance of a mission to the asteroid belt, and the distribution of asteroids in the Solar System, there are different types of missions that could be devised. For instance, the Near Earth Asteroid Mission (NEAR) first flew by asteroid Mathilde on its way to Eros which is an S-class asteroid, while Rosetta performed successful flybys of Steins and Lutetia asteroids during its trajectory to comet 67P. The Dawn mission too visited two asteroids, Vesta and Ceres to study the internal structure and composition of the asteroids, which ultimately found evidence of organic material on Ceres. Although, there have been many successful missions as listed in Table 1.1, there has not been a successful attempt to survey the asteroid belt with a spacecraft traversing the asteroid belt [? ]. In fact, there have only been a handful of asteroids that have been visited, and a majority of the missions have been restricted to one or two flybys. Thus, evidently, there is an opportunity to assess the possibility of a mission that surveys the asteroid belt, and performs a series of flybys of a greater number of asteroids, in what is termed as an Asteroid Tour Mission. This thesis deals primarily with trajectory design for an asteroid tour mission focusing on developing a framework for exploring possible trajectories from Earth to the asteroid belt.

Specification of thesis objectives necessitates a discussion on the nature of the problem as below.

The trajectory of an asteroid tour involves many variables. Imagining a hypothetical mission, a spacecraft may depart from the Earth at td, with its first encounter with a near earth asteroid at t1. Let this asteroid be called A1. Thereafter the spacecraft may perform a swingby of Mars at tm, thereby taking advantage of Mars’ gravity to propel further into the asteroid belt as shown in figure 1.2. The spacecraft then performs successive flybys of m − 1 asteroids such that the total number of flybys is finally m. After, all of these flybys, the spacecraft continues to traverse the asteroid belt until the end of mission duration. As is the case with space missions, the spacecraft can only carry a certain amount of propellant, and can only be given a limited ∆V depending on the launcher. Thus, the mission has to be completed within these constraints.

From the above example, we understand that the spacecraft needs to know how many and which asteroids to flyby. There are about 700,000 asteroids for which the orbital elements are known to a certain degree of accuracy.

This means that the trajectory design method would need to identify the right set of asteroids to be visited given the constraints, and in the right order. This set of asteroids combined with Earth and Mars is hereby defined as the integer solution set, Sint, as expressed by equation 1.2. Thus, the first

(20)

Asteroid Tour Mission 6

Mission Year Target Asteroids Outcome

OSIRIS- Rex

2016 NEA Bennu will reach Bennu

in December 2018

Hayabusa 2 2014 Ryugu (1999 JU3) landed on Ryugu (June 2018)

PROCYON 2014 2000 DP107 flyby failed

Dawn 2007 4 Vesta and 1 Ceres reached Vesta on 16 July, Ceres on 6 March, in orbit around Ceres Rosetta 2004 Steins and Lutetia successful flybys

of Steins and Lutetia

Hayabusa 2003 25143 Itokawa (1998 SF36)

successful sam- ple return in June 2010

Deep Space 1

1998 asteroid Braille 26 km flyby in July 1999

NEAR 1996 S-class 433 Eros, 253 Mathilde

orbit insertion around Eros in February 2000 Galileo 1989 951 Gaspra and 243

Ida

successful flyby of both asteroids Table 1.1: Missions to asteroids

(21)

Asteroid Tour Mission 7

Figure 1.2: Asteroid tour mission trajectory

optimization step in trajectory design for an asteroid tour mission involves identifying the best set of asteroids, and their optimal sequential arrangement in the set, Sint.

Sint= [Earth, M ars, A1, A2, ..., Ai]; i = [1, m] (1.2) As can be noted from the above discussion, in order to estimate the optimal sequential arrangement, the time of each flyby is required. This can be obtained using the minimum orbital intersection distance (MOID) of the trajectory of the spacecraft and the asteroids. A full treatment of this optimization procedure is described in Chapter 3. Once the time of the closest approach is determined from the MOID for each flyby, there is a need to optimize the flyby times around these values. This is evidently a continuous time domain optimization. The set of times of flyby, Scsis expressed in equation 1.3.

Scs = [td, tm, t1, t2, ..., ti]; i = [1, m] (1.3) In effect, trajectory design for an asteroid tour mission involves integer op- timization for asteroid selection and continuous optimization for flyby times,

(22)

Aim of the Thesis 8

earth departure and mars swingby time, thus making it a mixed integer pro- gramming (MIP) problem, which is discussed in more detail in Chapter 2.

The full description of how the trajectory is optimized is presented in Chap- ter 3. The present discussion provides the context for formulation of the objectives of this thesis in a specific and quantitative manner as follows.

1.4 Aim of the Thesis

The aim of this thesis is to develop an efficient and robust algorithmic frame- work for finding optimal trajectories of an asteroid tour mission with at least ten asteroid flybys, in accordance with ESA’s medium mission guidelines.

1.4.1 Objectives of Thesis

1. Selecting a global optimization technique for mixed integer program- ming of the asteroid tour mission trajectory design problem.

2. Implementation of an extended mixed integer programming solver that can search trajectories with ten asteroid flybys out of which at least five asteroids are bigger than 10 km in diameter.

3. Defining and implementing a photometric filtering technique for reduc- ing the solution space of the asteroid tour mission problem, and assess- ing its performance. The filter must be capable of rejecting asteroids that are invisible at an early stage of the algorithm.

4. Obtaining the best trajectories and analyzing them for navigation, and communication feasibility.

5. Analyze the statistical nature of the solutions and use it to set the parameters of the algorithmic framework.

6. Investigate computational improvement of the algorithmic framework.

In order to achieve the aforementioned aim and objectives, a thorough assessment of interplanetary astrodynamics is presented in the following section. Thereafter, Chapter 2, explores a variety of optimization tech- niques, and develops the rationale for selecting and implementing a par- ticular method for the algorithmic framework. The algorithmic framework, later defined as CASTPath MIP Solver, is presented in Chapter 3.

(23)

Interplanetary Trajectory Design 9

1.5 Interplanetary Trajectory Design

An asteroid tour mission is, in effect, an interplanetary mission with the spacecraft in a heliocentric orbit performing several maneuvers, i.e. swing- bys and flybys. Therefore, standard techniques in interplanetary trajectory design apply to the problem. Departing from the Earth, each transfer can essentially be thought of as a two body orbital boundary value problem, such that the Earth departure date can be determined using porkchop plot. How- ever, there is the added complexity of swingby at Mars and flybys which can be treated with the patched conics methodology. These concepts are visited in an explanatory manner in the following subsections, providing insights with an example transfer from Earth to Mars.

1.5.1 Two Body Orbital Boundary Value Problem

The two body orbital boundary value problem is one of the fundamental problems in astrodynamics. It has been and continues to be a subject of analytic, numerical and theoretical research. In simpler terms, the problem can be formulated as finding the transfer orbit between two locations given their position vectors and time of flight between the two points. Consider points P1 and P2 in figure 1.3, with radius vector ~r1 and ~r2 respectively from the focus F. The line joining the two points is called the chord, and has a length c. The included angle at focus F in the triangle P1F P2 is denoted by θ, and the terminal velocities at the two points are ~v1 and ~v2.

Many elliptical transfer orbits are possible between the two points, how- ever, it is geometrically easy to prove that there exists a minimum value of semi-major axis of the transfer orbit below which no elliptical transfer can take place from P1 to P2 as shown by Prussing [? ]. This minimum value of the semi-major axis, am can be calculated to be

am = |~r1| + |~r2| + c

4 (1.4)

Since the energy of an orbit can be expressed by equation 1.5, the mini- mum semi-major axis also corresponds to the minimum energy transfer. For instance, equation 1.4 reduces to Hohmann transfer orbit for c = |~r1|+|~2 r2|, which is a minimum energy transfer and is discussed in 2.3.

E = −µ

2a (1.5)

Furthermore, there exists an ellipse with a minimum eccentricity, termed as the fundamental ellipse by Battin [? ]. However, the fundamental ellipse

(24)

Interplanetary Trajectory Design 10

Figure 1.3: Geometry of Lambert’s problem

is not relevant to this discussion, it is one of the many interesting properties of this problem which is also usually known as Lambert’s problem [? ]. One of the first explorers of the problem, Johann Heinrich Lambert (1728-1779) stated that the time required to traverse an elliptical arc between two points depends only on the semi-major axis of the ellipse, the sum of the distances of the two points from the focus and the chord length giving us equation 1.6. The solution to transfer time is given by what is known as Lambert’s equation and can be found in standard textbooks on astrodynamics.

∆t = t2− t1 = F (a, c, r1 + r2) (1.6) The minimum energy transfer as described above plays a major role in determining the transfer orbit between two points given the time duration of the transfer. Mathematically, Lambert’s problem can be formulated as in equation 1.7.

[~v1, ~v2] = f (~r1, ~r2, ∆t) (1.7) Various solutions to the problem have been proposed since the original method proposed by Carl Friedrich Gauss in 1809. The algorithm imple- mented in the CASTPath MIP Solver, discussed in Chapter 3, follows Bat- tin’s solution based on Universal variables mainly due to its robustness [? ].

(25)

Interplanetary Trajectory Design 11

-4 -3 -2 -1 0 1 2 3 4

108 -2

-1.5 -1 -0.5 0 0.5 1 1.5 2

108

Earth Mars Lambert Transfers

Figure 1.4: Lambert arcs : Earth to Mars

The mathematical formulation can be found in the original work of Battin [? ], and so the present section is restricted to its implementation for an Earth-Mars transfer.

Figure 1.4 shows trajectories called Lambert arcs for different arrival dates on Mars from the same Earth departure date. The cyan colored inner circle represents Earth’s orbit, and the outer red orbit is Mars’. Notice how the eccentricity of the orbit increases as the arrival position moves away from 180° transfer. In case of the asteroid tour mission trajectory design, Lambert arcs would be used for transfer from one asteroid to the next, the solution to which will then be used to compute velocity changes required to effect the transfers.

1.5.2 Porkchop Plot

Solution to Lambert’s problem provides the departure and arrival velocities for a given set of departure and arrival time, thus enabling expression of total cost of the transfer in terms of change in velocity, or delta-V (∆V ). However, these delta-V values can be quickly plotted for more than one departure and arrival date. Thus, contour plots of the velocities for different combinations of departure and arrival can be easily plotted and these plots are known as porkchop plots for their resemblance to pork chops. Other performance parameters like dry mass and semi-major axis can also be plotted giving a quick visualization of the best launch periods and time of transfer. Consid- ering the Earth-Mars transfer discussed in section 1.5.1, a typical porkchop plot is presented in figure 1.5. It can be noticed that the porkchop repeats

(26)

Interplanetary Trajectory Design 12

Figure 1.5: Porkchop : Earth - Mars transfer

itself in nearly two years due to the relative motion of Earth and Mars. Such a porkchop plot be used to assess which departure and arrival dates are the most efficient ones, i.e. the lowest delta-V requirement. This technique would be employed to select suitable Earth departure date and time of flight for the asteroid tour mission in section 4.2.4.

1.5.3 Excess Velocity and Declination

Once a spacecraft has left the sphere of influence (SOI) of Earth, its trajectory can be described in spherical coordinates with respect to the Earth’s equator in terms of the excess velocity V, the declination angle δ, and the right ascension (α). Figure 1.6 shows that declination is the angle between the equatorial plane of the planet (Earth) and the plane of the orbit of the spacecraft as it travels away from the planet at V. The excess velocity refers to the velocity of the spacecraft with respect to the planet at an infinite distance from the planet, and is the difference between the velocity of the spacecraft and the velocity of the planet. The excess velocity is a measure of the launcher’s performance which is generally expressed as square of the excess velocity and called C3 (equation 1.8). It must be noted that ~Vs/c and V~planet are expressed with respect to the Sun.

V~= ~Vs/c− ~Vplanet C3 =

V~

2 (1.8)

(27)

Interplanetary Trajectory Design 13

Figure 1.6: Declination angle, δ

On the other hand, declination depends on the launch location’s latitude, if there are no declination changing maneuvers. If the launch is eastward the declination angle is equal to the latitude of the launch site. A useful fact to visualize declination angle is that interplanetary trajectories in the ecliptic plane are at a declination of 23° as the ecliptic is tilted by the same angle with respect to earth’s equatorial plane. The specific values used for declination and C3 used for the asteroid tour mission are described in section 4.2.

1.5.4 Patched Conics and Planetary Flyby

An asteroid tour mission with flyby of ten asteroids involves a total of 12 celestial bodies, including the Earth and Mars. However, the gravity fields of the asteroids would be unknown and hence, the flybys are merely deep space maneuvers. On the other hand, in the case of Mars swingby a patched conics approach is necessitated.

Patched conics is the method of computing trajectories that involve more than one central body. To put into perspective, ESA’s Rosetta mission was launched in March 2004, and performed a total of three flyby maneuvers around the Earth. A year after its launch, Rosetta encountered Earth again in March 2015 at a height of 300 - 14000 km, and then headed towards Mars.

The spacecraft returned to Earth for two more flybys, first in November 2007, and again in November 2009. It also made two flybys in the asteroid

(28)

Interplanetary Trajectory Design 14

belt observing asteroids Steins and Lutetia [? ]. Considering Rosetta’s flyby of Mars, its trajectory before entering Mars’ sphere of influence (SOI) can be described as elliptical with the Sun at its focus, and hyperbolic as it enters the SOI of Mars with its focal point at Mars. After departing the Mars’ SOI, it can again be considered elliptical. Thus, the total trajectory becomes a combination of conic sections, and is generally an efficient method of exploring the initial and boundary conditions of the mission trajectory [?

]. Even though the patched conics approximation is not precise, making some adjustments to the approximations can provide precise orbit determination.

Flyby or Swingby is an orbital maneuver that uses the gravity of a celestial body to change the velocity of a spacecraft. The change could be positive or negative depending on the relative velocities of the two bodies at the time of flyby. If the magnitude of velocity of the spacecraft is increased the energy of the planet decreases and vice versa. For the purpose of analysis, the following terms are defined, and equation 1.9 defines the relationship between them.

V~ps : velocity of target planet w.r.t Sun

V~ts1 : velocity of spacecraft w.r.t Sun before SOI V~ts2 : velocity of spacecraft w.r.t Sun leaving SOI V~tp1 : velocity of spacecraft w.r.t planet entering SOI V~tp2 : velocity of spacecraft w.r.t planet leaving SOI

∆~Vtp : change in velocity of spacecraft w.r.t Sun

V~ts1= ~Vtp1+ ~Vps V~ts2= ~Vtp2+ ~Vps

∆~Vtp= ~Vts2− ~Vts1= ~Vtp2− ~Vtp1

(1.9)

The change in magnitude of velocity, ∆~Vtp

of the spacecraft depends on the direction of approach to the planet as shown by figure 1.7 and 1.8. A decrease in the velocity magnitude of the spacecraft will lower its heliocen- tric orbit, and the opposite will raise it. Hence, planetary flybys can be used for escaping the Solar System, like Voyager I and II, or slow down for earth re-entry. The terms flyby and swingby are generally interchangeably used, however in this document, swingby is used only if the spacecraft uses the gravity of a planet, for example, Earth or Mars, while flyby is used with reference to the spacecraft flying by an asteroid which is imple- mented as a maneuver in space without taking advantage of the asteroid’s

(29)

Interplanetary Trajectory Design 15

Figure 1.7: Mars swingby : deceleration

gravity. In context of the asteroid tour mission, the spacecraft can use the gravity of Mars to increase its heliocentric velocity, that is its energy, and go further into the asteroid main belt. Tisserand Graphs are an essential tool to assess the usefulness of gravity assist maneuvers or swingbys as shown in figure 1.9. In general, for transfers with a planet in between the target body and launch body (for example, Earth to Mercury with gravity assist at Venus), gravity assist can decrease the departure and arrival C3.

1.5.5 Deep Space Maneuver

Deep space maneuvers (DSM) are called so because they are performed dur- ing the trajectory of the spacecraft without using the gravity of other planets or bodies. DSMs can be useful in between two gravity assists as required by the mission. It must be noted that DSM requires additional storage of pro- pellant. As already mentioned, the flybys of all the asteroids in an asteroid tour mission can be implemented as deep space maneuvers.

1.5.6 Minimum Orbital Intersection Distance

Hereafter referred to as MOID, minimum orbital intersection distance repre- sents the proximity of two confocal elliptical orbits. For two given elliptical orbits, there can be up to four proximities, and at least one proximity is

(30)

Interplanetary Trajectory Design 16

Figure 1.8: Mars swingby : acceleration

0 200 400 600 800 1000 1200 1400

Orbital Period (Earth Days)

0.2 0.4 0.6 0.8 1 1.2 1.4 1.6

Perihelion (AU)

Venus

Earth

Mars

Venus

Earth

Mars 1 km/s

2 km/s 3 km/s 4 km/s 5 km/s 6 km/s 7 km/s 8 km/s 9 km/s 10 km/s V

Figure 1.9: Tisserand graph : Venus, Earth, Mars

(31)

Interplanetary Trajectory Design 17

always guaranteed [? ]. MOID is a useful parameter for approximating any potential collision of an asteroid with another planet, and the MOID of NEAs is generally tracked as they can end up colliding with the Earth . For the asteroid tour trajectory design problem, MOID serves an important role of finding asteroids that are nearest to the spacecraft’s orbit. This is further discussed in section 3.1.

(32)

Chapter 2

Mathematical Optimization

People in all walks of life come across problems that require them to make a decision based on certain constraints. Consider the simple example of a finance engineer going to work from home. There are different modes of transport that the engineer can use each having a different cost and travel time. If the person chooses to take a particular mode, there are generally more than one route to reach the destination. In normal life, these decisions are taken arbitrarily or based on intuition or experience. However, the choice may or may not be the optimal one in that case. On the other hand, in the professional world, people seek more optimal answers to problems. For instance, a ride-sharing company would prefer to know the shortest routes to save fuel costs or find a distribution of driver locations to maximize revenue.

Such solutions are optimal for the given constraints of time, manpower, fuel, etc., and the process of finding these solutions is termed optimization. In mathematical terms, an optimization problem can be defined as

find min

x∈RF (x) (2.1)

subject to Gi(x) ≤ ci, i = 1, 2, .., n (2.2) Equation 2.1 refers to what is typically known as the cost function or the objective function and is formulated based on the problem under consider- ation. In case of the example of ride sharing, the objective function could be a weighted function of fuel, time of travel, distance etc. Similarly, in optimal control theory, the cost function of an LQR problem, described by equation 2.3, is minimized subject to a feedback controller action. For the purpose of maximizing a certain function, its negative is minimized and thus the definition remains consistent.

(33)

Mixed Integer Non-linear Programming 19

J =

Z

0

[xT(t)Qx(t) + uT(t)Ru(t)] dt (2.3) Equation 2.2 defines the constraints on the problem, both equality and inequality constraints. As in the example of optimal control, the feedback control action, u(t) = −Kx(t) is an equality constraint, and restricts the solution space of the problem to only those values of control input and state variables that satisfy the above condition. Furthermore, there can be inequal- ity constraints on the tracking error. These constraints will finally decide the optimal solution to the problem.

Returning to the optimization of trajectories for an asteroid tour mission, a pertinent question is how this general framework of optimization can be utilized and cast into the earlier mentioned combinatorial and continuous domain optimization. The next chapter deals with the explanation of this implementation, while in this chapter, methods of global optimization are evaluated keeping in mind that an asteroid tour mission trajectory design requires global optimization owing to the sheer large number of possible solutions.

2.1 Mixed Integer Non-linear Programming

As we have already observed for the asteroid tour mission, some problems may involve certain decision variables being constrained to integer numbers while others are continuous real numbers. Furthermore, if the cost function in such a problem is non-linear, it is termed as a Mixed Integer Non-linear Programming problem (MINLP). In a generic mathematical sense, MINLP can be described as below [? ].

find min

x,y F (x, y) (2.4)

subject to Gi(x) ≤ ci, i = 1, 2, .., n (2.5)

x ∈ R (2.6)

y ∈ I (2.7)

In such a formulation, the complexity of the problem increases as the number of integer variables increase. For instance, in a ride sharing applica- tion, the number of ride sharing customers and the route length may decide the revenue generation. However, the number of customers can only be inte- ger values while route length can be continuous. Any solver that implements

(34)

Global Optimization 20

revenue maximization for the above car sharing decision variables will need to consider the distinction. In some applications, there may be more than one integer decision variables and a change in their combination could also change the value of the objective function. This adds a combinatorial com- plexity to the problem, and hence may require large computational resources to find the optimal solution if there are non-linearities and a large num- ber of possible combinations. Thus, mixed integer non-linear programming problems are NP-complete from the point of view of complexity [? ].

2.2 Global Optimization

The general form of the optimization problem as in equations 2.1 and 2.2 seems rather simple, but can be very challenging to solve even for a small number of variables. If the cost function or any of the constraints are non- linear, the problem is said to be nonlinear optimization or non-linear pro- gramming problem [? ]. There are various methods for nonlinear optimiza- tion, including global optimization which is the focal point of the asteroid tour trajectory optimization problem. As the name suggests, global opti- mization techniques aim to find the global or the overall minimum of the objective function given the constraints. There are several deterministic and heuristic methods for optimization, and some of them are discussed briefly in the following sections. A detailed discussion on genetic algorithms and evolutionary computing follows in sub-section 2.2.2. The major challenge in developing a global optimization solver is the ability to find the optimal solu- tion within a reasonable amount of time, the definition of which may depend on the problem. Deterministic methods guarantee a solution in limited time, but are prone to get stuck at local optima. Heuristic methods are capable of searching the solution space although there are concerns around the time taken to arrive at the optimal solution, and so they do not necessarily guar- antee the optimal solution. Regardless, they are widely used and are very effective methods in global optimization.

2.2.1 Deterministic Algorithms

Branch and bound method

Branch and bound algorithm is used to solve combinatorial optimization problems. It relies on the depth search technique and using upper or lower bounds to prune feasible solutions. The idea is to search for solutions in areas of solution space where the optimal solution may lie. Searching from an initial point and using the bounding function to estimate an optimistic solution,

(35)

Global Optimization 21

branches are created in depth, and the objective function is evaluated at each step (node). Also, the current best value is updated using the bounding function. At the step where the objective function is worse than the current best or constraints are violated, the branch is cut off, and newer branches are created. Thus, the search is driven by the bounds which play a major role in avoiding an exhaustive search of the complete solution space as it may be very large and computationally non-feasible. The bounds can be found by relaxation of the constraints on the problem. However, if the search space is very large as in the case of visiting ten asteroids from a possible 700000, branch and bound quickly becomes a non-feasible method in terms of its ability to find the global optimal solution within acceptable computational time. The number of branches to be explored would surely be very large, and even though branch and bound is able to reject many branches at an early stage, it would not suffice the demands of the problem in a reasonable amount of time. The pseudo-code for branch and bound is presented below for completeness.

Algorithm 1 Branch and Bound

1: procedure BranchBound(f, g) . f: objective function; g: bounding function

2: q ← ∅ . queue, q holds children of node

3: estimate heuristic based best value, fb

4: initialize current best, b

5: while q 6= empty do

6: select branching node, n ∈ q

7: remove n from q

8: generate child nodes, ci of n, s.t. i = 1 : m

9: for i=1 to m do

10: estimate optimistic bound gi

11: if gi worse than best value,fb then

12: discard ci

13: else if ci completes the solution then

14: b = ci ; fb = gi,

15: else

16: add ci to q

17: end if

18: end for

19: end while

20: end procedure

(36)

Global Optimization 22

2.2.2 Metaheuristic Algorithms

Metaheuristics are non-deterministic procedures that aid the exploration of the solution or search space to find near-optimal solutions. Many techniques employ metaheuristics and some of these are reviewed below.

Ant Colony Optimization (ACO)

In 1990, Marco Dorigo presented the Ant System calling it a new computa- tional paradigm [? ]. The system mimics the shortest path finding ability of a group of ants. Ants communicate through pheromones, which they leave in their trails. When a randomly searching ant finds food, it may choose to follow a different trail to the colony, and deposit more pheromone on its way back. Other ants find this trail and are able to find food. In the event of an obstruction to the path, the ants tend to follow the shorter path around the object. This is because given equal probability due to random choice in the beginning, more ants will be able to cross the shorter path in the same time, which will deposit more pheromone on the shorter path. It must also be noted that the pheromone evaporates with time and thus the longer path becomes less and less attractive. ACO is attractive for problems reducible to route finding on graphs, and its implementation to the Traveling Salesman Problem is well described by the inventor. The reader is thus directed to the original paper [? ]. The integer optimization in asteroid tour trajec- tory design problem is combinatorial in nature, and resembles in many ways the Traveling Salesman Problem. Thus, Ant Colony Optimization is a suit- able method for selecting which asteroids to visit out of the possible list of asteroids.

Simulated Annealing

Annealing is the process of cooling down a metal from a high temperature in a controlled manner. This technique has been exploited to develop a meta-heuristic optimization method called simulated annealing, and is an attractive method for problems with large solution space. This is because of one important feature of the method which is to accept bad solutions as well to allow search of the solution space and avoid getting stuck in a lo- cal minimum. The algorithm, as described by the pseudo-code (2), starts with defining a temperature variable set to a high value, and an initial guess which could be random if there is no knowledge of the solution landscape.

Then, the temperature is reduced and at each reducing step, a new solution is searched for by a predefined move, for instance, a random radial search near the current solution. If the new guess gives a better value of the objective

(37)

Global Optimization 23

function, it is accepted, and if it is a worse solution, it is accepted based on a probabilistic function. The probability of selection depends on the change in the energy (objective function) of the two solutions and the temperature.

The probability function is defined such that it decreases with decreasing temperature i.e., the algorithm is more acceptable of bad solutions at higher temperature than at lower temperature to expand search scope. Due to its ability to search large solution spaces, simulated annealing is an attractive method for combinatorial optimization of the asteroid tour mission trajec- tory.

Algorithm 2 Simulated Annealing

1: procedure SimAnneal(Tmax, iter) → xbest

2: initial guess, xc= x0 . guess can be random

3: xbest ← xc

4: for i = 1 to iter do

5: calculate current temperature, Tc← T (i, Tmax)

6: compute current energy, Ec← E(xc)

7: create neighbouring solution, xi

8: compute energy for neighbouring solution,Ei ← E(xi)

9: ∆E ← Ei− Ec

10: if ∆E > 0 then . for maximization

11: xbest = xc

12: else if Probability(∆E, Ti, ) > rand(0, 1) then

13: xbest = xc

14: end if

15: end for

16: end procedure

Particle Swarm Optimization

Swarm-based optimization is similar to evolutionary algorithms (sec 2.2.2) from the point of view that both techniques use populations to search for optimal solutions which are evaluated by a fitness function. However, the difference lies in the fact that particle swarm optimization implements con- cepts of self-awareness and social interaction. Each particle is moved towards its previous best and global best at each iteration in a stochastic manner.

First introduced by Kennedy et al., PSO is found to be faster than evo- lutionary algorithms on standard problems and easily avoids local minima or maxima [? ]. Particle swarm optimization is, therefore, a very relevant

(38)

Global Optimization 24

algorithm to the asteroid tour trajectory design problem discussed in this document.

Genetic Algorithms

Genetic algorithms (GA) are an important tool for global optimization. First introduced by John Holland, they are a subset of evolutionary algorithms and based on the idea of natural selection [? ]. The main tenets of the algorithm are :

ˆ Fitness function is the equivalent of cost function as previously de- scribed, and is always specific to the problem.

ˆ Initial population provides the initial search points to the algorithm.

An ideal initial population would be diverse and wide ranging such that it covers the whole solution spectrum. It may also be seeded with a few near optimal solutions to drive convergence.

ˆ Selection is the basis for propagating the population forward as is done by nature through its refined process of natural selection. Only the fittest of individuals survive to the next generation. However, in computer applications, there are different methods that can be applied to perform selection, like random pairing, weighted random pairing, tournament selection, roulette wheel, etc.

ˆ Mating resembles the natural process of combination of genes to pro- duce phenotype, for instance, mitosis and meiosis. Furthermore, in the latter stages of meiosis, crossover takes place between the chromosomes of the mother and father, thus providing a process of bringing about change in the species.

ˆ Mutation can be thought of as random changes in gene characteristics.

A certain portion of the population may undergo mutation, however, the best solutions are allowed to carry forward unaltered in a GA.

This is termed as elitism. Some robust mutation methods are uniform, Gaussian, and adaptive mutation.

The aforementioned processes are combined together to create what is known as a genetic algorithm, the flowchart of which is shown in figure 2.1.

Genetic algorithms have been very popular since their inception, and have been tested to solve a variety of problems. They are reliable and are use- ful for both integer and continuous domain optimization. They have also been extensively used in optimizing interplanetary trajectories over the last

(39)

Global Optimization 25

Dene decision variables xi

Exit Formulate cost function

f(xi)

Set GA parameters: population size, generations (n), mating

and mutation functions

Generate initial population

Compute cost for each chromosome

Mate selection

Mating

Mutation Begin GA

∀j ∈ N : j=[1,n]

Is j = n Or Solution converged?

Yes No

Figure 2.1: Continuous genetic algorithm

(40)

Genetic Algorithms for Interplanetary Trajectory Design 26

two decades, and a short review of their application to trajectory design is provided in the next section 2.3. Thus, it is hereby considered that genetic algorithms provide the most complete methodology to tackle the problem of finding near-optimal trajectories for an asteroid tour mission. In the next section, a standard problem in astrodynamics is defined which is then solved using a continuous genetic algorithm, thereby exploring the various control parameters of the algorithm in the context of interplanetary trajectory de- sign.

2.3 Genetic Algorithms for Interplanetary Tra- jectory Design

Genetic algorithms have been a favored method to optimize interplanetary trajectories for a long time, partly due to the large number of possible so- lutions for accessing two or more bodies in the Solar System. Gage et al.

concluded the superiority of GA over traditional grid search techniques for interplanetary trajectory optimization [? ]. Hartmann et al. presented the application of genetic algorithms to low-thrust trajectory optimization back in 1998 [? ]. These were notable early developments in the field of trajectory optimization using evolutionary algorithms which has grown considerably since then. More recent developments include implementations of hybrid ge- netic algorithms for multiple gravity assist problems similar to the asteroid tour mission. This makes genetic algorithms a suitable method for optimizing the trajectory of an asteroid tour mission. Before presenting the implemen- tation of GA for the asteroid tour mission, a test simulation is developed to explore the various parameters like population size, mutation rate, etc. The simulation concerns with developing a genetic algorithm to optimize two im- pulse transfer between co-planar orbits. This problem was chosen due to the fact that an analytical solution is available that proves Hohmann transfer to be fuel optimal for a two impulse co-planar transfer [? ]. This implies that the effect of various GA parameters like mutation, selection and crossover can be studied in the context of a standard astrodynamics problem.

Problem definition

Referring to figure 2.2, consider an inner circular orbit A of radius, r1, and an outer circular orbit B of radius r2, both centered at the Sun, S. A transfer orbit from orbit A to orbit B is desired such that the total velocity change, delta-V (∆V ) is minimum. Thus, the challenge is to find optimal values for the eccentricity and semi-major axis of the transfer orbit. The spacecraft has

(41)

Genetic Algorithms for Interplanetary Trajectory Design 27

Figure 2.2: Generic two impulse transfer

a velocity, V1 in the inner circular orbit, V2 in the transfer orbit at the time of departing from the inner orbit, V3 at the intersection with the outer orbit and its target velocity in the outer orbit is V4. The mathematical formulation of the problem, and implementation of the genetic algorithm is as follows.

Assumptions

1. The initial orbit is circular and has a radius r1. 2. The target orbit is circular and has a radius r2. 3. The transfer orbit is a conic section.

Definition of terms

1. e : eccentricity of transfer orbit 2. a : semi-major axis of transfer orbit

3. ∆va : velocity change required at the inner orbit to enter the transfer orbit

(42)

Genetic Algorithms for Interplanetary Trajectory Design 28

4. ∆vb : velocity change required at the transfer orbit to enter the target outer orbit

5. h : Angular Momentum

6. µs: Gravitational Parameter of Sun 7. p : semi-latus rectum

Equations

1. ∆va for first impulse

∆va = q

v21+ v22− 2v1v2cos(θ1) (2.8) where,

v2cos(θ1) = h

r1 = pµsa(1 − e2)

r1 (2.9)

v1 =r µs

r1 (2.10)

v2 =r 2µs r1 −µs

a (2.11)

2. ∆vb for second impulse

∆vb = q

v32+ v42− 2v3v4cos(θ2) (2.12) where,

v3cos(θ2) = h

r2 = pµsa(1 − e2)

r2 (2.13)

v4 =r µs

r2 (2.14)

v3 =r 2µs

r2 −µs

a (2.15)

3. Total ∆v

∆v = |∆va| + |∆vb| (2.16)

(43)

Genetic Algorithms for Interplanetary Trajectory Design 29

Fitness function and constraints

The total ∆v = |∆va| + |∆vb| has been defined as the fitness function subject to the following constraints :

ˆ Since the perigee distance of the object has to be smaller than r1 and the apogee distance greater than r2 in order to intersect both the inner as well as the outer orbit, the following equations apply [? ].

rp = p

1 + e = a(1 − e) ≤ r1 (2.17) ra= p

1 − e = a(1 + e) ≥ r2 (2.18) GA methods and parameters

Once, the constraints are decided, the algorithm is set up with the following parameters and methodologies.

ˆ Initial Population is generated with random numbers within the bounds of the decision variables, and subject to the aforementioned constraints. Random numbers are generated until all the chromosomes meet the constraint criteria. This is a simple and effective strategy as the problem is relatively simple to solve, and does not demand any complicated methods for population seeding.

ˆ Parent selection for mating is done using roulette selection, a popular and simple method for random selection. Algorithm 3 describes the routine succinctly.

Algorithm 3 Roulette Selection

1: procedure Roulette(p, n) . p: population; n : number of parents to select → parents

2: initialize mating pool

3: rank chromosomes in order of better objective function value

4: create mating pool proportional to rank of chromosome

5: select n chromosomes randomly from the pool

6: end procedure

ˆ Mutation is performed through replacement of a randomly chosen gene of the chromosome by a randomly generated value within the bounds of the gene and global constraints. Mutation rate defines the portion of the population that will undergo this random change.

(44)

Genetic Algorithms for Interplanetary Trajectory Design 30

Algorithm 4 Mating and Crossover

1: procedure MatCross(parents, r1, r2) → children . r1: inner orbit radius, r2: outer orbit radius

2: randomly select which gene to edit in each parent

3: while constraints not met do

4: define mating factor randomly, β ← rand(0, 1)

5: create new genes by mating, gn1 ← g1− β ∗ (g1− g2)

6: gn2 ← g2+ β ∗ (g1− g2)

7: end while

8: create new parents with g1 and g2

9: define crossover point c ← int(n ∗ rand(0, 1))

10: split new parents at crossover point and exchange genes → children

11: if constraint satisfied then

12: keep children

13: else

14: discard children

15: end if

16: end procedure

2.3.1 Analysis of GA Parameters

Having set up the genetic algorithm, it is generally good practice to tune the parameters to satisfy the convergence criteria of the problem at hand. In this study, the effect of population size, bounds and mutation rate is studied on the convergence of the solution.

Effect of population size

In order to assess the effect of population size on the solution, three different population sizes were tested viz., a very low value of 10 individuals, followed by a population of 100 and a relatively large population of 1000. Other pa- rameters common to all three cases are listed in table 2.1. It can be observed in figure 2.3 and 2.5, that a very small population size does not converge to the solution, while a population of 1000 individuals is good enough to con- verge to the solution. The semi-major axis and eccentricity predictions are found to be correct to two decimal places for initial population of 1000. This trend can be attributed to the fact that if the population size is too small there is not enough diversity and breadth in the population to be able to close in or converge to the optimal solution. This can also lead to the algo- rithm getting stuck in a local minimum. The following figures demonstrate

References

Related documents

In this work, we address the problem of design and devel- opment of an UAV and its manipulation mechanism along with object detection and localization, coverage planning and

These values shown in Table 3.4 utilised equations 3.19, 3.20 and 3.31, to calculate the minimum mean output temperature(X 1 ), the required power per radiator area (Q/A) needed

The implementation consisted of obtaining and adapting the Yamazaki- Pertoft Genetic Algorithm UCTP solver, re-implementing the Renman- Fristedt Tabu Search UCTP solver, and

This value is attained by subtracting the sum of the length of all complete solutions of the tabu search algorithm from the sum of the length of all complete solutions of the

Program directors Karin Fernler and Fredrik Lange   President of SASSE Arian

The aim of this research has been to view how practitioners of contemporary shamanism relate to the sacred, and how the meaning of the sacred is expressing itself in

Slutsatsen som dras är att en variant av discrete differential evolution presterar betydligt bättre än en generisk genetisk algoritm på detta problem, men inget generellt antagande

Representation-based hardness results are interesting for a number of rea- sons, two of which we have already mentioned: they can be used to give formal veri cation to the importance