Guidance strategies for the boosted landing of reusable launch vehicles

IN

DEGREE PROJECT MECHANICAL ENGINEERING, SECOND CYCLE, 30 CREDITS

,

STOCKHOLM SWEDEN 2019

Guidance strategies for the boosted landing of reusable launch vehicles

AGATHE CARPENTIER

KTH ROYAL INSTITUTE OF TECHNOLOGY

Page 1/137

Master thesis report: Guidance strategies for the boosted landing of reusable launch vehicles.

First name - NAME – Department Date Signature

Redacted by : Agathe Carpentier – DLA/SDT/SPC 11/10/2019

Approved by : Eric Bourgeois – DLA/SDT/SPC

Authorized by: Benjamin Carpentier – DLA/SDT/SPC

Siège Direction des lanceurs Centre spatial de Toulouse Centre spatial guyanais

Page 1/137

Programme Supervisor

Aerospace Engineering, TAEEM, 120 ECTS Gunnar Tibert (KTH)

Systems Engineering Eric Bourgeois (CNES)

Title

Master thesis report: Guidance strategies for the boosted landing of reusable launch vehicles.

Abstract

This document presents the results of the master thesis conducted from April 2019 to October 2019 under the direction of CNES engineer Eric Bourgeois, as part of the KTH Master of Science in Aerospace Engineering curriculum.

Within the framework of development studies for the Callisto demonstrator, this master thesis aims at studying and developing possible guidance strategies for the boosted landing. Two main approaches are described in this document :

• Adaptive pseudo-spectral interpolation

• Convex optimization

The satisfying results yielded give strong arguments for choosing the latter as part of the Callisto GNC systems and describe possible implementation strategies as well as complementary analyses that could be conducted.

Sammanfattning

Denna rapport presenterar resultaten av ett examensarbete som genomfördes från april till oktober 2019 under ledning av CNES-ingenjören Eric Bourgeois, som en del av en masterexamen i flyg- och rymdteknik från KTH, Kungliga tekniska högskolan.

Inom ramen för utvecklingsstudier för bärraketen Callisto syftar detta arbete att studera och utveckla möjliga reglerstrategier för Callistos landing som kontrolleras med raketer. Två huvudsakliga metoder beskrivs:

• Adaptiv pseudospektral interpolering

• Konvex optimering

Resultaten ger starka argument för att välja den senare av dessa två metoder för Callistos reglersystem och beskriver möjliga implementeringsstrategier samt vilka kompletterande analyser som bör genomföras Keywords

guidance, Callisto, convex optimization, adaptative pseudo-spectral interpolation, GNC, thesis

Page 2/137 Acknowledgments

All work conducted during this thesis has been performed at the Launchers Directorate of CNES, the French Space Agency.

I am extremely grateful to my mentor Eric Bourgeois for his continued and valuable insight during this 6 months experience. His support, patience and expertise have allowed me to carry out this project in a great work environment and to learn more than I could have hoped for about guidance and trajectory generation for reusable space launchers.

I also want to thank my examiner Gunnar Tibert for his help and support in carrying this degree project, especially in a foreign space agency with the confidentiality constraints it came with. As he was the teacher responsible for the very first space dynamics course taught at the beginning of KTH Master’s, I am very happy to present this thesis under his authority a couple of years later.

Revisions history

Edition Revision Date Modification subject

1 0 11/10/2019 First edition

Page 3/137

Page 4/137

Contents

Table of contents 5

1. Acronyms and notations 10

1.1 Acronyms . . . . 10

1.2 Notations . . . . 10

2. Introduction 12 3. Mission presentation 14 3.1 Vehicle specifications . . . . 14

3.2 Introduction to guidance and its role in the GNC system . . . . 15

3.3 General optimisation problem for the boosted phase . . . . 17

3.4 Nominal boosted trajectory . . . . 18

3.5 Overview of already existing guidance algorithms and their driving principle . . . 19

3.5.1 Tangent linear steering law . . . . 19

3.5.2 Predictor-corrector guidance law . . . . 22

4. Adaptive pseudo-spectral interpolation method for the boosted reentry 26 4.1 Driving principles of the methodology and existing results . . . . 26

4.1.1 Mathematical principles . . . . 26

4.1.2 Previous evaluations: DLR study . . . . 35

4.2 Implementation of the methodology . . . . 36

4.2.1 Flight envelope characterization . . . . 36

4.2.2 APSI algorithm structure . . . . 40

4.3 Results and implemented method adjustments . . . . 43

4.3.1 Evaluating error induced by tools and basic operations. . . . . 43

4.3.2 First-hand evaluation for 3-D dispersions . . . . 43

4.3.3 Refinements implemented to reduce simulation errors . . . . 46

4.3.4 6-D dispersions implementation . . . . 54

4.4 Conclusion . . . . 57

5. Convex optimization based guidance for the boosted re-entry 62 5.1 Introduction to convex optimization . . . . 62

5.1.1 Mathematical theory background of convex optimization . . . . 62

5.1.2 A few examples of convex functions and subsets . . . . 63

5.1.3 Convex problems classification and dedicated solvers . . . . 64

5.2 Formulation of the guidance problem . . . . 66

5.3 Convexification of the optimization problem . . . . 69

5.3.1 Convexification of the feasible set . . . . 69

5.3.2 Convexification and linearisation of the dynamics . . . . 70

5.3.3 Reformulation of the optimization problem as a convex problem . . . . . 72

73subsubsection.5.3.4 5.3.5 Implementation of the convex problem to an optimal trajectory generation 74 5.3.6 Preliminary results . . . . 75

5.4 Specific adjustments implemented for a trajectory optimization problem . . . . . 75

Page 5/137

5.4.1 Trust region . . . . 75

