Trajectory design for a solar-sail mission to asteroid 2016 HO3

This paper proposes the use of solar-sail technology currently under development at NASA Langley Research Center for a CubeSat rendezvous mission with asteroid 2016 HO3, a quasi-satellite of Earth. Time-optimal trajectories are sought for within a 2022–2023 launch window, starting from an assumed launcher ejection condition in the Earth-Moon system. The optimal control problem is solved through a particular implementation of a direct pseudo-spectral method for which initial guesses are generated through a relatively simple and straightforward genetic algorithm search on the optimal launch date and sail attitude. The results show that the trajectories take 2.16–4.21 years to complete, depending on the assumed solar-sail reflectance model and solar-sail technology. To assess the performance of solar-sail propulsion for this mission, the trajectory is also designed assuming the use of solar electric propulsion. The resulting fuel-optimal trajectories take longer to complete than the solar-sail trajectories and require a propellant consumption that exceeds the expected propellant capacity onboard the CubeSat. This comparison demonstrates the superior performance of solar-sail technology for this mission.


Introduction
On 27 April 2016, the Pan-STARRS 1 asteroid survey telescope on Haleakalā, Hawaii, USA, detected a remarkable asteroid: 2016 HO 3 x . Its orbit is extremely similar to that of Earth, to the extent that asteroid 2016 HO 3 appears to orbit around our planet. The asteroid is therefore considered a near-Earth companion or a quasisatellite of Earth and is expected to accompany the Earth for hundreds of years. This unique characteristic, together with the already significant increase in interest in small-body research over recent years, makes 2016 HO 3 an interesting mission target.
Even though the orbit of 2016 HO 3 is very similar to that of Earth, its phasing and the 7.8 deg inclination with respect to the ecliptic makes it a difficult target x Jet Propulsion Laboratory. Small asteroid is Earth's constant companion. Available at https://www.jpl.nasa.gov/news/news.php? feature=6537. [ Accessed 25 February 2019] m.j.heiligers@tudelft.nl to reach. Low-thrust propulsion, either in the form of solar electric propulsion (SEP) or solar sailing [1,2], has been proven to enable high-energy missions. Examples for the use of SEP include JAXA's asteroid sample return mission Hayabusa [3], NASA's Dawn mission that visited the two largest bodies in the asteroid belt [4], and ESA's BepiColombo mission to Mercury [5]. Examples of proposed high-energy solar-sail missions include NASA's NEA Scout mission [6], as well as a range of theoretical mission concepts such as the Solar Polar Orbiter [7], the Geostorm mission concept [8,9], the Interstellar Heliopause Probe [10], asteroid rendezvous missions [11,12], and, more generally, a wide range of highly non-Keplerian orbits for novel space applications [1,[13][14][15][16][17][18].
This paper investigates the use of solar-sail propulsion to rendezvous with asteroid 2016 HO 3 . In particular, the solar-sail technology currently under development at NASA Langley Research Center is considered [19]. The assumed mission configuration is that of a CubeSat platform and a launch as secondary payload (e.g., onboard one of the Exploration Missions of the SLS launch vehicle) within a wide 2022-2023 launch window.
The assumption of a CubeSat-sized platform drove the choice for solar-sail propulsion as mass and dimension constraints limit the available space for SEP propellant as well as solar arrays to provide power to the SEP system. To confirm this choice, the trajectory will not only be designed for the use of solar-sail propulsion, but also for the use of solar electric propulsion.
The objective of the work in this paper is to find time-optimal solar-sail-or alternatively fuel-optimal SEP-trajectories from an assumed launcher ejection condition in the Earth-Moon circular restricted threebody problem (CR3BP) to 2016 HO 3 . An initial, ballistic trajectory up to the sphere of influence of the Earth is assumed to allow spacecraft testing and verification. Once at the sphere of influence, the lowthrust propulsion system is activated and the modelling of the trajectory continues in the Sun-Earth CR3BP. Time-and fuel-optimal trajectories are found through the application of a specific direct pseudospectral optimal control solver, PSOPT [20]. Initial guesses for the optimal control solver are obtained through a relatively straightforward genetic algorithm routine that finds the optimal launch date and constant direction of the low-thrust acceleration vector with respect to the direction of sunlight to minimize the miss-distance and miss-velocity at the asteroid. By subsequently feeding this initial guess to the optimal control solver, where the control is allowed to vary over time, these errors are overcome and the time of flight-or alternatively the fuel consumption-is minimized.

