diff --git a/LaTeX/fig/laguerre_plot.png b/LaTeX/fig/laguerre_plot.png new file mode 100644 index 0000000..96961e0 Binary files /dev/null and b/LaTeX/fig/laguerre_plot.png differ diff --git a/LaTeX/fig/spiral_plot.png b/LaTeX/fig/spiral_plot.png new file mode 100644 index 0000000..8671cfb Binary files /dev/null and b/LaTeX/fig/spiral_plot.png differ diff --git a/LaTeX/thesis.bib b/LaTeX/thesis.bib index f2cce9c..38be63d 100644 --- a/LaTeX/thesis.bib +++ b/LaTeX/thesis.bib @@ -8,3 +8,19 @@ year={1986}, publisher={Springer} } + +@software{snow, + author = {Andrew Ning}, + title = {SNOW.jl}, + url = {https://github.com/byuflowlab/SNOW.jl}, + version = {0.2.0}, + date = {2022-02-10}, +} + +@article{RevelsLubinPapamarkou2016, + title = {Forward-Mode Automatic Differentiation in {J}ulia}, + author = {{Revels}, J. and {Lubin}, M. and {Papamarkou}, T.}, + journal = {arXiv:1607.07892 [cs.MS]}, + year = {2016}, + url = {https://arxiv.org/abs/1607.07892} +} diff --git a/LaTeX/thesis.tex b/LaTeX/thesis.tex index d4e8f3d..87dc440 100644 --- a/LaTeX/thesis.tex +++ b/LaTeX/thesis.tex @@ -653,12 +653,11 @@ \begin{figure} \centering - \includegraphics[width=\textwidth]{fig/kepler} + \includegraphics[width=\textwidth]{fig/laguerre_plot} \caption{Example of a natural trajectory propagated via the Laguerre-Conway approach to solving Kepler's Problem} \label{laguerre_plot} \end{figure} - % TODO: Generate an orbit figure for here % TODO: Consider adding a paragraph about the improvements in processor time \subsection{Sims-Flanagan Propagator} @@ -681,20 +680,85 @@ \begin{figure} \centering - \includegraphics[width=\textwidth]{fig/kepler} + \includegraphics[width=\textwidth]{fig/spiral_plot} \caption{An example trajectory showing that classic continuous-thrust orbit shapes, such as this orbit spiral, are easily achievable using a Sims-Flanagan model} \label{sft_plot} \end{figure} - % TODO: Generate an orbit figure for here Figure~\ref{sft_plot} shows that the Sims-Flanagan transcription model can be used - to effectively model these types of orbit trajectories. + to effectively model these types of orbit trajectories. In fact, the Sims-Flanagan + model is capable of modeling nearly any low-thrust trajectory with a sufficiently + high number of $n$ samples. \subsection{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. + + Now that we have the basic building blocks of a continuous-thrust trajectory, we can + leverage one of the many non-linear optimization packages to find solutions near to + a (proposed) trajectory. This trajectory need not be valid. + + For the purposes of discussion in this Section, we will assume that the inner-loop + algorithm starts with just such a ''Mission Guess``, which represents the proposed + trajectory. However, we'll briefly mention what quantities are needed for this + input: + + A Mission Guess object contains: + \begin{singlespacing} + \begin{itemize} + \item The spacecraft and thruster parameters for the mission + \item A launch date + \item A launch $v_\infty$ vector representing excess Earth velocity + \item For each phase of the mission: + \begin{itemize} + \item The planet that the spacecraft will encounter (either flyby or + complete the mission) at the end of the phase + \item The $v_{\infty,out}$ vector representing excess velocity at the + planetary flyby (or launch if phase 1) at the beginning of the phase + \item The $v_{\infty,in}$ vector representing excess velocity at the + planetary flyby (or completion of mission) at the end of the phase + \item The time of flight for the phase + \item The unit-thrust profile in a sun-fixed frame represented by a + series of vectors with each element ranging from 0 to 1. + \end{itemize} + \end{itemize} + \end{singlespacing} + + From this information, as can be seen in Figure~\ref{nlp}, we can formulate the + mission in terms of a non-linear problem. Specifically, the Mission Guess object can + be represented as a vector, $x$, the propagation function as a function $F$, and the + constraints as another function $G$ such that $G(x) = \vec{0}$. + + This is a format that we can apply directly to the IPOPT solver, which Julia (the + programming language used) can utilize via bindings supplied by the SNOW.jl + package\cite{snow}. + + IPOPT also requires the derivatives of both the $F$ and $G$ functions in the + formulation above. Generally speaking, a project designer has two options for + determining derivatives. The first option is to analytically determine the + derivatives, guaranteeing accuracy, but requiring processor time if determined + algorithmically and sometimes simply impossible or mathematically very rigorous to + determine manually. The second option is to numerically derive the derivatives, + using a technique such as finite differencing. This limits the accuracy, but can be + faster than algorithmic symbolic manipulation and doesn't require rigorous manual + derivations. + + However, the Julia language has an excellent interface to a new technique, known as + automatic differentiation\cite{RevelsLubinPapamarkou2016}. Automatic differentiation + takes a slightly different approach to numerical derivation. It takes advantage of + the fact that any algorithmic function, no matter how complicated, can be broken + down into a series of smaller arithmetic functions, down to the level of simple + arithmetic. Since all of these simple arithmetic functions have a known derivative, + we can define a new datatype that carries through the function both the float and a + second number representing the derivative. Then, by applying (to the derivative) the + chain rule for every minute arithmetic function derivative as that arithmetic + function is applied to the main float value, the derivative can be determined, + accurate to the machine precision of the float type being used, with a processing + equivalent of two function calls (this of course depends on the simplicity of the + chained derivatives compared to the function pieces themselves). Generally speaking + this is much faster than the three or more function calls necessary for accurate + finite differencing and removes the need for the designer to tweak the epsilon value + in order to achieve maximum precision. \section{Outer Loop Implementation} Overview the outer loop. This may require a final flowchart, but might potentially be too diff --git a/julia/Manifest.toml b/julia/Manifest.toml index 86b5088..5a9fc7d 100644 --- a/julia/Manifest.toml +++ b/julia/Manifest.toml @@ -378,7 +378,7 @@ uuid = "093fc24a-ae57-5d10-9952-331d41423f4d" version = "1.3.5" [[LinearAlgebra]] -deps = ["Libdl"] +deps = ["Libdl", "libblastrampoline_jll"] uuid = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" [[LogExpFunctions]] @@ -477,6 +477,10 @@ git-tree-sha1 = "ba4a8f683303c9082e84afba96f25af3c7fb2436" uuid = "656ef2d0-ae68-5445-9ca0-591084a874a2" version = "0.3.12+1" +[[OpenBLAS_jll]] +deps = ["Artifacts", "CompilerSupportLibraries_jll", "Libdl"] +uuid = "4536629a-c528-5b80-bd46-f80d51c5b363" + [[OpenSpecFun_jll]] deps = ["Artifacts", "CompilerSupportLibraries_jll", "JLLWrappers", "Libdl", "Pkg"] git-tree-sha1 = "13652491f6856acfd2db29360e1bbcd4565d04f1" @@ -535,7 +539,7 @@ deps = ["InteractiveUtils", "Markdown", "Sockets", "Unicode"] uuid = "3fa0cd96-eef1-5676-8a61-b3b8758bbffb" [[Random]] -deps = ["Serialization"] +deps = ["SHA", "Serialization"] uuid = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" [[Reexport]] @@ -723,6 +727,10 @@ git-tree-sha1 = "9e7a1e8ca60b742e508a315c17eef5211e7fbfd7" uuid = "700de1a5-db45-46bc-99cf-38207098b444" version = "0.2.1" +[[libblastrampoline_jll]] +deps = ["Artifacts", "Libdl", "OpenBLAS_jll"] +uuid = "8e850b90-86db-534c-a0d3-1478176c7d93" + [[nghttp2_jll]] deps = ["Artifacts", "Libdl"] uuid = "8e850ede-7688-5339-a07c-302acd2aaf8d" diff --git a/julia/make_plots.jl b/julia/make_plots.jl new file mode 100644 index 0000000..7bc4409 --- /dev/null +++ b/julia/make_plots.jl @@ -0,0 +1,20 @@ +using Thesis +using PlotlyJS: savefig + +a = ∛(Sun.μ * ( year/2π )^2 ) +start = [ oe_to_xyz([ a, 0., 0., 0., 0., 0.5 ], Sun.μ); 10_000. ] +path = Thesis.prop(zeros(100,3), start, bepi, year)[1] +p = Thesis.plot([path], + colors=["#0FF"], + title="Single LaGuerre-Conway-Propagated Orbit") +savefig(p, "LaTeX/fig/laguerre_plot.png") + +t = 5year +spiral_sc = Sc("test", 1000., 2000., 0.0005, 50, 0.9) +start = [ oe_to_xyz([ 0.7a, 0.05, 0.1, 0., 0., 0.5 ], Sun.μ); 100_000. ] +spiral_thrust = spiral(0.2, 100_000, start, spiral_sc, t) +spiral_path = Thesis.prop(spiral_thrust, start, spiral_sc, t)[1] +p = Thesis.plot([spiral_path], + colors=["#0FF"], + title="Spiral Orbit") +savefig(p, "LaTeX/fig/spiral_plot.png") diff --git a/julia/src/utilities/conversions.jl b/julia/src/utilities/conversions.jl index ee389db..a545cb3 100644 --- a/julia/src/utilities/conversions.jl +++ b/julia/src/utilities/conversions.jl @@ -76,14 +76,14 @@ Output: a random reasonable orbit """ function gen_orbit(T::Float64, mass::Float64, primary::Body=Sun) μ = primary.μ - i = rand(0.0:0.01:0.4999π) + inc = rand(0.0:0.01:0.4999π) θ = rand(0.0:0.01:2π) i = 0 while true i += 1 e = rand(0.0:0.01:0.5) a = ∛(μ * ( T/2π )^2 ) - a*(1-e) < 1.1primary.r || return [ oe_to_xyz([ a, e, i, 0., 0., θ ], μ); mass ] + a*(1-e) < 1.1primary.r || return [ oe_to_xyz([ a, e, inc, 0., 0., θ ], μ); mass ] i < 100 || throw(GenOrbit_Error) end end diff --git a/julia/src/utilities/plotting.jl b/julia/src/utilities/plotting.jl index a0fd403..ea2f6e2 100644 --- a/julia/src/utilities/plotting.jl +++ b/julia/src/utilities/plotting.jl @@ -44,14 +44,14 @@ Generates a layout that works for most plots. Requires a title and a limit. """ function standard_layout(limit::Float64, title::AbstractString) limit *= 1.1 - Layout(title=title, - # width=1400, - # height=800, + Layout(title=attr(font=attr(color="rgb(250,250,250)"), yanchor="top", y=0.95, text=title), + textfont = attr(color="rgb(250,250,250)"), + margin=attr(b=20, t=0, l=0, r=0), paper_bgcolor="rgba(5,10,40,1.0)", - plot_bgcolor="rgba(100,100,100,0.01)", - scene = attr(xaxis = attr(autorange = false,range=[-limit,limit]), - yaxis = attr(autorange = false,range=[-limit,limit]), - zaxis = attr(autorange = false,range=[-limit,limit]), + plot_bgcolor="rgba(200,200,200,0.01)", + scene = attr(xaxis = attr(autorange = false,range=[-limit,limit],color="rgb(255,255,255)"), + yaxis = attr(autorange = false,range=[-limit,limit],color="rgb(255,255,255)"), + zaxis = attr(autorange = false,range=[-limit,limit],color="rgb(255,255,255)"), aspectratio=attr(x=1,y=1,z=1), aspectmode="manual")) end