5.4.2 Angle of attack . . . . 76

5.4.3 Flight time and ignition criterion . . . . 77

5.5 Final convex formulation of the convex problem . . . . 82

5.6 Guidance system implementation . . . . 85

5.6.1 Adapted guidance algorithm for each flight phase . . . . 85

5.6.2 Shut-down criterion . . . . 86

5.6.3 Guidance call frequency . . . . 89

5.6.4 Dispersions . . . . 90

5.7 First-hand performance assessment . . . . 90

5.7.1 Manual backstop in the time-to-go estimation for cases where the vehicle starts with a low altitude . . . . 91

5.7.2 Shut-down criterion backstops . . . . 92

5.7.3 Presentation of consolidated results . . . . 95

5.8 Consolidated performance assessment, coupled with results from glided phase . . 97

5.8.1 Idle mode phase . . . . 97

5.8.2 New dispersions considered . . . . 99

5.8.3 Presentation of the glided descent trajectories . . . 105

5.8.4 Definition of a new range in accordance with the nominal trajectory from the glided descent phase . . . 108

5.8.5 Results . . . 110

5.8.6 Allowance of an extra-iteration to ensure the guidance convergence . . . . 111

112subsubsection.5.8.7 5.8.8 Redefinition of the norms in the cost function . . . 116

5.8.9 Average fuel consumption analysis . . . 122

5.8.10 Consumption sensitivity analysis . . . 125

126subsubsection.5.8.11 5.9 Conclusion . . . 127

6. Conclusion 132

Page 6/137

List of Figures

2.1 Callisto’s boosted landing . . . . 13

3.1 Overall aspect of Callisto with fins and pods . . . . 14

3.2 Overview of the overall Callisto trajectory . . . . 15

3.3 Logic of a guidance algorithm during landing . . . . 16

3.4 Guidance role in the GNC system . . . . 17

3.5 Overview of the global optimization problem . . . . 18

3.6 Nominal trajectory with way-point . . . . 19

3.7 Nominal trajectory with way-point flight parameters . . . . 20

3.8 Simplified guidance dynamics model: optimal angle to join a circular orbit . . . . 20

3.9 Predictor-Corrector guidance strategy schematics. . . . . 23

4.1 Examples of dispersions in n dimensions . . . . 27

4.2 Reference subspace selection for 3-D dispersions . . . . 28

4.3 Adaptive Pseudo-spectral interpolation principle- step 1 . . . . 29

4.4 Adaptive Pseudo-spectral interpolation principle- step 2 . . . . 30

4.5 Legendre polynomials on [−1, 1] for different n ∈ N . . . . 30

4.6 Adaptive Pseudo-spectral interpolation principle- step 3 . . . . 31

4.7 Example of a multivariate interpolation with a decomposition in 1-D problems . 32 4.8 Pseudo-spectral approximation of the Runge function with evenly distributed collocation nodes . . . . 33

4.9 Pseudo-spectral approximation of the Runge function with collocation nodes cho- sen to mitigate the Runge phenomenon . . . . 34

4.10 Adaptive Pseudo-spectral interpolation principle- step 4 . . . . 35

4.11 Geodesic coordinates . . . . 36

4.12 East-North-Up speed coordinates . . . . 37

4.13 Bank angle . . . . 37

4.14 Simulation loopback test scheme . . . . 42

4.15 Pseudo-spectral and simulation-induced errors . . . . 44

4.16 First-hand evaluation of APSI with 3-D longitudinal dispersions . . . . 45

4.17 Reference subspace with smoothed controls . . . . 46

4.18 Evaluation of APSI with 3-D longitudinal dispersions and smoothed controls . . 47

48figure.caption.35 4.20 Evaluation of APSI with 3-D longitudinal dispersions, smoothed controls and refined reference subspace . . . . 49

4.21 Monte Carlo campaign for longitudinal dispersions . . . . 50

4.22 Monte Carlo campaign for lateral dispersions . . . . 51

4.23 In flight re-interpolation scheme . . . . 52

4.24 Monte Carlo campaign for longitudinal dispersions with re-interpolation at 10s . 53 4.25 Monte Carlo campaign for lateral dispersions with re-interpolation at 10s . . . . 54

4.26 Monte Carlo campaign for longitudinal dispersions with linear interpolation. . . . 55

4.27 Paired controls scheme . . . . 56

4.28 Example of trajectory yielded with paired controls method . . . . 57

4.29 Monte Carlo campaign with classic APSI paired controls . . . . 58

4.30 Monte Carlo campaign with re-interpolated paired controls . . . . 59

4.31 Monte Carlo campaign with linear interpolation paired controls . . . . 59

Page 7/137

4.32 Monte Carlo campaign with 6 dimensions dispersions . . . . 60

5.1 Classification of convex problems . . . . 65

5.2 Launcher state coordinates . . . . 67

5.3 Launcher angles definitions . . . . 68

5.4 Example on constraints defining a non convex set. . . . . 70

71figure.caption.53 5.6 Implementation of the convex optimization problem for a variable x . . . . 74

5.7 Guidance algorithm principle description . . . . 75

5.8 Example of the solution convergence for the nominal trajectory after one iteration 76 5.9 Effect of penalizing the distance of discretised trajectory state to final state on the trajectory: without constraint (left) and with constraint (right) . . . . 77

5.10 Effect of penalizing a non-zero angle of attack on the smoothness of a trajectory 78 79figure.caption.59 81figure.caption.60 81figure.caption.61 83figure.caption.62 5.15 NewFrontier, taking into account the optimal behaviour, in the (v

, h) plane . . . 83

5.16 Example of a trajectory in the (v

, h) plane with a time-to-go and ignition criterion derived from the NewFrontier . . . . 84