Orbit of 2016 HO 3
The orbital elements of asteroid 2016 HO 3 are provided in Table 1. With a semi-major axis very close to that of Earth and only a small eccentricity, its orbit resembles that of Earth (only inclined at 7.8 deg with respect to the ecliptic). Graphical representations of the asteroid's motion in the timeframe 1960-2020 appear in Fig. 1, where the Earth is assumed to be on a Keplerian, 1 AU (astronomical unit) circular orbit. Figure 1(a) shows the orbit in an inertial frame, A(X I , Y I , Z I ), where the X Iaxis points towards the vernal equinox, the Z I -axis is oriented perpendicular to the ecliptic, and the Y I -axis  LightSail A) and LightSail 2 solar sails [22]. TRAC booms have been problematic for larger solar sails due to their high coefficient of thermal expansion (CTE), very low torsional stiffness, and low deployed precision [23,24]. An improved version of the composite-based EDU solar sail-the Advanced Composites-Based Solar Sail System (ACS3)-is now under development by NASA LaRC and NASA Ames Research Center for a low Earth orbit (LEO) solar-sail technology risk reduction mission in the 2021 timeframe. The 12U ACS3 flight experiment is intended as a technology development pathfinder for a future, larger composite-based small satellite solar-sail system suitable for 12U to 27U CubeSat class spacecraft. For purposes of this study, a lightness number range bounding the anticipated solar-sail performance of a notional 12U-27U CubeSat-class spacecraft using the ACS3 solar-sail technology is assumed.

Mission assumptions
To design the trajectory from launch ejection to 2016 HO 3 , a set of assumptions are made: • Launch is assumed to take place in 2022-2023. • The trajectory is assumed to start from the following launcher ejection conditions in a synodic Earth-Moon frame B EM (x EM , y EM , z EM ) with r 0 the initial position vector (km) and v 0 the initial velocity vector (km/s):  (7)). • It is assumed that the solar-sail system can be replaced by a solar electric propulsion system (and power system) with a performance that is based on the following assumptions: An initial spacecraft mass of 14-21 kg. A specific impulse of 1600 s [25]. A maximum thrust of 0.9 mN [25]. A maximum propellant mass capacity of 1.5 kg [25]. A maximum thruster operation duration of 20,000 hours (2.3 years) [25].

Dynamics
As outlined in the previous section, two different sets of dynamical frameworks are adopted to design the trajectory to 2016 HO 3 : the Earth-Moon CR3BP from the initial condition in Eq. (1) up to the Earth's SOI (hereafter referred to as the Earth-Moon ballistic phase) and the Sun-Earth CR3BP from the Earth's SOI up to rendezvous with the asteroid (hereafter referred to as the interplanetary phase). The dynamics in both phases will be detailed in the following subsections, starting with the interplanetary phase.

Interplanetary phase
Considering the relatively close proximity of the asteroid to the Earth and its semi-bounded motion around the Earth as shown in Fig. 1(b), the trajectory to 2016 HO 3 is designed in the framework of the Sun-Earth CR3BP.
In the CR3BP, the motion of an infinitely small mass, m (i.e., the spacecraft), is described under the influence of the gravitational attraction of two much larger primary masses, m 1 (here, the Sun) and m 2 (here, the Earth).
The gravitational influence of the small mass on the primary masses is neglected and the primary masses are assumed to move in circular orbits about their common center-of-mass. The reference frame employed to define the spacecraft's dynamics in the Sun-Earth CR3BP is that of Fig. 1 New units are introduced: the sum of the two primary masses is taken as the unit of mass, i.e., m 1 + m 2 = 1. Then, with the mass ratio μ = m 2 /(m 1 + m 2 ), the masses of these primary bodies become m 1 = 1 − μ and m 2 = μ. As unit of length, the distance between the primary bodies is selected, and 1/ω is chosen as the unit of time, yielding ω = 1. Then, one revolution of the reference frame (i.e., one year for the Sun-Earth CR3BP) is represented by 2π. In this framework, the motion of a low-thrust propelled spacecraft is described byr The terms on the left-hand side of Eq. (2) are the kinematic, coriolis, and centripetal accelerations, respectively, while the terms on the right-hand side are the low-thrust acceleration and the gravitational acceleration exerted by the primary masses. In frame given by where r 1 and r 2 are the magnitudes of the vectors respectively; see Fig. 2. Finally, the term a is the low-thrust acceleration vector which in this paper is either provided by a solar sail, a = a s , or a solar electric propulsion system, a = a T . Its definition therefore depends on the type of low-thrust propulsion system employed which will be discussed separately in the following two subsections.

