From 03825ea5584d9ef97c1b9e8c35847d9f0175a8ed Mon Sep 17 00:00:00 2001 From: Connor Date: Sun, 23 Jan 2022 16:13:38 -0700 Subject: [PATCH] Section 3 --- LaTeX/thesis.tex | 996 +++++++++++++++++++++++++---------------------- 1 file changed, 527 insertions(+), 469 deletions(-) diff --git a/LaTeX/thesis.tex b/LaTeX/thesis.tex index 8df1e3a..88369fc 100644 --- a/LaTeX/thesis.tex +++ b/LaTeX/thesis.tex @@ -17,20 +17,20 @@ \newcommand{\sectionbreak}{\clearpage} \titleformat{\section} - {\bfseries\fontspec{Roboto}\LARGE} - {\thesection} - {0.5em} - {} + {\bfseries\fontspec{Roboto}\LARGE} + {\thesection} + {0.5em} + {} \titleformat{\subsection} - {\fontspec{Roboto}\Large} - {\thesubsection} - {0.5em} - {} + {\fontspec{Roboto}\Large} + {\thesubsection} + {0.5em} + {} \titleformat{\subsubsection} - {\fontspec{Roboto}\large} - {\thesubsubsection} - {0.5em} - {} + {\fontspec{Roboto}\large} + {\thesubsubsection} + {0.5em} + {} \newcommand{\thesisTitle}{Designing Optimal Low-Thrust Interplanetary Trajectories Utilizing Monotonic Basin Hopping} @@ -40,495 +40,553 @@ Monotonic Basin Hopping} \begin{document} - \title{\thesisTitle} - \author{Richard Connor Johnstone \\ - B.S., University of Kentucky, 2016} + \title{\thesisTitle} + \author{Richard Connor Johnstone \\ + B.S., University of Kentucky, 2016} - \maketitle + \maketitle - \vspace{3.5in} + \vspace{3.5in} - \begin{adjustwidth}{100pt}{100pt} - \begin{center} - A thesis submitted to the Faculty of the Graduate School of the University of Colorado in - partial fulfillment of the requirements for the degree of Master of Science Department of - Aerospace Engineering Sciences \\ 2022 - \end{center} - \end{adjustwidth} + \begin{adjustwidth}{100pt}{100pt} + \begin{center} + A thesis submitted to the Faculty of the Graduate School of the University of Colorado in + partial fulfillment of the requirements for the degree of Master of Science Department of + Aerospace Engineering Sciences \\ 2022 + \end{center} + \end{adjustwidth} - \newpage + \newpage - This will be the copyright page. + This will be the copyright page. - \newpage + \newpage - \begin{center} - This thesis entitled: + \begin{center} + This thesis entitled: - \thesisTitle + \thesisTitle - has been approved for the Department of Aerospace Engineering Sciences + has been approved for the Department of Aerospace Engineering Sciences - \end{center} + \end{center} - \centerline{\begin{minipage}{4in} - \begin{center} + \centerline{\begin{minipage}{4in} + \begin{center} - \vspace{1.5in} - \hrule - Dr. Natasha Bosanac + \vspace{1.5in} + \hrule + Dr. Natasha Bosanac - \vspace{1.5in} - \hrule - Dr. Daniel J. Scheeres + \vspace{1.5in} + \hrule + Dr. Daniel J. Scheeres - \vspace{1.5in} - \hrule - Dr. Kathryn Davis + \vspace{1.5in} + \hrule + Dr. Kathryn Davis - \vspace{1.5in} + \vspace{1.5in} - \end{center} - \end{minipage}} + \end{center} + \end{minipage}} - \hspace{4.25in} Date: \hrulefill + \hspace{4.25in} Date: \hrulefill - The final copy of this thesis has been examined by the signatories, and we find that both the - content and the form meet acceptable presentation standards of scholarly work in the above - mentioned discipline. + The final copy of this thesis has been examined by the signatories, and we find that both the + content and the form meet acceptable presentation standards of scholarly work in the above + mentioned discipline. - \newpage + \newpage - %TODO: This should be better - Dedicated to some people. + %TODO: This should be better + Dedicated to some people. - \newpage + \newpage - Richard Connor Johnstone + Richard Connor Johnstone - \thesisTitle - - %TODO: Don't directly copy Bryce's formatting - Thesis directed by Dr. Natasha Bosanac - - %TODO: This is just a quick abstract for now. Rewrite more towards the end. - \begin{abstract} - - There are a variety of approaches to finding and optimizing low-thrust trajectories in - interplanetary space. This thesis analyzes one such approach, Sims-Flanagan transcriptions, and - its applications in a multiple-shooting non-linear solver for the purpose of finding valid - low-thrust trajectory arcs between planets given poor initial conditions. These valid arcs are - then fed into a Monotonic Basin Hopping (MBH) algorithm, which combines these arcs in order to - find and optimize interplanetary trajectories, given a set of flyby planets. This allows for a - fairly rapid searching of a very large solution space of low-thrust profiles via a medium - fidelity inner-loop solver and a well-suited optimization routine. The trajectories found by - this method can then be optimized further by feeding the solutions back, once again, into the - non-linear solver, this time allowing the solver to perform optimization. - - \end{abstract} - - \newpage - - \tableofcontents - - \listoffigures - - \newpage + \thesisTitle + + %TODO: Don't directly copy Bryce's formatting + Thesis directed by Dr. Natasha Bosanac + + %TODO: This is just a quick abstract for now. Rewrite more towards the end. + \begin{abstract} + + There are a variety of approaches to finding and optimizing low-thrust trajectories in + interplanetary space. This thesis analyzes one such approach, Sims-Flanagan transcriptions, and + its applications in a multiple-shooting non-linear solver for the purpose of finding valid + low-thrust trajectory arcs between planets given poor initial conditions. These valid arcs are + then fed into a Monotonic Basin Hopping (MBH) algorithm, which combines these arcs in order to + find and optimize interplanetary trajectories, given a set of flyby planets. This allows for a + fairly rapid searching of a very large solution space of low-thrust profiles via a medium + fidelity inner-loop solver and a well-suited optimization routine. The trajectories found by + this method can then be optimized further by feeding the solutions back, once again, into the + non-linear solver, this time allowing the solver to perform optimization. + + \end{abstract} + + \newpage + + \tableofcontents + + \listoffigures + + \newpage - \section{Introduction} + \section{Introduction} - Continuous low-thrust arcs utilizing technologies such as Ion propulsion, Hall thrusters, and - others can be a powerful tool in the design of interplanetary space missions. They tend to be - particularly suited to missions which require very high total change in velocity or $\Delta V$ - values and take place over a particularly long duration. Traditional impulsive thrusting - techniques can achieve these changes in velocity, but they typically have a far lower specific - impulse and, as such, are much less efficient and use more fuel, costing the mission valuable - financial resources that could instead be used for science. Because of their inherently high - specific impulse (and thus efficiency), low-thrust fuels are well-suited to interplanetary - missions. - - For instance, low thrust ion propulsion was used on the Bepi-Colombo, Dawn, and Deep - Space 1 missions. In general, anytime an interplanetary trajectory is posed, it is advisable to - first explore the possibility of low-thrust technologies. In an interplanetary mission, the - primary downside to low-thrust orbits (that they require significant time to achieve large - $\Delta V$ changes) is made irrelevant by the fact that interplanetary trajectories take such a - long time as a matter of course. - - Another technique often leveraged by interplanetary trajectory designers is the gravity assist. - Gravity assists cleverly utilize the inertia of a large planetary body to ''slingshot`` a - spacecraft, modifying the direction of its velocity with respect to the central body, the Sun. - This technique lends itself very well to impulsive trajectories. The gravity assist maneuver - itself can be modeled very effectively by an impulsive maneuver with certain constraints, placed - right at the moment of closest approach to the (flyby) target body. Because of this, - optimization with impulsive trajectories and gravity assists are common. - - However, there is no physical reason why low-thrust trajectories can't also incorporate gravity - assists. The optimization problem becomes much more complicated. The separate problems of - optimizing flyby parameters (planet, flyby date, etc.) and optimizing the low-thrust control - arcs don't combine very easily. In this paper, a technique is explored by setting the - dual-problem up as a Hybrid Optimal Control Problem (HOCP). - - This thesis will explore these concepts in a number of different sections. Section - \ref{traj_opt} will explore the basic principles of trajectory optimization in a manner agnostic - to the differences between continuous low-thrust and impulsive high-thrust techniques. Section - \ref{low_thrust} will then delve into the different aspects to consider when optimizing a low - thrust mission profile over an impulsive one. Section \ref{interplanetary} provides more detail - on the interplanetary considerations, including force models and gravity assists. Section - \ref{algorithm} will cover the implementation details of the HOCP optimization algorithm - developed for this paper. Finally, section \ref{results} will explore the results of some - hypothetical missions to Saturn. - - \section{Trajectory Optimization} \label{traj_opt} - - Trajectory optimization is concerned with a narrow problem (namely, optimizing a spaceflight - trajectory to an end state) with a wide range of possible techniques, approaches, and even - solutions. In this section, the foundations for direct optimization of these sorts of problems - will be explored by first introducing the Two-Body Problem, then an algorithm for directly - solving for states in that system, then exploring approaches to Non-Linear Problem (NLP) solving - in general and how they apply to spaceflight trajectories. - - \subsection{The Two-Body Problem} - The motion of a spacecraft in space is governed by a large number of forces. When planning and - designing a spacecraft trajectory, we often want to use the most complete (and often complex) - model of these forces that is available. However, in the process of designing these - trajectories, we often have to compute the path of the spacecraft many hundreds, thousands, or - even millions of times. Utilizing very high-fidelity force models that account for aerodynamic - pressures, solar radiation pressures, multi-body effects, and many others may be infeasible - for the method being used if the computations take too long. - - Therefore, a common approach (and the one utilized in this implementation) is to first look - simply at the single largest force governing the spacecraft in motion, the gravitational force - due to the primary body around which it is orbiting. This can provide an excellent - low-to-medium fidelity model that can be extremely useful in categorizing the optimization - space as quickly as possible. In many cases, including the algorithm used in this paper, it is - unlikely that local cost-function minima would be missed due to the lack of fidelity of the - Two Body Problem. - - In order to explore the Two Body Problem, we must first examine the full set of assumptions - associated with the force model. Firstly, we are only concerned with the nominative two - bodies: the spacecraft and the planetary body around which it is orbiting. Secondly, both of - these bodies are modeled as simple point masses. This removes the need to account for - non-uniform densities and asymmetry. The third assumption is that the mass of the spacecraft - ($m_2$) is much much smaller than the mass of the planetary body ($m_1$) and enough so as to be - considered negligible. The only force acting on this system is then the force of gravity that - the primary body enacts upon the secondary. Lastly, we'll assume a fixed inertial frame. This - isn't necessary for the formulation of a solution, but will simplify the derivation. - - Reducing the system to two point masses with a single gravitational force acting between them - (and only in one direction) we can model the force on the secondary body as: - - \begin{equation} - \ddot{\vec{r}} = - \frac{G \left( m_1 + m_2 \right)}{r^2} \frac{\vec{r}}{\left| r \right|} - \end{equation} - - Where $\vec{r}$ is the position of the spacecraft, $G$ is the universal gravitational - parameter, $m_1$ is the mass of the planetary body, and $m_2$ is the mass of the spacecraft. - Due to our assumption that the mass of the spacecraft is significantly smaller than the mass - of the primary body ($m_1 >> m_2$) we can reduce that formulation to simply: - - \begin{equation} - \ddot{\vec{r}} = - \frac{\mu}{r^2} \hat{r} - \end{equation} - - Where $\mu = G m_1$ is the specific gravitational parameter for our primary body of interest. - - \subsubsection{Kepler's Laws and Equations} - - % TODO: Can I segue better from 2BP to Keplerian geometry? - - Now that we've fully qualified the forces acting within the Two Body Problem, we can concern - ourselves with more practical applications of this as a force model. It should be noted, - firstly, that the spacecraft's position and velocity (given an initial position and velocity - and of course the $\mu$ value of the primary body) is actually analytically solvable for all - future points in time. This can be easily observed by noting that there are three - one-dimensional equations (one for each component of the three-dimensional position) and - three unknowns (the three components of the second derivative of the position). - - In the early 1600s, Johannes Kepler produced just such a solution. By taking advantages of - what is also known as ``Kepler's Laws'' which are: - - \begin{enumerate} - \item Each planet's orbit is an ellipse with the Sun at one of the foci. This can be - expanded to any orbit by re-wording as ``all orbital paths follow a conic section - (circle, ellipse, parabola, or hyperbola) with a primary body at one of the foci''. - \item The area swept out by the imaginary line connecting the primary and secondary - bodies increases linearly with respect to time. This implies that the magnitude of the - orbital speed is not constant. - \item The square of the orbital period is proportional to the cube of the semi-major - axis of the orbit, regardless of eccentricity. Specifically, the relationship is: $T = 2 - \pi \sqrt{\frac{a^3}{\mu}}$ where $T$ is the period and $a$ is the semi-major axis. - \end{enumerate} - - \subsection{Analytical Solutions to Kepler's Equations} - - Kepler was able to produce an equation to represent the angular displacement of an orbiting - body around a primary body as a function of time, which we'll derive now for the elliptical - case. Since the total area of an ellipse is the product of $\pi$, the semi-major axis, and - the semi-minor axis ($\pi a b$), we can relate (by Kepler's second law) the area swept out - by an orbit as a function of time: - - \begin{equation}\label{swept} - \frac{\Delta t}{T} = \frac{k}{\pi a b} - \end{equation} - - This leaves just one unknown variable $k$, which we can determine through use of the - geometric auxiliary circle, which is a circle with radius equal to the ellipse's semi-major - axis and center directly between the two foci, as in Figure~\ref{aux_circ}. - - \begin{figure} - \centering - \includegraphics[width=0.8\textwidth]{LaTeX/fig/kepler} - \caption{Geometric Representation of Auxiliary Circle}\label{aux_circ} - \end{figure} - - In order to find the area swept by the spacecraft, $k$, we can take advantage of the fact - that that area is the triangle $k_1$ subtracted from the elliptical segment $PCB$: - - \begin{equation}\label{areas_eq} - k = area(seg_{PCB}) - area(k_1) - \end{equation} - - Where the area of the triangle $k_1$ can be found easily using geometric formulae: - - \begin{align} - area(k_1) &= \frac{1}{2} \left( ae - a \cos E \right) \left( \frac{b}{a} a \sin E \right) \\ - &= \frac{ab}{2} \left(e \sin E - \cos E \sin E \right) - \end{align} - - Now we can find the area for the elliptical segment $PCB$ by first finding the circular - segment $POB'$, subtracting the triangle $OB'C$, then applying the fact that an ellipse is - merely a vertical scaling of a circle by the amount $\frac{b}{a}$. - - \begin{align} - area(PCB) &= \frac{b}{a} \left( \frac{a^2 E}{2} - \frac{1}{2} \left( a \cos E \right) - \left( a \sin E \right) \right) \\ - &= \frac{abE}{2} - \frac{ab}{2} \left( \cos E \sin E \right) \\ - &= \frac{ab}{2} \left( E - \cos E \sin E \right) - \end{align} - - By substituting the two areas back into Equation~\ref{areas_eq} we can get the $k$ area - swept out by the spacecraft: - - \begin{equation} - k = \frac{ab}{2} \left( E - e \sin E \right) - \end{equation} - - Which we can then substitute back into the equation for the swept area as a function of - time (Equation~\ref{swept}): - - \begin{equation} - \frac{\Delta t}{T} = \frac{E - e \sin E}{2 \pi} - \end{equation} - - Which is, effectively, Kepler's equation. It is commonly known by a different form: - - \begin{align} - M &= E - e \sin E \\ - &= \sqrt{\frac{\mu}{a^3}} \Delta t - \end{align} - - Where we've defined the mean anomaly as $M$ and used the fact that $T = - \sqrt{\frac{a^3}{\mu}}$. This provides us a useful relationship between Eccentric Anomaly - ($E$) which can be related to spacecraft position, and time, but we still need a useful - algorithm for solving this equation. - - \subsubsection{LaGuerre-Conway Algorithm} - For this application, I used an algorithm known as the LaGuerre-Conway algorithm, which was - presented in 1986 as a faster algorithm for directly solving Kepler's equation and has been - in use in many applications since. This algorithm is known for its convergence robustness - and also its speed of convergence when compared to higher order Newton methods. - - This thesis will omit a step-through of the algorithm itself, but the code will be present - in the Appendix. - - \subsection{Non-Linear Problem Optimization} - - Now we can consider the formulation of the problem in a more useful way. For instance, given a - desired final state in position and velocity we can relatively easily determine the initial - state necessary to end up at that desired state over a pre-defined period of time by solving - Kepler's equation. In fact, this is often how impulsive trajectories are calculated since, - other than the impulsive thrusting event itself, the trajectory is entirely natural. - - However, often in trajectory design we want to consider a number of other inputs. For - instance, a low thrust profile, a planetary flyby, the effects of rotating a solar panel on - solar radiation pressure, etc. Once these inputs have been accepted as part of the model, the - system is generally no longer analytically solvable, or, if it is, is too complex to calculate - directly. - - Therefore an approach is needed, in trajectory optimization and many other fields, to optimize - highly non-linear, unpredictable systems such as this. The field that developed to approach - this problem is known as Non-Linear Problem (NLP) Optimization. In these cases, generally - speaking, direct methods are impossible, so a number of tools have been developed to optimize - NLPs in the general case, via indirect, iterative methods. - - \subsubsection{Non-Linear Solvers} - For these types of non-linear, constrained problems, a number of tools have been developed - that act as frameworks for applying a large number of different algorithms. This allows for - simple testing of many different algorithms to find what works best for the nuances of the - problem in question. - - One of the most common of these NLP optimizers is SNOPT, which is a proprietary package - written primarily using a number of Fortran libraries by the Systems Optimization Laboratory - at Stanford University. It uses a sparse sequential quadratic programming approach. - - Another common NLP optimization packages (and the one used in this implementation) is the - Interior Point Optimizer or IPOPT. It can be used in much the same way as SNOPT and uses an - Interior Point Linesearch Filter Method and was developed as an open-source project by the - organization COIN-OR under the Eclipse Public License. - - Both of these methods utilize similar approaches to solve general constrained non-linear - problems iteratively. Both of them can make heavy use of derivative Jacobians and Hessians - to improve the convergence speed and both have been ported for use in a number of - programming languages, including in Julia, which was used for this project. - - This is by no means an exhaustive list, as there are a number of other optimization - libraries that utilize a massive number of different algorithms. For the most part, the - libraries that port these are quite modular in the sense that multiple algorithms can be - tested without changing much source code. - - \subsubsection{Linesearch Method} - As mentioned above, this project utilized IPOPT which leveraged an Interior Point Linesearch - method. A linesearch algorithm is one which attempts to find the optimum of a non-linear - problem by first taking an initial guess $x_k$. The algorithm then determines a step - direction (in this case through the use of automatic differentiation to calculate the - derivatives of the non-linear problem) and a step length. The linesearch algorithm then - continues to step the initial guess, now labeled $x_{k+1}$ after the addition of the - ``step'' vector and iterates this process until predefined termination conditions are met. - - In this case, the IPOPT algorithm was used, not as an optimizer, but as a solver. For - reasons that will be explained in the algorithm description in Section~\ref{algorithm} it - was sufficient merely that the non-linear constraints were met, therefore optimization (in - the particular step in which IPOPT was used) was unnecessary. - - \section{Low-Thrust Considerations} \label{low_thrust} - Highlight the differences between high and low-thrust mission profiles. - - \subsection{Low Thrust Overview} - Dive deeper into the concept of low thrust trajectories, highlighting the added trouble with - propagation/integration. - - \subsection{Sims-Flanagan Transcription} - Reveal the advantages of Sims-Flanagan transcription as an alternative to higher-fidelity - propagation models. Be sure to mention its uses in many legitimate places. - - \section{Interplanetary Trajectory Considerations} \label{interplanetary} - Highlight the problems with the 2BP in co-ordinating influences of extra bodies over an - interplanetary journey. - - \subsection{Patched Conics} - Describe the method of patched conics. - - \subsection{Gravity Assist Maneuvers} - Describe how a gravity assist maneuver would work in the framework of patched conics. Also - discuss the advantages of such a maneuver. - - \subsection{Multiple Gravity Assist Techniques} - Discuss the advantages of chaining together multiple gravity assists and highlight the - difficulties in choosing these assists. Here I can mention porkchop plots, Lambert's problem, - etc. Here I can also talk about Hybrid Optimal Control Problems. - - \subsection{Ephemeris Considerations} - I can quickly mention SPICE here and talk a bit about validation. - - % \section{Genetic Algorithms} - % I will probably give only a brief overview of genetic algorithms here. I don't personally know - % that much about them. Then in the following subsections I can discuss the parts that are - % relevant to the specific algorithm that I'm using. - - % \subsection{Decision Vectors} - % Discuss what a decision vector is in the context of an optimization problem. - - % \subsection{Selection and Fitness Evaluation} - % Discuss the costing being used as well as the different types of fitness evaluation that are - % common. Also discuss the concept of generations and ``survival''. - - % \subsubsection{Tournament Selection} - % Dive deeper into the specific selection algorithm being used here. - - % \subsection{Crossover} - % Discuss the concept of crossover and procreation in a genetic algorithm. - - % \subsubsection{Binary Crossover} - % Discuss specific crossover algorithm used here. - - % \subsubsection{Mutation} - % Discuss both the necessity for mutation and the mutation algorithm being used. - - \section{Algorithm Overview} \label{algorithm} - Highlight the algorithm at a high-level. This is likely where flowcharts and diagrams will go to - give a high-level overview. - - \subsection{Trajectory Composition} - Discuss briefly the nomenclature used in defining these trajectories. Currently this isn't - ``baked in'' to the code, so I have some freedom to adopt Englander's notation or use my own - (since my intended use case is a little simpler). - - \subsection{Inner Loop Implementation} - Give a better overview of the inner loop specifically. Probably this section will have a more - in-depth flowchart. - - \subsubsection{LaGuerre-Conway Kepler Solver} - Discuss how the LaGuerre-Conway algorithm is used in the code to provide a fundamental - ``natural trajectory'' between two quantized, but not necessarily close points. Mention - validation. - - \subsubsection{Sims-Flanagan Propagator} - Discuss how this algorithm can then be expanded by using SFT to propagate any number of - low-thrust steps over a specific arc. Mention validation. Here I can also mention the ``Sc'' - object and talk about how those parameters were chosen and effected the propagator. - - \subsubsection{Non-Linear Problem Solver} - Mention the package being used to solve NLPs and how it works, highlighting the trust region - method used and error-handling. Mention validation. + Continuous low-thrust arcs utilizing technologies such as Ion propulsion, Hall thrusters, and + others can be a powerful tool in the design of interplanetary space missions. They tend to be + particularly suited to missions which require very high total change in velocity or $\Delta V$ + values and take place over a particularly long duration. Traditional impulsive thrusting + techniques can achieve these changes in velocity, but they typically have a far lower specific + impulse and, as such, are much less efficient and use more fuel, costing the mission valuable + financial resources that could instead be used for science. Because of their inherently high + specific impulse (and thus efficiency), low-thrust fuels are well-suited to interplanetary + missions. + + For instance, low thrust ion propulsion was used on the Bepi-Colombo, Dawn, and Deep + Space 1 missions. In general, anytime an interplanetary trajectory is posed, it is advisable to + first explore the possibility of low-thrust technologies. In an interplanetary mission, the + primary downside to low-thrust orbits (that they require significant time to achieve large + $\Delta V$ changes) is made irrelevant by the fact that interplanetary trajectories take such a + long time as a matter of course. + + Another technique often leveraged by interplanetary trajectory designers is the gravity assist. + Gravity assists cleverly utilize the inertia of a large planetary body to ''slingshot`` a + spacecraft, modifying the direction of its velocity with respect to the central body, the Sun. + This technique lends itself very well to impulsive trajectories. The gravity assist maneuver + itself can be modeled very effectively by an impulsive maneuver with certain constraints, placed + right at the moment of closest approach to the (flyby) target body. Because of this, + optimization with impulsive trajectories and gravity assists are common. + + However, there is no physical reason why low-thrust trajectories can't also incorporate gravity + assists. The optimization problem becomes much more complicated. The separate problems of + optimizing flyby parameters (planet, flyby date, etc.) and optimizing the low-thrust control + arcs don't combine very easily. In this paper, a technique is explored by setting the + dual-problem up as a Hybrid Optimal Control Problem (HOCP). + + This thesis will explore these concepts in a number of different sections. Section + \ref{traj_opt} will explore the basic principles of trajectory optimization in a manner agnostic + to the differences between continuous low-thrust and impulsive high-thrust techniques. Section + \ref{low_thrust} will then delve into the different aspects to consider when optimizing a low + thrust mission profile over an impulsive one. Section \ref{interplanetary} provides more detail + on the interplanetary considerations, including force models and gravity assists. Section + \ref{algorithm} will cover the implementation details of the HOCP optimization algorithm + developed for this paper. Finally, section \ref{results} will explore the results of some + hypothetical missions to Saturn. + + \section{Trajectory Optimization} \label{traj_opt} + + Trajectory optimization is concerned with a narrow problem (namely, optimizing a spaceflight + trajectory to an end state) with a wide range of possible techniques, approaches, and even + solutions. In this section, the foundations for direct optimization of these sorts of problems + will be explored by first introducing the Two-Body Problem, then an algorithm for directly + solving for states in that system, then exploring approaches to Non-Linear Problem (NLP) solving + in general and how they apply to spaceflight trajectories. + + \subsection{The Two-Body Problem} + The motion of a spacecraft in space is governed by a large number of forces. When planning and + designing a spacecraft trajectory, we often want to use the most complete (and often complex) + model of these forces that is available. However, in the process of designing these + trajectories, we often have to compute the path of the spacecraft many hundreds, thousands, or + even millions of times. Utilizing very high-fidelity force models that account for aerodynamic + pressures, solar radiation pressures, multi-body effects, and many others may be infeasible + for the method being used if the computations take too long. + + Therefore, a common approach (and the one utilized in this implementation) is to first look + simply at the single largest force governing the spacecraft in motion, the gravitational force + due to the primary body around which it is orbiting. This can provide an excellent + low-to-medium fidelity model that can be extremely useful in categorizing the optimization + space as quickly as possible. In many cases, including the algorithm used in this paper, it is + unlikely that local cost-function minima would be missed due to the lack of fidelity of the + Two Body Problem. + + In order to explore the Two Body Problem, we must first examine the full set of assumptions + associated with the force model. Firstly, we are only concerned with the nominative two + bodies: the spacecraft and the planetary body around which it is orbiting. Secondly, both of + these bodies are modeled as simple point masses. This removes the need to account for + non-uniform densities and asymmetry. The third assumption is that the mass of the spacecraft + ($m_2$) is much much smaller than the mass of the planetary body ($m_1$) and enough so as to be + considered negligible. The only force acting on this system is then the force of gravity that + the primary body enacts upon the secondary. Lastly, we'll assume a fixed inertial frame. This + isn't necessary for the formulation of a solution, but will simplify the derivation. + + Reducing the system to two point masses with a single gravitational force acting between them + (and only in one direction) we can model the force on the secondary body as: + + \begin{equation} + \ddot{\vec{r}} = - \frac{G \left( m_1 + m_2 \right)}{r^2} \frac{\vec{r}}{\left| r \right|} + \end{equation} + + Where $\vec{r}$ is the position of the spacecraft, $G$ is the universal gravitational + parameter, $m_1$ is the mass of the planetary body, and $m_2$ is the mass of the spacecraft. + Due to our assumption that the mass of the spacecraft is significantly smaller than the mass + of the primary body ($m_1 >> m_2$) we can reduce that formulation to simply: + + \begin{equation} + \ddot{\vec{r}} = - \frac{\mu}{r^2} \hat{r} + \end{equation} + + Where $\mu = G m_1$ is the specific gravitational parameter for our primary body of interest. + + \subsubsection{Kepler's Laws and Equations} + + % TODO: Can I segue better from 2BP to Keplerian geometry? + + Now that we've fully qualified the forces acting within the Two Body Problem, we can concern + ourselves with more practical applications of this as a force model. It should be noted, + firstly, that the spacecraft's position and velocity (given an initial position and velocity + and of course the $\mu$ value of the primary body) is actually analytically solvable for all + future points in time. This can be easily observed by noting that there are three + one-dimensional equations (one for each component of the three-dimensional position) and + three unknowns (the three components of the second derivative of the position). + + In the early 1600s, Johannes Kepler produced just such a solution. By taking advantages of + what is also known as ``Kepler's Laws'' which are: + + \begin{enumerate} + \item Each planet's orbit is an ellipse with the Sun at one of the foci. This can be + expanded to any orbit by re-wording as ``all orbital paths follow a conic section + (circle, ellipse, parabola, or hyperbola) with a primary body at one of the foci''. + \item The area swept out by the imaginary line connecting the primary and secondary + bodies increases linearly with respect to time. This implies that the magnitude of the + orbital speed is not constant. + \item The square of the orbital period is proportional to the cube of the semi-major + axis of the orbit, regardless of eccentricity. Specifically, the relationship is: $T = 2 + \pi \sqrt{\frac{a^3}{\mu}}$ where $T$ is the period and $a$ is the semi-major axis. + \end{enumerate} + + \subsection{Analytical Solutions to Kepler's Equations} + + Kepler was able to produce an equation to represent the angular displacement of an orbiting + body around a primary body as a function of time, which we'll derive now for the elliptical + case. Since the total area of an ellipse is the product of $\pi$, the semi-major axis, and + the semi-minor axis ($\pi a b$), we can relate (by Kepler's second law) the area swept out + by an orbit as a function of time: + + \begin{equation}\label{swept} + \frac{\Delta t}{T} = \frac{k}{\pi a b} + \end{equation} + + This leaves just one unknown variable $k$, which we can determine through use of the + geometric auxiliary circle, which is a circle with radius equal to the ellipse's semi-major + axis and center directly between the two foci, as in Figure~\ref{aux_circ}. + + \begin{figure} + \centering + \includegraphics[width=0.8\textwidth]{LaTeX/fig/kepler} + \caption{Geometric Representation of Auxiliary Circle}\label{aux_circ} + \end{figure} + + In order to find the area swept by the spacecraft, $k$, we can take advantage of the fact + that that area is the triangle $k_1$ subtracted from the elliptical segment $PCB$: + + \begin{equation}\label{areas_eq} + k = area(seg_{PCB}) - area(k_1) + \end{equation} + + Where the area of the triangle $k_1$ can be found easily using geometric formulae: + + \begin{align} + area(k_1) &= \frac{1}{2} \left( ae - a \cos E \right) \left( \frac{b}{a} a \sin E \right) \\ + &= \frac{ab}{2} \left(e \sin E - \cos E \sin E \right) + \end{align} + + Now we can find the area for the elliptical segment $PCB$ by first finding the circular + segment $POB'$, subtracting the triangle $OB'C$, then applying the fact that an ellipse is + merely a vertical scaling of a circle by the amount $\frac{b}{a}$. + + \begin{align} + area(PCB) &= \frac{b}{a} \left( \frac{a^2 E}{2} - \frac{1}{2} \left( a \cos E \right) + \left( a \sin E \right) \right) \\ + &= \frac{abE}{2} - \frac{ab}{2} \left( \cos E \sin E \right) \\ + &= \frac{ab}{2} \left( E - \cos E \sin E \right) + \end{align} + + By substituting the two areas back into Equation~\ref{areas_eq} we can get the $k$ area + swept out by the spacecraft: + + \begin{equation} + k = \frac{ab}{2} \left( E - e \sin E \right) + \end{equation} + + Which we can then substitute back into the equation for the swept area as a function of + time (Equation~\ref{swept}): + + \begin{equation} + \frac{\Delta t}{T} = \frac{E - e \sin E}{2 \pi} + \end{equation} + + Which is, effectively, Kepler's equation. It is commonly known by a different form: + + \begin{align} + M &= E - e \sin E \\ + &= \sqrt{\frac{\mu}{a^3}} \Delta t + \end{align} + + Where we've defined the mean anomaly as $M$ and used the fact that $T = + \sqrt{\frac{a^3}{\mu}}$. This provides us a useful relationship between Eccentric Anomaly + ($E$) which can be related to spacecraft position, and time, but we still need a useful + algorithm for solving this equation. + + \subsubsection{LaGuerre-Conway Algorithm}\label{laguerre} + For this application, I used an algorithm known as the LaGuerre-Conway algorithm, which was + presented in 1986 as a faster algorithm for directly solving Kepler's equation and has been + in use in many applications since. This algorithm is known for its convergence robustness + and also its speed of convergence when compared to higher order Newton methods. + + This thesis will omit a step-through of the algorithm itself, but the code will be present + in the Appendix. + + \subsection{Non-Linear Problem Optimization} + + Now we can consider the formulation of the problem in a more useful way. For instance, given a + desired final state in position and velocity we can relatively easily determine the initial + state necessary to end up at that desired state over a pre-defined period of time by solving + Kepler's equation. In fact, this is often how impulsive trajectories are calculated since, + other than the impulsive thrusting event itself, the trajectory is entirely natural. + + However, often in trajectory design we want to consider a number of other inputs. For + instance, a low thrust profile, a planetary flyby, the effects of rotating a solar panel on + solar radiation pressure, etc. Once these inputs have been accepted as part of the model, the + system is generally no longer analytically solvable, or, if it is, is too complex to calculate + directly. + + Therefore an approach is needed, in trajectory optimization and many other fields, to optimize + highly non-linear, unpredictable systems such as this. The field that developed to approach + this problem is known as Non-Linear Problem (NLP) Optimization. + + There are, however, two categories of approaches to solving an NLP. The first category, + indirect methods, involve declaring a set of necessary and/or sufficient conditions for declaring + the problem optimal. These conditions then allow the non-linear problem (generally) to be + reformulated as a two point boundary value problem. Solving this boundary value problem can + provide a control law for the optimal path. Indirect approaches for spacecraft trajectory + optimization have given us the Primer Vector Theory. + + The other category is the direct methods. In a direct optimization problem, the cost function + itself is calculated to provide the optimal solution. The problem is usually thought of as a + collection of dynamics and controls. Then these controls can be modified to minimize the cost + function. A number of tools have been developed to optimize NLPs via this direct method in the + general case. For this particular problem, direct approaches were used as the low-thrust + system dynamics adds too much complexity to quickly optimize indirectly and the individual + optimization routines needed to proceed as quickly as possible. + + \subsubsection{Non-Linear Solvers} + For these types of non-linear, constrained problems, a number of tools have been developed + that act as frameworks for applying a large number of different algorithms. This allows for + simple testing of many different algorithms to find what works best for the nuances of the + problem in question. + + One of the most common of these NLP optimizers is SNOPT, which is a proprietary package + written primarily using a number of Fortran libraries by the Systems Optimization Laboratory + at Stanford University. It uses a sparse sequential quadratic programming approach. + + Another common NLP optimization packages (and the one used in this implementation) is the + Interior Point Optimizer or IPOPT. It can be used in much the same way as SNOPT and uses an + Interior Point Linesearch Filter Method and was developed as an open-source project by the + organization COIN-OR under the Eclipse Public License. + + Both of these methods utilize similar approaches to solve general constrained non-linear + problems iteratively. Both of them can make heavy use of derivative Jacobians and Hessians + to improve the convergence speed and both have been ported for use in a number of + programming languages, including in Julia, which was used for this project. + + This is by no means an exhaustive list, as there are a number of other optimization + libraries that utilize a massive number of different algorithms. For the most part, the + libraries that port these are quite modular in the sense that multiple algorithms can be + tested without changing much source code. + + \subsubsection{Linesearch Method} + As mentioned above, this project utilized IPOPT which leveraged an Interior Point Linesearch + method. A linesearch algorithm is one which attempts to find the optimum of a non-linear + problem by first taking an initial guess $x_k$. The algorithm then determines a step + direction (in this case through the use of automatic differentiation to calculate the + derivatives of the non-linear problem) and a step length. The linesearch algorithm then + continues to step the initial guess, now labeled $x_{k+1}$ after the addition of the + ``step'' vector and iterates this process until predefined termination conditions are met. + + In this case, the IPOPT algorithm was used, not as an optimizer, but as a solver. For + reasons that will be explained in the algorithm description in Section~\ref{algorithm} it + was sufficient merely that the non-linear constraints were met, therefore optimization (in + the particular step in which IPOPT was used) was unnecessary. + + \section{Low-Thrust Considerations} \label{low_thrust} + + Thus far, the techniques that have been discussed can be equally useful for both impulsive and + continuous thrust mission profiles. In this section, we'll discuss the intricacies of continuous + low-thrust trajectories in particular. There are many methods for optimizing such profiles and + we'll briefly discuss the difference between a direct and indirect optimization of a low-thrust + trajectory as well as introduce the concept of a control law and the notation used in this + thesis for modelling low-thrust trajectories more simply. + + \subsection{Low-Thrust Control Laws} + + In determining a low-thrust arc, a number of variables must be accounted for and, ideally, + optimized. + + Firstly, we must determine the presence or absence of thrust. Often, this is a question of + preference in the arsenal of the mission designer. Generally speaking, there are points along + an orbit at which thrusting in order to achieve the final orbit are more or less efficient. + For instance, in a classic orbit raising, if increasing the semi-major axis is the only goal, + then thrusting nearer to the periapsis is far more efficient than thrusting near the apoapsis. + For this reason, a mission designer may choose to reduce the thrust or turn it off altogether + during certain segments of the trajectory. + + Secondly, the direction of thrust must also be determined. The methods for determining this + direction varies greatly depending on the particular control law chosen for that mission. + Generally speaking, a control law determines these two parameters: thrust presence and thrust + direction, at each point along the arc. + + This is, of course, also true for impulsive trajectories. However, since the thrust presence + for those trajectories are generally taken to be impulse functions, the control laws can + afford to be much less complicated for a given mission goal, by simply thrusting only at the + moment on the orbit when the transition will be most efficient. For a low-thrust mission, + however, the control law must be continuous rather than discrete and therefore the control law + inherently gains a lot of complexity. + + \subsection{Sims-Flanagan Transcription} + + The major problem with optimizing low thrust paths is that the control law must necessarily be + continuous. Also, since indirect optimization approaches are quite difficult, the problem must + necessarily be reformulated as a discrete one in order to apply a direct approach. Therefore, + this thesis chose to use a model well suited for discretizing low-thrust paths: the + Sims-Flanagan transcription (SFT). + + The SFT is actually quite a simple method for discretizing low-thrust arcs. First the + continuous arc is subdivided into a number ($N$) of individual consistent timesteps of length + $\frac{tof}{N}$. The control thrust is then applied at the center of each of these time steps. + + Using the SFT, it is relatively straightforward to propagate a state (in the context of the + Two-Body Problem) that utilizes a continuous low-thrust control, without the need for + computationally expensive numeric integration algorithms, by simply solving Kepler's equation + (using the LaGuerre-Conway algorithm introduced in Section~\ref{laguerre}) $N$ times. This + greatly reduces the computation complexity, which is particularly useful for cases in which + low-thrust trajectories need to be calculated many millions of times, as is the case in this + thesis. The fidelity of the model can also be easily fine-tuned. By simply increasing the + number of sub-arcs, one can rapidly approach a fidelity equal to a continuous low-thrust + trajectory within the Two-Body Problem, with only linearly-increasing computation time. + + \section{Interplanetary Trajectory Considerations} \label{interplanetary} + Highlight the problems with the 2BP in co-ordinating influences of extra bodies over an + interplanetary journey. + + \subsection{Patched Conics} + Describe the method of patched conics. + + \subsection{Gravity Assist Maneuvers} + Describe how a gravity assist maneuver would work in the framework of patched conics. Also + discuss the advantages of such a maneuver. + + \subsection{Multiple Gravity Assist Techniques} + Discuss the advantages of chaining together multiple gravity assists and highlight the + difficulties in choosing these assists. Here I can mention porkchop plots, Lambert's problem, + etc. Here I can also talk about Hybrid Optimal Control Problems. + + \subsection{Ephemeris Considerations} + I can quickly mention SPICE here and talk a bit about validation. + + % \section{Genetic Algorithms} + % I will probably give only a brief overview of genetic algorithms here. I don't personally know + % that much about them. Then in the following subsections I can discuss the parts that are + % relevant to the specific algorithm that I'm using. + + % \subsection{Decision Vectors} + % Discuss what a decision vector is in the context of an optimization problem. + + % \subsection{Selection and Fitness Evaluation} + % Discuss the costing being used as well as the different types of fitness evaluation that are + % common. Also discuss the concept of generations and ``survival''. + + % \subsubsection{Tournament Selection} + % Dive deeper into the specific selection algorithm being used here. + + % \subsection{Crossover} + % Discuss the concept of crossover and procreation in a genetic algorithm. + + % \subsubsection{Binary Crossover} + % Discuss specific crossover algorithm used here. + + % \subsubsection{Mutation} + % Discuss both the necessity for mutation and the mutation algorithm being used. + + \section{Algorithm Overview} \label{algorithm} + Highlight the algorithm at a high-level. This is likely where flowcharts and diagrams will go to + give a high-level overview. + + \subsection{Trajectory Composition} + Discuss briefly the nomenclature used in defining these trajectories. Currently this isn't + ``baked in'' to the code, so I have some freedom to adopt Englander's notation or use my own + (since my intended use case is a little simpler). + + \subsection{Inner Loop Implementation} + Give a better overview of the inner loop specifically. Probably this section will have a more + in-depth flowchart. + + \subsubsection{LaGuerre-Conway Kepler Solver} + Discuss how the LaGuerre-Conway algorithm is used in the code to provide a fundamental + ``natural trajectory'' between two quantized, but not necessarily close points. Mention + validation. + + \subsubsection{Sims-Flanagan Propagator} + Discuss how this algorithm can then be expanded by using SFT to propagate any number of + low-thrust steps over a specific arc. Mention validation. Here I can also mention the ``Sc'' + object and talk about how those parameters were chosen and effected the propagator. + + \subsubsection{Non-Linear Problem Solver} + Mention the package being used to solve NLPs and how it works, highlighting the trust region + method used and error-handling. Mention validation. + + \subsubsection{Monotonic Basin Hopping} + Outline the MBH algorithm, going into detail at each step. Mention the long-tailed PDF being + used and go into quite a bit of detail. Englander's paper on the MBH algorithm specifically + should be a good guide. Mention validation. + + \subsection{Outer Loop Implementation} + Overview the outer loop. This may require a final flowchart, but might potentially be too + simple to lend itself to one. + + \subsubsection{Inner Loop Calling Function} + The primary reason for including this section is to discuss the error handling. + + \subsubsection{Genetic Algorithm Description} + Similar to the MBH section, there are a lot of implementation details to go over here. Many + will have already been discussed in the background sections above. But I can step through + each of the decisions, similar to Englander's paper on this. + + \section{Results Analysis} \label{results} + Simply highlight that the algorithm was tested on a sample trajectory to Saturn. - \subsubsection{Monotonic Basin Hopping} - Outline the MBH algorithm, going into detail at each step. Mention the long-tailed PDF being - used and go into quite a bit of detail. Englander's paper on the MBH algorithm specifically - should be a good guide. Mention validation. + \subsection{Sample Trajectory to Saturn} + Give an overview of the trajectory that was ultimately chosen. - \subsection{Outer Loop Implementation} - Overview the outer loop. This may require a final flowchart, but might potentially be too - simple to lend itself to one. + \subsubsection{Comparison to Less Optimal Solutions} + I should have a number of elite but less-optimal solutions. Honestly, I may write the + algorithm to keep all of the solutions to provide many points of comparison here. - \subsubsection{Inner Loop Calling Function} - The primary reason for including this section is to discuss the error handling. + \subsubsection{Cost Function Analysis} + Give some real-world context for the mass-use, time-of-flight, etc. - \subsubsection{Genetic Algorithm Description} - Similar to the MBH section, there are a lot of implementation details to go over here. Many - will have already been discussed in the background sections above. But I can step through - each of the decisions, similar to Englander's paper on this. + \subsubsection{Comparison to Impulsive Trajectories} + I may also remove this section. I could do a quick comparison (using porkchop plots) to + similar impulsive trajectories. Honestly, this is a lot of work for very little gain, + though, so probably the first place to chop if needed. - \section{Results Analysis} \label{results} - Simply highlight that the algorithm was tested on a sample trajectory to Saturn. + \section{Conclusion} \label{conclusion} + \subsection{Overview of Results} + Quick re-wording of the previous section in a paragraph or two for reader's convenience. - \subsection{Sample Trajectory to Saturn} - Give an overview of the trajectory that was ultimately chosen. + \subsection{Applications of Algorithm} + Talk a bit about why this work is valuable. Missions that could have benefited, missions that + this enables, etc. - \subsubsection{Comparison to Less Optimal Solutions} - I should have a number of elite but less-optimal solutions. Honestly, I may write the - algorithm to keep all of the solutions to provide many points of comparison here. + \subsection{Recommendations for Future Work} + Recommend future work, obviously. There are a \emph{ton} of opportunities for improvement + including parallelization, cluster computing, etc. - \subsubsection{Cost Function Analysis} - Give some real-world context for the mass-use, time-of-flight, etc. - - \subsubsection{Comparison to Impulsive Trajectories} - I may also remove this section. I could do a quick comparison (using porkchop plots) to - similar impulsive trajectories. Honestly, this is a lot of work for very little gain, - though, so probably the first place to chop if needed. - - \section{Conclusion} \label{conclusion} - \subsection{Overview of Results} - Quick re-wording of the previous section in a paragraph or two for reader's convenience. - - \subsection{Applications of Algorithm} - Talk a bit about why this work is valuable. Missions that could have benefited, missions that - this enables, etc. - - \subsection{Recommendations for Future Work} - Recommend future work, obviously. There are a \emph{ton} of opportunities for improvement - including parallelization, cluster computing, etc. - - % \bibliography{biblio}{} - % \bibliographystyle{plain} + % \bibliography{biblio}{} + % \bibliographystyle{plain} \end{document}