5.17 Cost function of the final convex guidance optimization problem . . . . 85

86figure.caption.66 87figure.caption.67 5.20 Logic for estimating equivalent mass flow rate and time-to-go during the whole trajectory . . . . 87

88figure.caption.69 5.22 FinalFrontier, triggering engine shut-down . . . . 89

5.23 Structure of the code triggering the shut-down depending on the scenario . . . . 93

5.24 Examples of loopholes, and the decision taken by the backstop function . . . . . 94

5.25 Monte Carlo trajectories in the (v

, h) plane . . . . 96

5.26 Monte Carlo positions at initial point, waypoint and touchdown . . . . 97

5.27 Zoom on Monte Carlo positions at initial point, waypoint and touchdown: 96 successes . . . . 98

5.28 Zoom on shut-down and touchdown: trajectories in the (v

, h) plane . . . . 99

5.29 Monte Carlo touchdown conditions . . . 100

5.30 Monte Carlo trajectories . . . 101

5.31 Dispersions for the Monte Carlo 4 failures . . . 102

5.32 Scale of dispersions on plots . . . 102

103figure.caption.82 5.34 Example of optimized trajectory with InitialFrontier and engine ramp-up . . . . 103

104figure.caption.84 104figure.caption.85 5.37 Way-point position precision at way-point with (empty) and without (filled) delay on the angle-of-attack dynamics. . . . 105

107figure.caption.88

108figure.caption.89

5.40 Monte Carlo trajectories in the (v

, h) plane . . . 109

Page 8/137

5.41 Monte Carlo trajectories in the (v

, h) plane . . . 110

5.42 Monte Carlo positions at initial point, waypoint and touchdown . . . 111

5.43 Zoom on shut-down and touchdown: trajectories in the (v

, h) plane . . . 112

5.44 Monte Carlo touchdown conditions . . . 113

5.45 Monte Carlo trajectories . . . 114

5.46 Failed case vs success similar case in the (v

, h) plane . . . 115

5.47 Failed case vs success similar case trajectories . . . 115

116figure.caption.100 117figure.caption.101 118figure.caption.102 119figure.caption.103 119figure.caption.104 120figure.caption.105 120figure.caption.106 121figure.caption.107 121figure.caption.108 123figure.caption.109 123figure.caption.110 5.59 Cumulative distribution function of fuel consumption on 1000 Monte Carlo cases 124 5.60 Propellant consumption analysis for a lower nominal consumption . . . 126

5.61 Final Monte Carlo trajectories . . . 127

5.62 Final Monte Carlo trajectories in the (v

, h) plane . . . 128

5.63 Final Monte Carlo positions at initial point, waypoint and touchdown . . . 129

5.64 Final Monte Carlo touchdown conditions . . . 130

5.65 Final Monte Carlo fuel consumption analysis . . . 131

5.66 Effect of wind speed on launcher dynamics . . . 131

Page 9/137

Page 10/137

1. Acronyms and notations

1.1 Acronyms

Acronyms

CNES Centre National d’Etudes Spatiales (French National Space Agency)

DLR Deutsches Zentrum für Luft- und Raumfahrt (German National Space Agency) JAXA Japan Aerospace eXploration Agency

GNC Guidance, Navigation and Control

LOX Liquid Oxygen

LH2 Liquid Hydrogen

APSI Adaptive Pseudo-Spectral Interpolation

LD Low Density

HD High Density

MC Monte Carlo

SOCP Second Order Cone Programming

LP Linear Programming

QP Quadratic Programming

ttg time-to-go

N/A Not Applicable

FE Flight envelope

1.2 Notations

Shared notations

n

Number of states

n

Number of controls

t

Time-node

x Vehicle state

u Controls

V Relative speed

t Time

t

Final time

h Altitude

α Angle of attack

m Vehicle mass

q Mass-flow-rate

g Gravity constant

I

Engine specific impulse

x

Initial state

x

Final state

Page 11/137 APSI specific notations

µ Bank angle

T

Matrix of states and positions at all low-density time-nodes T

Matrix of states and positions at all high-density time-nodes M

Low-density to high-density transition matrix

N

Number of low-density time nodes N

Number of high-density time nodes p Set of parameter space dimensions p

Parameter space dimension (APSI) f

Multivariate interpolation function

SV Supporting values vector for multivariate interpolation τ

Pseudo-spectral time-node

θ Latitude

φ Longitude

ψ Velocity heading angle (azimuth)

γ Flight path angle

Convex optimization specific notations

C Convex subset of R

v

Vertical speed

v

Lateral speed

r Range

q

Realised mass flow rate q

Controlled mass flow rate

M Vehicle mass

T Thrust

L Lift

D Drag

θ Launcher attitude

ν Speed angle with the local horizontal plane A

, B

, C

Discrete dynamics matrixes at time node t

∆t Time step

N

Number of sub time steps in integrating linear dynamics over one time-step ∆t C

Lift coefficient (vehicle referential)

C

Drag coefficient (vehicle referential) S

Reference surface for aerodynamics forces S

Engine reference surface for thrust forces Γ

Weight on element · · · in cost function

τ Time constant on 1

order dynamics on angle of attack and mass flow rate

Page 12/137

2. Introduction

Ever since the Apollo space flight program, re-entry guidance has been an important topic of aerospace research.

Landing a space launcher represents an important technological challenge for guidance algo- rithms and GNC systems in general. Fluctuating flight conditions, important dispersions and unknown parameters as well as strong control and structural constraints coupled with high ac- curacy demands impose stringent requirements on guidance techniques used for the landing of launchers.

This is especially true for the boosted landing phase. While it represents a very short phase of the overall descent trajectory, it is critical with regard to two main aspects:

• It is the last phase before touchdown, and the accuracy requirements on final state are therefore the most important.

• The flight phase only lasts less than a minute, however the atmosphere conditions and vehicle parameters vary drastically. The unsteadiness of forces represents a great challenge for the guidance robustness.

Various sorts of entry guidance strategies have been developed and studied, mostly in the literature, but also in-flight, to try and succeed in guiding the vehicle to land safely. In every case, the guidance model uses approximations and leads to errors. The landing accuracy re- quirements make use of a feedback controller essential to ensure the safe landing of the guided space vehicle. The guidance needs to be called frequently with navigation updates and therefore needs to be particularly fast in assessing an updated trajectory.

CNES, JAXA and DLR are currently developing a launcher demonstrator, Callisto, whose goal would be to perform a ballistic flight before returning to Earth and landing safely on a pad.

Guidance algorithms are currently being studied to determine possible strategies for the return phase.

The goal of the master thesis was therefore to try and elaborate a guidance strategy that would satisfy the demands of Callisto’s boosted landing. An artistic view of Callisto’s boosted landing can be seen on figure 2.1. Two main strategies were studied in order to assess the performance of each and get some insight into which algorithm would be best suited for Callisto’s flight.

The first approach that is presented is a method combining the use of a data-base of pre- computed optimal trajectories and multivariate interpolation to generate trajectories for off- nominal conditions. This approach was studied in the literature for Lunar Landing problems for example and yielded satisfying results [6]

Its main advantages are its capability to deal with asymmetric re-entry scenarios and its adapt-

Page 13/137

Figure 2.1: Callisto’s boosted landing

ability in case important modifications were brought to the nominal landing trajectory. It also shows the strong advantage of requiring very little computation power and therefore the possi- bility of frequent guidance calls to ensure robustness.

The second guidance strategy that is presented is a convex-optimisation based guidance algorithm.

This type of guidance method has gained a strong popularity in aerospace applications in recent years. It was applied to many landing cases, among others a 6-dof Mars landing guidance study [3], [4] and tested in real-flight conditions [19] and showed convincing results.

Indeed, the unique theoretical advantages of convex optimisation problems make them especially suited to guidance requirements. The challenge it poses though is the need to reformulate the trajectory optimisation problem as a convex one, as the vehicle dynamics and path constraints are generally not convex.

This master thesis proposes a convex guidance implementation, adapted to the Callisto boosted

landing, based on linearisation and . It also studies the possibility of

combining this strategy with engine ignition and extinction logics, to be able to study the

boosted landing as a whole.

Page 14/137

3. Mission presentation

3.1 Vehicle specifications

Callisto is a demonstrator co-developed by the French, German and Japanese space agencies (CNES, DLR and JAXA). It is a small-scale launcher propelled with a Japanese variable-thrust engine which aim is to test the recovery and reuse of the first-stage of a rocket. It is developed in preparation of the next demonstrator Themis, which will be equipped with a Prometheus engine. Callisto’s and Themis’ role is to pave the way for the development of Ariane NEXT (successor of Ariane6).

. Its overall aspect can be seen in Figure 3.1

Figure 3.1: Overall aspect of Callisto with fins and pods

Callisto’s nominal flight can be seen in Figure 3.2. After take-off, a ballistic phase will lead the vehicle to the apogee of its trajectory. It will there perform a flip manoeuvre and will start descending guided by its fins. About before touchdown, the engine will ignite again and the boosted landing will ensure that the vehicle lands safely on a pad on the ocean. On the trajectory Figure, the blue parts represent the maximum boost phases and the red parts the phases with the engine off.

One of the key features of Callisto is therefore its ability to land safely and its GNC systems are critical for the success of the mission.

The flip manoeuvre and descent are the flight phases requiring the most precision. They therefore

require important technical innovation as they have not yet been demonstrated on a launcher

Page 15/137

Figure 3.2: Overview of the overall Callisto trajectory

the size of Callisto.

3.2 Introduction to guidance and its role in the GNC system

This master thesis focuses on establishing a guidance logic for the boosted landing.

Guidance is part of the Guidance-Navigation-Control system of the vehicle. The GNC system role is to control the movement of the vehicle, as it cannot be performed manually. GNC stands for:

• Guidance, which the re-computation of the best (or at least feasible) trajectory during the flight given the position of the vehicle at a certain time and the desired target, while respecting control capacity and physical constraints.

• Navigation, which is the determination of the position and velocity of the vehicle at a given instant.

• Control, which is the steering of the controls (fins or thruster vectoring) and adjustment of the mass-flow-rate to follow the guidance-provided trajectory, while maintaining the stability of the vehicle.

Indeed, the vehicle follows a nominal trajectory computed before flight and stored in the launcher. However, during flight, dispersions will occur, such as deviations in mass, wind, atmospheric parameters, realised controls, etc.

Those dispersions will result in the vehicle not being on the theoretical trajectory that was

Page 16/137 planned and such dispersions will only increase as time goes by if no correction manoeuvre is performed. Guidance therefore has to take into account these dispersions and compute a new trajectory that will allow the vehicle to rejoin the nominal trajectory as fast as possible. This logic can be seen in Figure 3.3.

Figure 3.3: Logic of a guidance algorithm during landing

Guidance is continuously communicating with Navigation and Control in the GNC system:

it retrieves the vehicle position and velocity from Navigation and computes a trajectory from this data.

A trajectory is composed of a list of positions, velocities and all time instants, but also of controls, such as acceleration, angle of attack or mass flow rate. Those controls are then communicated to the pilot (Control), which will translate them into actual steering controls, such as fins angle or engine mass flow rate. Those orders will result in a modification of the vehicle position and velocity, which will then be assessed again by Navigation and the GNC loop can start again.