Solar-sail propulsion
To model the solar-sail acceleration, this paper will consider both an ideal and an optical solar-sail reflectance model. Considering both these models will allow to compute not only the theoretically fastest trajectory possible (for the ideal model) but also a more realistic trajectory (for the optical model). While the ideal model assumes the sail to be a perfect, specular reflector, the optical model also includes the effects of absorption, diffuse reflection, and thermal emission. Though different in performance, both solar-sail models can be captured in the mathematical definition provided below.
The solar-sail acceleration vector can be decomposed into a component normal to the sail, a nn , and a component tangential to the sail, a tt ; see Fig. 3: The normal to the sail,n, can be defined through two angles, the solar-sail pitch and clock angles, that define the solar sail's orientation with respect to the direction of sunlight; see Fig. 4. For this, a new reference frame S(r 1 ,θ 1 ,φ 1 ) is defined where ther 1 -unit vector is directed along the Sun-sail line (see also below Eq. (3)) and the two remaining axes are defined as in Fig. 4 (left schematic). The pitch angle, α, is then defined as the angle between the normal vector,n, andr 1 , and the clock angle, δ, is defined as the angle between the projection ofn onto the (θ 1 ,φ 1 )-plane andφ 1 . This gives the following definition of the normal vector with respect to frame S(r 1 ,θ 1 ,φ 1 ): Note that, due to the solar sail's inability to generate an acceleration component in the direction of the Sun, the normal vector always points away from the Sun. This can be reflected through appropriate bounds on the pitch and clock angles: The magnitudes of the solar-sail acceleration components along the normal and tangential directions in Eq. (4) are given by [1]: In Eq. (7), β is the solar-sail lightness number,r is the reflectivity coefficient that indicates the fraction of reflected photons, and s indicates the fraction of photons that are specularly reflected, while the term (1−s) indicates the fraction of photons that are diffusely reflected; B f and B b are the non-Lambertian coefficients of the front (subscript "f") and back (subscript "b") of the sail, and ε f and ε b are the corresponding emissivity coefficients. Values for these optical coefficients for both an ideal sail and an optical sail model appear in Table 2. The optical sail coefficients have recently been obtained for NASA's proposed Near Earth Asteroid (NEA) Scout mission [27]. Finally, note that, by substituting the values for the ideal sail model into Eq. (7), only an acceleration component normal to the sail remains, i.e., a t = 0 and therefore a s = a nn andm =n. From Eq. (7), the resulting acceleration direction,m, can be The direction ofm with respect to frame S(r 1 ,θ 1 ,φ 1 ) can then be defined aŝ Through a transformation matrix, this normal vector can be transformed to the synodic Sun-Earth frame B SE (x SE , y SE , z SE ) for use in Eq. (4):

Solar electric propulsion
In the case of employing solar electric propulsion, the low-thrust acceleration vector in Eq.
(2) is defined as where T = [T x,SE T y,SE T z,SE ] T is the Cartesian SEP thrust vector and m is the spacecraft mass. The SEP thrust vector direction can be defined in a similar way as the solar-sail normal vector using the pitch and clock angles, α T and δ T , respectively: with T the SEP thrust magnitude. Note that this time no restrictions need to be imposed on the pitch angle: Finally, due to the consumption of propellant, the spacecraft mass decreases over time according tȯ with I sp = 1600 s the assumed SEP thruster's specific impulse; see Section 4. g 0 is the Earth's standard free fall constant. The differential equation in Eq. (15) needs to be integrated simultaneously with the spacecraft dynamics in Eq. (2).

