Project Didymos
Earth to NEO 65803 Preliminary Analysis
Porkchop Plot: 1000-Day Launch Window Iteration
Objective
To find the optimal transfer time to reach 65803 Didymos, a Near-Earth Object (NEO), in a 1000 day interval.
The provided orbital elements of Earth and Didymos at some reference time (tref) are:
Orbital Elements {a, e, i, Ω, ω, θ(tref)}
- 1 AU = 148.6 * 106 km
- 1 day = 86400s
- μ = 132712 * 106 km3/s2
65803 Didymos (Credit: Wikipedia)
Original Positions & OrbitAssumptions
- Zero-sphere of influence at both Earth and Didymos
- Restricted two body equations are valid (Justified by the fact that MSun >> MDidymos and MSun >> MEarth)
- No perturbations (Negligible drag, J2 perturbations, 3rd bodies, etc)
- No orbital decay (The only Keplerian Element that changes over time is true anomaly)
Why MATLAB?
MATLAB provides the ode78 suite, which is essential for high-order propagation required in orbital mechanics. Implementing a comparable variable time-step ODE solver in a lower-level language is a significant undertaking on its own.
Furthermore, the native graphing and data visualization tools in MATLAB are significantly more streamlined than system-level alternatives like C++. This allows for rapid iteration of complex Porkchop plots and 3D trajectory renders, which were critical for the preliminary analysis phase of this project.
Restricted 2-Body Problem (R2BP)
Technically, every body in the universe is subjected to the gravitational pull of a near infinite amount of bodies. However, being technically correct, is of little to no use to anyone, if the difference between the technically right answer and an approximation is within reason (or possibly negligible in certain applications), especially since perturbations exist (even when the near infinite amount of bodies are accounted for). The technically correct way to analyze any orbit is to consider N-bodies and create a differential equation using Newton’s Law of Gravitation for each body. The N-body problem (for N>1) does not elicit an analytical solution, and is computationally expensive to compute. (Note that unless you consider every body that currently exists in the universe - which is not possible since humans have not identified every solar system, let alone star, let alone planet, etc - the N-body problem will not be the technically correct solution).
The next logical approximation (after ruling out 8-bodies) is to consider just 2-bodies. As previously stated the N-body problem doesn’t elicit an analytical solution, which includes the 2-body case. So, then the next logical approximation is to consider a restricted case of the 2-body problem, where the orbiting mass is so much smaller than the central mass that it can be reasoned that the orbiting mass does not influence the central mass. This is the R2BP, which does elicit a semi-analytical (Kepler’s Problem).
Situation we are applying R2BP logic to(Credit: Wiki)
Keplerian Element Pic(Credit: Sun, He & Bouman, Katherine & Tiede, Paul & Wang, et al (2022))
Orbital Elements
In a R2BP orbit, there are 6 degrees of freedom that are needed to define the state. The two most common ways to express a state are with a Cartesian three-dimensional position and velocity vector or with Keplerian elements.
Keplerian elements are defined in a perifocal frame and describe the orbit more intuitively. The Keplerian elements used in this report are {a, e, i, Ω, ω, ν}. Compared to a position and velocity vector, which change constantly throughout the orbit, the only Keplerian element that changes with time (assuming no perturbation like secular drift or oscillations), is true anomaly (ν). Converting between a common frame (like earth-centered-inertial) to the perifocal frame is a trivial calculation with the use of the rotation matrix. Furthermore, converting between Keplerian elements and Cartesian position and velocity vectors is trivial.
Kepler’s Problem
The main advantage of the R2BP compared to the N-body problem is that it has a semi-analytical solution. In other words, given a known state at an initial time, we are able to find the state after a given amount of time has elapsed. Rather, Kepler’s problem, which is this exact initial value problem, provides us with the semi-analytical solution. Kepler’s Problem relies on making use of a fictitious mean anomaly and eccentric anomaly in order to solve this problem. The reason it is semi-analytical is because it relies on solving a transcendental equation (M=E-esin(E)), which necessitates the use of numerical methods/techniques (like Newton-Raphons).
Having this semi-analytical solution is advantageous because it significantly decreases computational expenses. In other words, instead of needing to propagate the orbit using a time-step method (which scales with the change in time), we are instead able to simply solve a transcendental equation and do some algebra to get a more accurate result (time-step propagation is prone to inaccuracies that tend to increase as change in time increases).
Kep Problem Reference Circle
Position Boundary Value ProblemLambert Solver / Gauss Problem
In mission analysis, it is of utmost importance to calculate the expense, or cost, or a transfer/maneuver. This expense is often referred to as the ΔV of the transfer/maneuver, since ΔV is an independent metric of the geometry of a spacecraft and can easily be calculated from the Rocket Equation for any given spacecraft. Very intuitively, ΔV is literally the change in velocity required. For orbital transfers ΔVtot= |V|to leave your current orbit and enter the transfer orbit + |V|to leave the transfer orbit and get into the wanted final orbit. So, it is of utmost importance to know the initial and final velocities that the spacecraft needs. In other words, it is of utmost importance to be able to calculate the initial and final velocity vector given the initial position and a change in time.
The Gauss/Lambert Problem is this exact boundary value problem. Similar to Kepler’s Problem it also has a semi-analytical solution. A universal Lambert solver (“lambert_universal.m”) is used in this report. The implementation details can be found in the code, since the implementation is commonly accepted, and not very useful for this report. Although, it should be noted that issues can arise if the position vectors are linearly dependent.
Takeaways
The goal of this report is to perform preliminary mission analysis for a spacecraft transfer from Earth to Didymos. In order to do this we must be able to efficiently calculate certain things. First, we need to be able to convert between a Cartesian-defined state to a Keplerian-defined state and vice versa. This is accomplished analytically with trivial calculations in “kep2rv.m” and “rv2kep.m”. Additionally, we need to be able to calculate the position of an orbiting mass about the Sun after any given time past the reference time. This is accomplished by solving Kepler’s Problem (“propogate_elliptical_ta.m” or “propogate_elliptical.m”). Furthermore, we need to be able to calculate the initial and final velocities given two position vectors and a time-of-flight (tof). This is accomplished using a Lambert solver (“lambert_universal.m”).
With these scripts, it is possible to perform mission analysis given a departure and arrival time, fairly trivially. It is also possible to perform mission analysis given departure and arrival time range, which is best visualized using a Porkchop Plot. Furthermore, it is useful to visualize what these transfers will look like, which can be done by propagating the orbits of Earth and Didymos about the Sun, along with the transfer trajectory. This can be accomplished using ode78, a simplistic rhs_2body.m script and plotting commands.
All the code used in this analysis is used to accomplish the previous paragraph, however there are additional helper functions/scripts which decompose the process, allowing for increased readability, functionality, and a reduction in the amount of repeated code.
Porkchop Plot (Prograde) - [tref, tref + 1000 days]
Specific Time Interval
Analysis Window:
- t₀ ≥ t_ref
- t_f ≤ t_ref + 1000 days
- [t_ref, t_ref + 1000 days]
Mission Results
- Departure (t₀): t_ref + 600 days
- Arrival (t_f): t_ref + 1000 days
- TOF: 400 days
- Local Optimal ΔV 7.86 km/s
Constraints & Parameters:
- Any TOF < 100 days is ignored for feasibility.
- Data calculated at 5-day intervals for optimization.
Optimal Earth to Didymos Transfer Trajectory