The interfacing logic between the different GNC subsystems can be seen in Figure 3.4.

Guidance can therefore be summarized as the real-time re-evaluation of the trajectory based on observable parameters (speed and position) and non-observable ones (vehicle properties, en- vironment) and the determination of controls that are given to the pilot.

It needs to be precise and robust to face dispersions, to minimize efficiently fuel consumption.

It also needs to ensure that the controls it gives the pilot are applicable and that the trajectory is controllable. Moreover, its capacity to be taken on-board is critical as it needs to be fast and reliable, as well as compatible with the on-board systems.

More specifically, for the boosted landing, the specific challenges guidance faces are the unsteadi-

ness of the physical forces on the vehicle which vary quickly. The flight time is also very short,

Page 17/137

Figure 3.4: Guidance role in the GNC system

so the navigation parameters vary very quickly.

3.3 General optimisation problem for the boosted phase

Guidance role is to re-optimize a trajectory at a given time, given navigation parameters in a very short time.

This computation is based on a nominal trajectory that the vehicle needs to follow as close as possible. The nominal trajectory is computed off-line, using dedicated optimisation software.

An optimal trajectory is a trajectory that minimizes fuel consumption and allows the vehicle to land safely in a controllable way. The global optimization problem is described in Figure 3.5.

The criterion that needs to be minimized is the fuel consumption, while respecting some constraints (fins and engine steering constraints as well as structural and stability constraints through ). In order to do so, the optimization has some parameters that vary, which are the final time, the mass flow rate, the angle of attack and the bank angle. Those parameters impact the vehicle behaviour through the vehicle dynamics, consisting of forces (gravity, thrust and aerodynamics). The vehicle is entirely determined by 6 states: 3 positions and 3 velocities, and 3 controls are considered.

This optimization problem is solved using dedicated CNES tools, . It cannot be

Page 18/137

Figure 3.5: Overview of the global optimization problem

solved analytically as the dynamics equation are quite complex and depend on vehicle parame- ters that are not necessarily expressed in an analytical way. Direct methods are therefore used through the software .

.

.

. This allows the vehicle to focus entirely on verticality during the last seconds and helps the GNC system to control the vehicle.

3.4 Nominal boosted trajectory

The optimal nominal trajectory, solution of the problem described above, can be seen in Figure 3.6 and its main parameters in Figure 3.7

One can see that the nominal trajectory is very smooth. It was computed in such a way that . . The mass-flow rate is also optimized as the maximum thrust allowed is used during the whole trajectory until the way-point, and then

the robustness finish the landing with a minimum thrust

final boost, on a vertical trajectory.

Page 19/137

Figure 3.6: Nominal trajectory with way-point

3.5 Overview of already existing guidance algorithms and their driving prin- ciple

Many guidance laws have been developed and implemented on launchers, applied to ascent phases, rendez-vous manoeuvres and descent phases. A global overview of guidance strategies cannot be described in this report. However, two main strategies are briefly introduced, to present the classical logic of a guidance strategy.

The first one is analytical and was used among other launchers on the guidance system of SaturnV that sent men to the Moon. Despite its apparent simplicity, it therefore shows strong results [12].

3.5.1 Tangent linear steering law

Based on very simplified assumptions, the tangent linear steering law guidance can be solved

analytically, which ensures its efficiency in terms of calculation burden. The steering law ignores

the aerodynamics forces, assumes a constant thrust and only aims at optimising the thrust angle

Page 20/137

Figure 3.7: Nominal trajectory with way-point flight parameters

to reach a circular orbit at a certain altitude. The model seen by the guidance can be graphically represented as in Figure 3.8.

Figure 3.8: Simplified guidance dynamics model: optimal angle to join a circular orbit

The launcher is here considered only as a point. Drag and lift forces are ignored and the

Page 21/137 thrust is assumed constant. Its goal is to join a circular orbit using the least possible fuel.

Since the thrust is assumed constant, the optimisation problem is equivalent to minimising the final time and can be written as the following:

minimize

α

t

such that ∀t ∈ [t

, t

], X(t) = ˙



















z = v ˙

x = v ˙

v ˙

=

0+ ˙mt

cos α v ˙

=

0+ ˙mt

sin α − g t

= 0

x(0) = x

, x(t

) = h z(0) = z

v

(0) = v

, v

(t

) = 0 v

(0) = v

, v

(t

) = v

(1)

This optimisation problem can be solved analytically sing the Pontryagin maximum principle.

The Hamiltonian can be written as:

H = hλ, f i = λ

v

+ λ

v

+ λ

T

m cos α + λ



T

m sin α − g



(2) Using the costate equations, one gets:

∃c

, c

, c

, c

∈ R,



















λ

= c

λ

= c

λ

= −c

t + c

λ

= −c

t + c

(3)

The optimal control can then be computed using the derivative of the Hamiltonian as







cos α = √

(−c1t+c3)²+(−c2t+c4)²

sin α = √

(−c1t+c3)²+(−c2t+c4)²

(4)

Using the boundary conditions gives:

tan α = − c

c

t + c

c

(5)

Page 22/137 The only step left is to solve the two-point boundary-value problem (2pbvp). In this case it can be found analytically or numerically without representing an important computation burden.

This guidance strategy is extremely efficient and gives an analytical expression of the steer- ing law given navigation parameters. However, this is only possible because the dynamics and optimisation problem were simplified drastically thanks to some strong assumptions.

It may seem that this guidance strategy would be too simplistic for a launcher GNC system.

However, frequent calls to the guidance may bring enough robustness and the simplicity of the guidance model make it easy to implement, fast and reliable.

Simple models are therefore sometimes best when developing a guidance strategy.

3.5.2 Predictor-corrector guidance law