Ballistic Earth-Moon phase
As highlighted in the mission assumptions section, the first phase of the trajectory is assumed to be ballistic (no use of the solar sail or SEP thruster) and is modelled in the Earth-Moon CR3BP. The dynamics are then as defined in Eq. (2), only now in the B EM (x EM , y EM , z EM ) frame with the Earth and Moon as primaries (μ = 0.01215) and a = 0. When integrating the dynamics forward from Eq. (1) up to the sphere of influence of the Earth (at a distance of 1,496,513 km), the trajectory in Fig. 5 is obtained. It takes the spacecraft 9 days to reach the sphere of influence. Transforming the end of the trajectory to the inertial frame A(X I , Y I , Z I ) and the synodic Sun-Earth frame B SE (x SE , y SE , z SE ), results in the conditions as shown in Fig. 6. In the inertial frame, the conditions at the SOI oscillate around the orbit of the Earth and complete one revolution in one year, whereas in the synodic Sun-Earth frame, these conditions conduct one revolution per synodic lunar month.

Optimal control problem
Depending on the low-thrust propulsion system employed, the objective in this study is to either minimize the time of flight in the interplanetary part of the trajectory (for the solar-sail configuration) or the propellant consumption (for the SEP configuration). The objective, J, can thus be defined as where t 0 and t f are the initial and final time in the interplanetary phase, respectively, and m f is the final spacecraft mass. The goal then is to find the states, x(t), and controls, u(t), that minimize Eq. (16) and satisfy the dynamics in Eq.
(2) as well as a set of boundary and path constraints. The states are the position and velocity vectors in the CR3BP. For the SEP configuration, the spacecraft mass is added: where the initial state, where the last row in the vectors only applies to the SEP case. The controls are also defined differently for the two low-thrust propulsion configurations: where the following bounds are imposed: In Eq. (20), T max = 0.9 mN is the assumed maximum thrust magnitude; see Section 4. In addition, a path constraint is included to avoid close Earth approaches: r 2 (t) 800,000 km (21) Finally, bounds on the initial and final time need to be specified to ensure a launch in the assumed 2022-2023 launch window and to limit the search space on the final time: 1 January 2022 t 0 31 December 2023 31 December 2023 t f 1 January 2028 (22) Note that the time in the actual implementation of the optimal control problem is defined in non-dimensional units after 1 January 2022, i.e., 1 January 2022 is represented by t = 0, 1 January 2023 is represented by t = 2π, and so on. The optimal control problem defined in Eqs. (16)-(22) is solved with a particular implementation of a direct pseudospectral method in C ++ , PSOPT [20]. PSOPT is an open source tool developed by Victor M. Becerra of the University of Reading, UK. It can use both Legendre and Chebyshev polynomials to approximate and interpolate the dependent variables at the nodes. However, in this work, only the Legendre pseudospectral method is used and PSOPT is interfaced to the NLP solver IPOPT (Interior Point OPTimizer), an open source C ++ implementation of an interior point method for large scale problems [28]. Furthermore, a consecutive mesh refinement of [50, 75, 100] nodes is applied, a convergence tolerance of 10 −6 is used, and a maximum number of iterations per mesh refinement of 1000 is enforced.