The linear tangent law is a representative example of analytical guidance models. Other numerical models can be implemented, which allow for more complex dynamics equations. It is

the case of the predictor-corrector guidance algorithm, Its

driving principles are here briefly described. More details can be found in [11] and [10].

One is being given a vehicle with state variable X = [r, V, θ, φ, γ, ψ], with r radial distance to Earth centre, V velocity magnitude, θ longitude, φ latitude, γ flight path angle and ψ the velocity heading angle.

The vehicle dynamics are written in differential form as the following:

X = f (X) ˙ (6)

Throughout the flight, the vehicle needs to follow feasible trajectory constraints f

(X) ≤ f

, i = 1, · · · , N . Those constraints are then rewritten as velocity-dependent upper bounds on the bank angle magnitude (noted as σ). One therefore now has







|σ| ≤ σ

(V )

|∆α| ≤ ∆α

(7)

where ∆α is the angle of attack adjustment, bounded by an upper value.

There are also terminal constraints on the vehicle state, more specifically on position (alti-

tude, longitude, latitude) and velocity. Those constraints are then rewritten only as a function

of altitude and velocity.

Page 23/137 The longitude and latitude constraints are first rewritten as one:

s

(t

) = s

(8)

where s

is a function of longitude and latitude defining the range-to-go from the vehicle to the target.

An energy parameter is then defined as the following: e =

−

, so that the terminal constraints are rewritten as:

s

(e

) = s

h(e

) = h

(9)

A predictor-corrector can then be implemented.

For each control variable, the profile of a control variable (bank angle σ and angle of attack α) is defined as the following, as a function of the energy variable e:

u(e) = u

(e) + ∆u

e − e

e

− e

(10)

u

is the reference profile of the control evolution as a function of e:

u

(e) = u

+ u

− u

e

− e

(e − e

) (11)

The guidance strategy then determines the magnitude profile of a control variable u, which satisfies the terminal conditions, functions of e as well as the feasible path constraints. This is done based on the prediction of the final state if no correction is applied, and then determin- ing the correction needed to satisfy the terminal conditions. The schematic of this guidance algorithm can be described graphically as shown in Figure 3.9

Figure 3.9: Predictor-Corrector guidance strategy schematics.

The guidance first predicts the error at e

, ∆h

and ∆s

. The combined terminal error J is then defined as follows:

J =



∆h

k

2

+



∆s

k

2

(12)

Page 24/137 with k

and k

specific weights on terminal errors. The corrector is then brought to a zero- finding problem, which can be solved numerically.

Note: this is a generic presentation of the predictor-corrector guidance strategy. The control profile, parameters definitions and numerical method used to solve the zero-finding problem need to be adjusted to the problem at hand.

This other guidance strategy is considerably different from the linear tangent steering law. It takes into account complex dynamics, and rewrites them as a function of an energy parameter, to simplify the computations of the optimal control. The problem is then finally brought to a simple zero-finding problem, which can be solved efficiently numerically.

Those are only two guidance strategies among many others. They aim to show that many

approaches can be taken when establishing a guidance algorithm. The best suited strategy will

depend on the mission and the vehicle computation power, as well as the number of controls

that need to be optimized.

Page 25/137

Page 26/137

4. Adaptive pseudo-spectral interpolation method for the boosted reentry

4.1 Driving principles of the methodology and existing results

Currently existing guidance algorithms are not necessarily the best suited for a boosted land- ing flight phase. Other strategies exist that have shown promising results and were therefore adapted to our case. One of them is the Adaptive Pseudo-Spectral Interpolation method which is presented in this section.

The methodology presented here stems out from a DLR paper that was published in 2016 [5].

The mathematical method described in the study was adapted to our case and its results tested.

4.1.1 Mathematical principles

The idea behind this strategy is quite intuitive. Instead of re-optimizing in-flight a trajec- tory in the guidance algorithm, one uses support trajectories with dispersed initial conditions to re-interpolate the vehicle pseudo-optimal trajectory. Interpolation is much simpler than op- timization in terms of computation burden, and is therefore much faster.

The vehicles enters its boosted landing phase with observable dispersions in speed and po- sition compared to its nominal trajectory, stemming from previous flight phases. It therefore needs to reassess its trajectory very efficiently to compensate for those dispersions and land safely. This operation has to be performed in much less time than needed off-line to compute a whole new trajectory.

This is made possible through the computation of off-line nominal trajectories for dispersions in each dimension and stores it in the guidance system. When the guidance is called, it will therefore look for the nominal trajectories the closest to its actual position and compute a new trajectory using those neighbours.

One can decompose the APSI algorithm in five phases, described in detail in [5].

Parameter space discretisation

As the guidance algorithm relies on neighbour trajectories to re-evaluate controls in-flight,

one needs to discretise the flight envelope in each dimension. As our vehicle is considered as a set

of 6 states (3 positions and 3 velocities), there are therefore 6 dimensions to consider. There can

be 6 different off-nominal navigation conditions, and therefore there has to be a discretisation

Page 27/137 of the flight envelope in each of these dimensions.

In order to do so, one constructs a 6-D hypercube, containing the nominal conditions, thus representing the 6-D uncertainties.

For each variable in position or speed, its value will be contained in a compact interval, which can be geometrically represented by a straight line, with edges being the extremal possible values. If one considers 2-D uncertainties, the state will then be inside a rectangle. For 3-D uncertainties, it will be a cube. Therefore, for 6-D uncertainties, it will be a hypercube, called a herexact.

Such parametric spaces can be seen in Figure 4.1.

(a) Example of 1-D dispersion (b) Example of 2-D dispersion

(c) Example of 3-D dispersion (d) Example of 6-D dispersion

Figure 4.1: Examples of dispersions in n dimensions