Initial guess
In order to initiate the optimization process, PSOPT requires an initial guess of the states, controls, and time, which are constructed through the following approach: • First, the launch date, t L , and the pitch and clock angles, (α, δ) or (α T , δ T ), are fixed. Furthermore, for the solar-sail case, the sail lightness number, β, is fixed whereas for the SEP case the thrust magnitude is set to its maximum value, T max . • Subsequently, the initial condition in Eq. (1) is forward integrated from t L up to the Earth's sphere of influence to construct the ballistic Earth-Moon phase and the end of the trajectory is transformed to the B SE (x SE , y SE , z SE ) frame. • Then, the integration is continued to construct the interplanetary phase, where the low-thrust acceleration is defined by either (α, δ, β) or (α T , δ T , T max ). The integration is truncated after 5 years.
In addition, to aid the trajectory in increasing its inclination to that of the asteroid, a rudimentary out-of-plane steering law is adopted where the out-of-plane component of the acceleration takes the sign of the y SE -coordinate. For example, for the solar-sail case: and similar for the thrust vector, T for the SEP case. • Subsequently, at each time step in the propagated trajectory, t, that occurs after 3 years of flight, the dimensionless error in distance, Δr, and error in velocity, ΔV , between the spacecraft's state-vector and that of the asteroid is computed. • Finally, the trajectory is truncated at the point where the sum of these errors, Δr +ΔV , is minimal. Note that, in dimensionless units, a position error of 5000 km and a velocity error of 1 m/s are of the same order of magnitude. For a given performance of the solar-sail or SEP configuration, i.e., for a given value for β or T max the trajectory is fully defined by the following set of three parameters: which define the objective J(p) = Δr + ΔV (25) To find the values for the parameters in Eq. (24) that minimize the objective in Eq. (25), a genetic algorithm is employed. In particular, the MATLAB function ga.m is used with 1000 individuals and a maximum of 50 generations while enforcing the following bounds on the optimization parameters: 31-12-2023 1 2 π π , solar sail [31-12-2023 π π], SEP Finally, note that, to account for the inherent randomness of the genetic algorithm approach, each simulation case is run for five different seed values.

Solar-sail initial guesses
Solar-sail initial guesses are generated for both the ideal and optical sail models and for four different lightness numbers, β = [0.025, 0.03, 0.035, 0.04], to cover the assumed lightness number range; see Section 4. For each lightness number, the genetic algorithm is run five times (for the five different seed values), resulting in a total of 20 runs per sail model. The best trajectory for each lightness number, i.e., the trajectory with the smallest objective function value among the five runs, appears in Fig. 7 (for an ideal sail model) with numerical values in Table 3 (for both sail models). From these results, it can be deduced that, the larger the lightness number, the smaller the objective function value (the 8th column in Table 3), which indicates that the rendezvous conditions at the asteroid can be more easily met with better solar-sail technology. Furthermore, due to the limited number of controls (especially the constant pitch and clock angles throughout the trajectory), a mismatch in position and significant mismatch in velocity still exist between the spacecraft and asteroid at the end of the trajectory. These errors can be overcome by the  optimizer, as will be shown later.

SEP initial guesses
For the SEP case, initial guesses are generated for the extremes of the range in initial spacecraft mass: m 0 = 14 kg and m 0 = 21 kg; see Section 4. Similar to the initial guesses for the solar-sail case, five runs (for the five different seed values) are conducted for each initial mass. The best trajectories for each value for m 0 appear in Fig. 8 with numerical values in Table 4. From these results, it can be observed that, contrary to the solarsail case, out-of-plane thrusting is barely exploited: the pitch and clock angles are close to 90 deg, indicating that thrusting takes place entirely in the ecliptic plane and along the y SE -axis. As Fig. 8 and Table 4 clearly show, while this leads to relatively small errors on the position, large errors on the velocity remain. Furthermore, due to the assumption of continuous thrusting at the maximum thrust magnitude, the propellant consumption is large: 6.93 kg and 7.47 kg for m 0 = 14 kg and m 0 = 21 kg, respectively. This greatly exceeds the expected propellant budget of 1.5 kg; see Section 4.

Results
Due to the assumed wide launch window (2022-2023), the problem at hand contains many local minima. This already became apparent in the results for the initial guess in the previous section. It can therefore be expected that PSOPT will converge to a different local  minimum when provided with a different initial guess.
To therefore best detect the global minimum, PSOPT is run for each generated initial guess, i.e., not only for the best initial guesses that appear in Table 3, but for all 40 (solar-sail case) and 10 (SEP case) generated initial guesses. The results of the runs that generated the best trajectory in terms of objective function value (i.e., the total time of flight) are presented in this section.