In our case, the hypercube considered can be described as the following, where each p rep-

Page 28/137 resents a parameter (speed or position) and the vertices of the hypercube:

p

= [h

, h

] p

= [θ

, θ

] p

= [φ

, φ

] p

= [γ

, γ

] p

= [ψ

, ψ

] p

= [V

, V

]

(13)

The initial mass is not considered as a dispersion as it cannot be assessed during flight.

The initial conditions of the nominal trajectory will then be at the hypercentre of the hypercube described by those equations.

Trajectory generation

The second step of the algorithm, which is also performed off-line, is to compute an optimal trajectory for each point of the discretised parameter space. If one chooses to discretise the pa- rameter space in only two points in each dimension, then the initial conditions of the trajectories will be the 2

vertices of the hypercube. If one chooses to discretise the parameter space into a finer grid, with 3 points for each parameter for example, then one will have to compute 3

trajectories with initial conditions at each vertex of the hexeract and at the middle of each edge for example. Those trajectories corresponding to each set of initial conditions are then stored in the guidance algorithm. The next three steps are then performed on-line, during the flight.

Reference subspace selection

When the vehicle enters its boosted phase, the navigation gives the guidance a set of 6 initial conditions. The guidance then has to select the 2

closest trajectories in the discretised parameter space. This operation can be represented as in Figure 4.2 as an example in 3-D.

Figure 4.2: Reference subspace selection for 3-D dispersions

The closest trajectories are here defined as the ones with the least distance in the initial

Page 29/137 conditions for each dimension. If the parameter space is a hypercube, the reference subspace will then be a smaller hypercube, included in the parameter space hypercube.

This is the first step of the APSI guidance algorithm, graphically represented in Figure 4.3

Figure 4.3: Adaptive Pseudo-spectral interpolation principle- step 1 Low density interpolation

When the guidance has determined the 2

trajectories the closest to its initial conditions, it will interpolate the positions, controls and final time along those trajectories to determine the pseudo-optimal trajectory for the vehicles initial conditions inside the smaller hypercube.

One trajectory is made of 6 state variables (3 speeds, 3 positions) and 3 controls for each flight instant, as well as one final time. There are therefore 9 variables that need to be interpolated for each flight instant and one variable that needs only to be determined at one instant. This step can be represented in Figure 4.4.

However, if one chose to interpolate those 9 variables at each flight instant, the computing power needed to perform this operation would be much too important. Therefore, those 9 vari- ables are only interpolated at specific flight times along the trajectories, and will then be used to reconstitute a continuous trajectory.

The choice of the N flight instants at which one determines the 9 variables is then crucial. In- deed, they need to be adequately chosen so that the reconstituted continuous trajectory is close enough to the interpolated one.

To determine those flight instants, one chooses to take the flight times that correspond to the

Radau-Legendre roots. The Radau-Legendre polynomials are defined on [−1, 1] as the following:

Page 30/137

Figure 4.4: Adaptive Pseudo-spectral interpolation principle- step 2

∀N ≥ 0, ∀τ ∈ [−1, 1], L

(τ ) = 1 2

N !

d

dτ



(τ

− 1)

∀N ≥ 1, ∀τ ∈ [−1, 1]; R

(τ ) = L

(τ ) − L

(τ )

(14)

where the L

are the Legendre polynomials of order N and R

the Radau-Legendre poly- nomials. The roots of the Radau-Legendre polynomials are therefore the points τ ∈ [−1, 1] for which L

(τ ) = L

(τ ) for a given N . The behaviour of the Legendre polynomials can be seen in Figure 4.5.

Figure 4.5: Legendre polynomials on [−1, 1] for different n ∈ N

Page 31/137 For each τ such that R

(τ ) = 0, a multivariate interpolation is then conducted to determine the 9 variables of the interpolated trajectory. This step can be represented schematically as in Figure 4.6.

Figure 4.6: Adaptive Pseudo-spectral interpolation principle- step 3

The multivariate interpolation function is therefore a function taking in 9 · 2

variables corresponding to the controls and states of each reference subspace trajectory at a specific instant, yielding 9 interpolated variables corresponding to the interpolated trajectory.

If one notes x

∈ [p

, p

], i = 1, ..., 6 each of the off-nominal initial conditions, one can denote for each state or control j as f

the function

f

:[−1, 1] × R

→ R

(τ

, p) 7→ SV

(15)

where SV is the 2

supporting values of the interpolation, by oppostion to the p

, which are the supporting points. The multivariate interpolation, noted as f

for each variable j, is therefore

f

:[−1, 1] × R

→ R (τ

, x) 7→ (SV, x) 7→ s

(16)

where s is the state or control of each of the 9 variables j that need to be interpolated.

The multivariate interpolation method can be chosen freely. In this case, one will perform

an interpolation via a tensor spline. This operation can be written as the following:

Page 32/137

f

(τ

, x) = s

(x) =

m1

X

i1=1

. . .

m6

X

i6=1

c

B

(x

) . . . B

(x

) (17)

where B

is the i

B-spline of order k for a knot vector t = (t

)

. The coefficients c

is chosen such that

∀i ∈ [1, , 6], ∀k ∈ [1, . . . , N ], s

(p

) = f

(τ

, p

) (18) As one can see, since our support points constitute a grid, each B-spline only depends on a sin- gle variable, the 6-D interpolation problem can be decomposed in 6 1-D interpolation problems.

An efficient algorithm for solving each problem is the de Boor algorithm, with k

= 2, i = 1, . . . , 6 and t

= (p

, p

, p

, p

).

A graphical representation of such an interpolation method can be seen in Figure 4.7 for a 2-D problem. The support points are the grid, the supporting values the green dots. Along the first dimension p