Solar-sail optimal results
The results for each considered value for the solar-sail lightness number appear in Table 5 for both the ideal and optical sail models with details for a subset of the results in Fig. 9 and Fig. 10 for β = 0.025 and β = 0.04, respectively. Note that the clock angle profiles in Fig. 9 and Fig. 10 may appear erratic, but a clock angle switch from −π to π does not require an actual physical change in the solar-sail attitude.
From the results in Table 5, Fig. 9, and Fig. 10, the following observations can be made: • PSOPT is able to overcome the discontinuities in the state vector of the initial guess at the end of the trajectory between the sailcraft and the asteroid, while decreasing the time of flight compared to the initial guess by (on average) 1.4 and 0.7 years for the ideal and optical sail models, respectively. • The best initial guess as in Table 3 does not necessarily lead to the best optimized result. In fact, only for an ideal sail model and a lightness number of β = 0.03 did the best initial guess provide the best optimized result. • The difference in launch date between the initial guess and the corresponding optimized trajectory is, on average, only 2 days. This implies that the reduction in time of flight is only due to an advancement of the arrival date, not due to a    Table 3 and Table 5.
The underlying reason can be found in the fact that the launch conditions change much more rapidly over time than the arrival conditions (see above

SEP optimal results
Following the same approach as for the solar-sail optimal results, minimum-propellant trajectories can be generated for the nominal SEP case, which appear in the first two rows of Table 6 with details for m 0 = 14 kg in Fig. 11. The table shows that two additional cases have been considered with larger specific impulses and maximum thrust magnitudes. These results will be discussed later. From the results in the first two rows of Table 6 and   • Similar conclusions as for the solar-sail case can be drawn regarding the fact that the best initial guess as in Table 4 does not necessarily lead to the best optimized result and that the launch date does not change much (if at all) with respect to the initial guess. From the observations listed above, it can be concluded that solar electric propulsion (under the assumptions of Section 4) does not seem a viable propulsion method for this mission from a propellant consumption point of view, thruster operating time and flight time. To further support this conclusion, additional simulations have been conducted for an even better performing SEP system (see the results in the last two rows in Table 6). Here, the specific impulse has been increased to I sp = 2100 s and the maximum thrust magnitude to T max = 1.15 mN, the maximum value as specified in Ref. [25]. The results in Table 6 show that these improvements lower the propellant consumption and the thruster operation duration. However, the propellant consumption is still larger than the expected budget of 1.5 kg and the thruster operation duration is still too long for m 0 = 21 kg. Finally, while the times of flight are also reduced, they are still longer than those for the solar-sail trajectories for β > 0.025 (ideal sail) and β > 0.03 (optical sail).

Conclusions
This paper has presented minimum-time solar-sail trajectories and minimum-propellant solar electric propulsion (SEP) trajectories for a CubeSat mission to asteroid 2016 HO 3 within a set of mission assumptions. In particular, the trajectory is assumed to start from a fixed state in the Earth-Moon system, followed by a ballistic arc up to the Earth's sphere of influence, after which the low-thrust propulsion system is activated to rendezvous with 2016 HO 3 . For the solar-sail configuration, trajectories take 2.16-3.51 years for lightness numbers in the range 0025-0.04 (with shorter flight times for larger lightness numbers) if the solar sail acts as a perfectly reflecting mirror. If slight optical imperfections are included in the solar-sail reflectance model, the time of flight increases by approximately 20 percent to 2.5-4.21 years. For the assumed nominal SEP thruster performance, times of flight of 4.38 years and 4.62 years are obtained for spacecraft initial masses of 14 kg and 21 kg, respectively, which is longer than any solar-sail case considered. Furthermore, the trajectories require propellant consumptions of 3.60 kg and 5.90 kg, which exceed the expected available onboard propellant mass capacity of 1.5 kg. Finally, the required thruster operation duration for an initial mass of 14 kg exceeds the expected maximum operating duration of 20,000 hours. When pushing the SEP technology to its maximum (assuming a larger specific impulse and maximum thrust magnitude), smaller propellant consumptions are obtained (though still beyond the available 1.5 kg) as well as shorter times of flight, though still not as fast as the solar-sail trajectories with lightness number larger than 0.025 (ideal sail) or 0.03 (optical sail). From these analyses the superior capabilities of solar-sail technology for this mission seem to be evident.