, univariate interpolation is made along each edge. Along the second dimension p

, a second univariate interpolation is made, using the results from the precedent operation.

Figure 4.7: Example of a multivariate interpolation with a decomposition in 1-D problems LD-HD Pseudo-spectral conversion

The multivariate interpolation applied at the precedent step yields 6 positions, 3 controls at N

specific flight instants. One now needs to convert those discrete variables into continuous functions of time defined on [0, t

].

The first step is to compute the final time of the interpolated trajectory. This is done using the multivariate function defined in the precedent section, at t = 0, using the final time of the 2

support trajectories.

Pseudo-spectral approximation allows one to approximate polynomial expressions of contin- uous functions sampled in N + 1 points at τ

, k ∈ [0, . . . , N ]:

∀τ, f (τ ) ≈

N

X

i=0

f

P

(τ ) (19)

Page 33/137 with

∀τ, ∀i ∈ [0, . . . , N ], P

(τ ) =

N

Y

k=0,k6=i

τ − τ

τ

− τ

(20)

The advantage of having taken the roots of Radau-Legendre polynomials as collocation nodes for the multivariate interpolation, is that pseudo-spectral approximation using nodes sampled at those roots guarantees a uniform convergence, and avoids the Runge phenomenon that can occur when approximating functions with unbounded derivatives.

A graphical representation of Runge phenomenon with evenly distributed collocation nodes can be seen in Figure 4.8. How approximating the function using Radau-Legendre roots smooths the result can be seen in Figure 4.9 for the same function. More details about the Runge phenomenon and how it can be mitigated can be found in [17].

Figure 4.8: Pseudo-spectral approximation of the Runge function with evenly distributed collo- cation nodes

Using those roots is therefore a good strategy for rebuilding our continuous trajectory. How- ever, the Radau-Legendre roots are in the [−1, 1] interval, and one needs to yield a trajectory in the [0, t

] interval. As simple affine correspondence can be built as the following:

t = t

2 (τ + 1) (21)

We therefore have a bijective correspondence between the spectral time τ and the real time t.

Starting from the 9 variables computed in each N

+ 1 nodes, one needs to approximate a function for each variable, so that one has the value of this variable for a much finer time grid.

Noting X

and U

the states and controls vectors and having n

states and n

controls, one

Page 34/137

Figure 4.9: Pseudo-spectral approximation of the Runge function with collocation nodes chosen to mitigate the Runge phenomenon

has

T

= X

U

!

= X

, ..., X

U

, ..., U

!

∈ R

(22) Our goal is therefore to have

T

= X

U

!

= X

, ..., X

U

, ..., U

!

∈ R

(23)

representing the values of the 9 variables (states and controls) at each of the N

+1 time nodes.

The idea is then to take N

large enough to have trajectory points in many time instants, that is to discretise the time vector to be sufficiently fine. One also need to convert τ

, the pseudo-spectral time vector, into t

, the real HD time vector.

This operation can be done in a matrix form. Indeed, one has

∀τ ∈ [−1, 1], f (τ ) ≈

NLD

X

i=0

f

Y

k=0,k6=i

τ − τ

τ

− τ

(24)

Here, f

can be seen as the p

row of T

. Moreover, one can compute f in the N

+ 1 high density pseudo-spectral time nodes ˜ τ

. One can therefore write

∀m ∈ [0, N

], ∀p ∈ [1, n

+ n

], T ˜

(˜ τ

) ≈

NLD

X

i=0

T

NLD

Y

k=0,k6=i

τ ˜

− τ

τ

− τ

(25)

Page 35/137 This equation can be expressed for all p in matrix form:

T ˜

= T

M

(26)

with

M

=











QNLD

k=1

˜ τ0−τ_k

τ0−τ_k

· · ·

0−τ_k

.. . . .. .. .

QNLD−1 k=0

˜ τ0−τ_k

τ_nLD−τ_k

· · ·

nLD−τ_k











(27)

Furthermore, the HD time-vector can be computed from the HD pseudo-spectral time vector τ ˜

in the following way:

˜ t

= t

2 (˜ τ

+ 1), ∀m ∈ [0, N

] (28)

A graphical representation of the LD-HD conversion can be seen in Figure 4.10.

Figure 4.10: Adaptive Pseudo-spectral interpolation principle- step 4

4.1.2 Previous evaluations: DLR study

This method was adapted from a DLR study focusing on the unpowered glided descent. It showed very promising results on reducing the error to the nominal trajectory at the end of the glided descent.

Indeed, out of the 1000 dispersed cases that were presented in their study, the error correction performed by the guidance seemed to yield precision levels satisfying enough to make all trajec- tories acceptable. More details about their results can be found in [5].

However, it was applied to a glided descent case, that is without the engine on. That means,

Page 36/137 first of all, that only two controls were considered. Moreover, the overall glided descent trajec- tory lasts also much longer and is smoother. The specific constraints on guidance for the boosted landing (stringent accuracy requirements at touchdown, forces unsteadiness, strong variation of flight parameters) do not apply to the glided descent, so it is highly possible that this strategy could be implemented for steady flight phases but prove themselves to be inefficient for the final boosted landing.

It is therefore important to adapt it and implement it for our flight phase to see if it could yield also satisfying results to our flight phase.

4.2 Implementation of the methodology 4.2.1 Flight envelope characterization

Optimal trajectory generation

The launch vehicle dynamics are here modelled in 3 dimensions. The state of the launcher consists therefore of 6 variables, characterizing its position and velocity.

The position is characterized using the geodesic coordinates, and the speed using its norm, slope with earth surface and its azimuth, as shown in Figures 4.11 and 4.12.

X =

























h θ φ γ V ψ

























(29)

Figure 4.11: Geodesic coordinates