Capturing small asteroids into a Sun–Earth Lagrangian point

63rd International Astronautical Congress, Naples, Italy. Copyright 2012 by the International Astronautical Federation. All rights reserved. IAC-12. C1.6.7 CAPTURING SMALL ASTEROIDS INTO SUN-EARTH LAGRANGIAN POINTS FOR MINING PURPOSES Neus Lladó∗ , Yuan Ren† , Josep J. Masdemont‡ , Gerard Gómez§ Abstract The aim of this paper is to study the capture of small Near Earth Objects (NEOs) into the Sun-Earth L2 using low-thrust propulsion for mining or science purposes. As it is well known, the vicinity of these points is inside a net of dynamical channels suitable for the transport in the Earth-Moon neighborhood, so different final destinations from here could be easily considered. Asteroids with very small mass and not representing a potential hazard are analyzed. An initial pruning of asteroids is made, considering NEOs with stellar magnitude bigger than 28, which are the smallest available, and NEOs close to the Earth orbit with semi-major axis between 0.85-1.15. Due to the difficult determination of their physical properties, two methods to estimate the asteroid masses are conducted. A procedure to find the low-thrust optimization trajectories has been implemented. The initial seed is obtained integrating forward the equations of motion plus its conjugated equations expressed in cartesian coordinates and applying the Pontryagins maximum principle to obtain the optimal control with a switching function for the thrust. To refine the trajectory a 4 order Runge-Kutta shooting method has been used. The objective function in this study is the fuel consumption. Finally, the capable asteroids to get captured by a low-thrust engine have been listed indicating the main parameters. 1 Introduction est NEAs closest to the Earth. The objective is to find a feasible trajectory with a technology already demonstrated, i.e. Variable Specific Impulse Magnetoplasma Rocket (VASIMR)3 engines. The mass of the asteroids is a parameter which is very difficult to know. However, this parameter is essential for the trajectory computation, so the first section discusses different solutions to estimate this value. Next section is centered to get a rough departure time and time of flight with a global optimization method, propagating with a two-body model forward and backward the trajectory to a mid-point. The third section describes the core of the paper, the trajectory optimization applying the 4th order Runge-Kutta shooting method. It has been noticed that it is very sensible to converge, thus an accurate initial guess for the optimization problem is needed. The last section lists the results of the trajectory optimization for the pruned asteroids of section 2. Asteroid mining1 will play a key role in providing the future resources for the exploration of the Solar System. A rough spectral taxonomy of asteroid types separates them in three types: C-type (carbonaceous), S-type (stony) and M-type (metallic). Type C asteroids comprise more than 70 % of all asteroids. Eventhough Near Earth Asteroids (NEAs) are potentially the most hazardous objects in space, they are the objects that could be easier to exploit for their raw materials. The current paper targets the transfer of asteroids to the Sun-Earth L2 lagrangian point geometrically defined. Another interesting approach could be to insert them into the stable manifold of a libration point orbit.2 However, due to the high number of possibilities and the number of asteroids considered in the study this has been left for further research. Capturing an asteroid near the Earth would make easier to mine it, as well as to exploit it in terms of studying its behaviour and physical properties. The problem of interest is to find the fuel optimal lowthrust capture trajectory from the original asteroid orbit to the Sun-Earth L2 libration point. The departure time has been set between 2456000.5 and 2460000.5 JD and the maximum time of flight allowed is 1800 days. We have followed a methodology which in this paper is divided in four main sections. The first section of this paper is dedicated to prune the asteroids from the Near Earth Asteroids database and select the small- 2 Asteroids Database Selection In this research, JPL’s Solar System Dynamics Group small-body database (SBDB) has been used to select the asteroids and get the orbital elements and stellar magnitude data for each one. As of 5th of September of 2012, 9049 NEAs (Near Earth Asteroids) have been identified. An initial pruning of the asteroids has been made with the criteria to select the smallest ones within the Earth’s neighborhood4 , in order to be capable to move them ∗ Elecnor Deimos, Spain. [email protected] University, Canada. [email protected] ‡ IEEC & Universitat Politècnica de Catalunya, Spain. [email protected] § IEEC & Universitat de Barcelona, Spain, [email protected] † York 1 63rd International Astronautical Congress, Naples, Italy. Copyright 2012 by the International Astronautical Federation. All rights reserved. with the current technology of a Solar Electric Propulsion (SEP) system. Then, the constraints in the list of the database asteroids are a semi-major axis range between 0.85-1.15 AU and a stellar magnitude bigger than 28. The stellar magnitude represents the brightness of the object, a direct relationship with its size. Applying these constraints, we get a final selection of 40 asteroids with the parameters in Table 1. The asteroid 2004 UH1, has been discarded because it has a very high eccentricity, plus it is a type of asteroid whose orbit intersects with the Earth orbit, but the relative speed is very high. The combination of these facts would need a very big effort to change the orbit. Name Class a [AU] e i [deg] Ω [deg] ω [deg] M [deg] n [deg/day] H 1991 VG 2000 LG6 2003 SW130 2003 WT153 2006 BV39 2006 JY26 2006 RH120 2007 EK 2007 UN12 2008 CM74 2008 GM2 2008 HU4 2008 JL24 2008 KT 2008 LD 2008 UA202 2008 UC202 2008 WO2 2009 BD 2009 WQ6 2009 WW7 2009 WR52 2009 YR 2010 JW34 2010 RF12 2010 UY7 2010 UE51 2010 VL65 2010 VQ98 2011 AM37 2011 BQ50 2011 CA7 2011 CH22 2011 JV10 2011 MD 2011 UD21 2012 AQ 2012 EP10 2012 FS35 Apollo Aten Aten Aten Apollo Apollo Apollo Apollo Apollo Apollo Apollo Apollo Apollo Apollo Aten Apollo Apollo Apollo Apollo Aten Apollo Apollo Aten Aten Apollo Aten Apollo Apollo Apollo Apollo Aten Apollo Aten Apollo Amor Aten Apollo Apollo Apollo 1.03 0.92 0.88 0.89 1.15 1.01 1.03 1.13 1.05 1.09 1.05 1.09 1.04 1.01 0.89 1.03 1.01 1.03 1.06 0.87 1.09 1.03 0.94 0.98 1.06 0.90 1.06 1.07 1.02 1.10 0.95 1.08 0.88 1.14 1.06 0.98 1.07 1.05 1.10 0.0491 0.1109 0.3043 0.1777 0.2714 0.0830 0.0245 0.2724 0.0605 0.1469 0.1572 0.0733 0.1066 0.0848 0.1547 0.0686 0.0685 0.1882 0.0516 0.4087 0.2618 0.1551 0.1102 0.0548 0.1882 0.1499 0.0597 0.1440 0.0271 0.1473 0.0982 0.2888 0.2358 0.2020 0.0371 0.0302 0.1038 0.1160 0.1185 1.45 2.83 3.67 0.37 0.74 1.44 0.60 1.21 0.24 0.86 4.10 1.26 0.55 1.98 6.54 0.26 7.46 2.01 1.27 5.82 2.53 4.24 0.70 2.26 0.88 0.46 0.62 4.40 1.48 2.63 0.36 0.12 0.13 1.40 2.45 1.06 2.86 1.03 2.34 73.98 72.55 176.45 55.61 127.09 43.50 51.14 168.58 216.11 321.58 195.11 221.34 225.82 240.66 250.90 21.06 37.43 238.15 253.33 55.68 57.18 61.03 86.95 49.81 163.85 39.95 32.29 223.12 46.17 291.28 281.01 311.00 334.67 221.39 271.63 22.52 97.32 348.04 186.57 24.51 8.19 47.80 148.91 74.96 273.45 10.14 83.26 134.34 242.73 278.25 341.50 281.97 101.86 201.42 300.89 91.24 85.70 316.73 227.27 273.71 269.88 127.87 43.61 267.56 210.44 47.25 253.97 341.60 129.20 1.27 278.61 27.59 297.52 5.84 208.45 316.09 105.73 42.23 340.17 185.75 49.55 55.61 116.33 29.55 221.25 181.71 238.24 339.86 121.93 327.11 124.19 7.44 202.43 330.08 230.71 331.19 115.11 288.88 241.39 329.31 257.46 294.67 254.92 228.37 239.36 203.68 316.69 213.92 141.56 108.09 115.28 101.68 56.38 144.91 280.99 249.95 126.35 0.95 1.12 1.19 1.17 0.80 0.97 0.94 0.82 0.91 0.87 0.91 0.86 0.93 0.97 1.17 0.94 0.97 0.95 0.90 1.22 0.87 0.94 1.08 1.01 0.90 1.16 0.91 0.90 0.95 0.85 1.06 0.88 1.20 0.81 0.91 1.02 0.89 0.92 0.86 28.39 29.019 29.117 28.048 28.984 28.349 29.527 29.258 28.741 28.043 28.356 28.223 29.572 28.215 28.864 29.44 28.242 29.779 28.236 29.186 28.894 28.32 28.003 28.148 28.369 28.527 28.311 28.423 28.2 29.69 28.341 30.319 28.961 29.706 28.073 28.483 30.698 29.165 30.286 Table 1: Orbital Elements, stellar magnitude and class of the Asteroids selected 2 63rd International Astronautical Congress, Naples, Italy. Copyright 2012 by the International Astronautical Federation. All rights reserved. where pv is the geometric albedo and H is the asteroid absolute magnitude from the SBDB database. In Table 3, we have represented the diameter and the corresponding mass that would have an asteroid with the correspondent stellar magnitude indicated, supposing that it is calculated with the Method 2. This is useful to be able to compare the diameter in both methods. In the upper limit differs up to 37 %, though the lower limit difference is only 6.5 %. If we compare the masses in both methods, we can notice that the upper limit of mass is very different in both methods because in the first method we have 2.1 Mass Estimation considered an albedo equal to 0.25 and in method 2 equal The mass of an asteroid is very difficult to obtain as it can to 0.5. only be determined by in-situ measurements or from the H D [m] Mass [Kg] observed dynamics (i.e., spacecraft tracking during en5 counters, natural satellites of asteroids). In this research, 28 4.72 − 14.93 143249 - 4529934 two methods to estimate the mass of the selected asteroids have been conducted both assuming an spherical shape of 28.5 3.75 − 11.86 71794 - 2270345 the asteroid and a bulk density6 of 2.6 g/cm3 . The mass estimation in the first method has been made 29 2.98 − 9.42 35982 - 1137868 from the relationship stated in JPL7 between the absolute magnitude H and the diameter of the asteroid, which as29.5 2.37 − 7.48 18033 - 570284 sumes an albedo ranging from 0.25 to 0.05. This method grouped the asteroids based on its absolute magnitude. 30 1.88 − 5.94 9038 - 285819 H D [m] Mass [Kg] The vast majority of the asteroids selected are from the Apollo family, which have a semi-major axis greater than 1 AU and a perihelion distance q < 1.017AU . In fact, 38 out of the 39 asteroids selected cross the Earth’s orbit and they are rather Apollo or Aten; except one that belongs to the Amor family (1 < q < 1.3 AU and a > 1 AU). The highest inclination is 7.45 deg from the ecliptic plane and the range of eccentricity among the asteroids selected is between 0.0244 and 0.4086. 28 7 − 15 466945 - 4594579 28.5 5 − 12 170169 - 2352424 4−9 87126 - 992429 Table 3: Method 2 relationship between H, D and Mass considering an albedo range of 0.05 - 0.5 Therefore, two absolute limits for each mass method have been estimated for asteroids, reflecting the uncertainty of the mass determination and being able to give 29.5 3 − 7 36756 - 466945 a range of the asteroid mass that the mission could find. In the following sections the procedure to be able to 30 3 − 6 36756 - 294053 capture asteroids is presented. The procedure is not trivial and requires a serial of steps to find an optimal trajectory to capture the asteroids. The objective of these previous steps from the actual refined trajectory optimization is to Table 2: Relationship from JPL between H, D and Mass dispose of suitable good initial guessed parameters such In practice, as we are studying the asteroids with as the departure time, the time of flight and the history H > 28, there will be 5 groups as it can be seen in Table of control of the trajectory to be able to converge to the 2. We have classified the asteroids which have a stellar optimal solution. magnitude between 28-28.3 in the group H=28, the asteroids with an H between 28.3-28.7 in the H=28 group, the 3 Transfer opportunity search asteroids with an H between 28.7-29.2 in the H=29, and so on. In the second method the diameter has been derived The goal of the initial orbit search step is to find a suitfrom the Eq. 1, based on the absolute magnitude H, ac- able departure time and time of flight. In this research, cording to Fowler and Chillemi8 and assuming an albedo we concluded that it is very important to have a good first estimation of them because they are very sensible in the range2 between 0.05 to 0.5. trajectory optimization. A low-thrust trajectory modeled as a series of impulses connected by conic arcs proposed 1329 −0.2H D = √ 10 , (1) by Sims and Flanagan9 has been applied. pv 29 3 63rd International Astronautical Congress, Naples, Italy. Copyright 2012 by the International Astronautical Federation. All rights reserved. N *3, the number of nodes of the ballistic segment multiplied by the magnitude plus the angles of thrust. Upper and lower bounds are set for the variables between 0 and 1, except for the departure time and time of flight. The departure time in our problem has been set between 2456000.5 − 2460000.5 JD and the time of flight between 200 − 1800 days. The matching conditions of the global optimization consists in the position and velocity of the encounter point. The positions should be concurrent and the velocities mismatch should be less than a tolerance error. The performance index of the global optimization problem includes the constraints in the objective function. Thus, the objective function to be minimized is the sum of the impulses in each node and the state vector mismatches: The model consists on a low-thrust trajectory divided into two segments and modeled as a N series of impulsive maneuvers connected by conic Lambert arcs, separated in an equal time step nodes. The trajectory is propagated (two-body regime) forward in time from the asteroid to the connection point and backward from the Sun-Earth Lagrangian Point L2 to the connection point. The Lagrangian Point has been geometrically defined using JPL ephemerides DE405. The trajectory is shown in Figure 1. The propagation between impulses is according to a sequential Kepler model, avoiding in this case the numerical integration, which is the most time consuming part of the trajectory generation. J= n X i=1 ki + λr ||Rconn ||2 + λv ||Vconn ||2 where λr and λv are the weighting factors, Rconn is the position error at the connection point, and Vconn is the velocity error at the connection point. The variables are expressed in canonical units, where the distance unit for the position mismatches is in AU . The weighting factors affect the feasibility and optimality of the trajectory found. The higher the weighting parameters are, the optimizer focuses more on the constraints. In our program, they have been set to 105 . By using the Differential Evolution algorithm, a solution can be found very fast. However, the position and the velocity mismatches at the connection point are rather big. In fact, the velocity error is bigger than the capability of low-thrust engine. Therefore, the trajectory obtained is not a feasible solution and only the rough departure time and time of flight have been used to find the optimal trajectory. The option to improve the results using a local optimizer has also been studied. The time history of control gets smoother using this method, but the improvement in terms of matching is insignificant. Therefore, we have decided to skip this step and go directly to the next one because the outputs from the global optimization method are already acceptable initial seeds for the trajectory optimization. Fig. 1: Trajectory model10 The impulsive maneuver ∆V on each node is defined by three parameters α, β and k, which represent the magnitude and the direction of the thrust. ∆V = ∆Vmax k[cos β cos α, cos β sin α, sin β]T where α ∈ [0, 2π), β ∈ [− π2 , π2 ], k ∈ [0, 1]. The ∆V at each node should not exceed maximum magnitude ∆Vmax = atstep , where a is the acceleration offered by the low-thrust engine when operated at full thrust and tstep is the time span between the two nodes. In our computations, the maximum magnitud of the acceleration has been set to 10−4 m/s2 . The global search of the trajectory10 ,11 is based on Differential Evolution algorithm. The version of Differential Evolution algorithm implemented is from Storn and Price12 , with a weighting factor F = 0.8 and a crossover constant CR = 0.8. For each run, the number of maximum iterations has been set to 20000 and the population size to 50. It is worth to mention that the sampling has been made uniformly distributed over the surface of the sphere. The number of optimization parameters is the sum of the Time of Flight (TOF), Departure Time (DT) and 4 Trajectory Generation By using the method introduced in section 3, the rough transfer trajectory, which is approximated by a series of Lambert arcs, is obtained. However, this transfer trajectory is inaccurate and its optimality also cannot be guaran4 63rd International Astronautical Congress, Naples, Italy. Copyright 2012 by the International Astronautical Federation. All rights reserved. teed (The general performance index contains the penalty function. The weight of the penalty function changes the characteristics of the optimization problem). In this section, the accurate trajectory optimization technique will be introduced. The sequential keplerian approximation is replaced by the integration of the ordinary differential equations of motion, which includes a more accurate lowthrust orbital dynamic model. A series of piecewise continuous control histories is generated by using the conjugated equations in calculus of variation. Then, these control histories are optimized by Runge-Kutta 4th order shooting method, and a series of local optimal solutions are obtained (The Runge-Kutta 4th order shooting method is a kind of gradient based optimization method. Therefore, like all the rest of gradient based algorithms, it can only obtain the local optimal solution). The solution which has the best performance index will be chosen, and considered as the global optimal solution. 4.1 HT = − T =0 T = Tmax 0 < T < Tmax v̇ = ṁ = λ˙r = HT > 0 HT < 0 HT = 0 1. We generate the initial values of the conjugated states [λTr , λTv , λm ] randomly, setting the limits of the values between −100 to 100 to have the same order of magnitude. The typical dynamical equations, expressed in cartesian coordinates, the conjugated equations and the optimal control equations to obtain the initial seed of the optimization problem are presented below. In this case the mass is propagated along the trajectory, while in the global optimization it was assumed to keep it constant throughout the trajectory. = if if if (4) where g = 9.8m/s2 and the engine parameters are Isp = 3000s and Tmax is set to when the engines operate at full thrust, which will depend on the case studied. The direct method does not need the conjugated equations in Eq. 2, but in this research the conjugated equations have been used to generate the initial guesses of the piecewise continuous time history of control. The procedure is listed as follows: Initial guess generation ṙ ||λv || λm − m gIsp 2. Then, we use these initial values, the departure time and the time of flight, which are obtained from the global search step, to propagate the set of Eq. 2 ; 3. If the error of the final states is smaller than a given tolerance (position error less than 0.5 AU and velocity less than 0.2 VU), the set of initial conjugated states [λTr , λTv , λm ] is accepted, otherwise it is rejected and step 1-3 are repeated, until a set of acceptable conjugated states is obtained; v µ T r+ γ r3 m T − g0 Isp − 4. Generate the discrete time history of control by using the control equation and switching function of Eq. (4) (In all computations the number of nodes that has been used is 81). 3λTv r µ − r r3 r5 −λ T −||λv || 2 m λv Using this procedure, a series of trajectories whose final states are not very far from the target state can be obλ̇m = (2) tained, and meanwhile the time histories of control are also obtained. In this research, 50 sets of time history of where, T and γ are the control variables. The perfor- control are obtained in this step. mance index is defined as mf 4.2 Local optimization −→ max (3) J= m0 The next step is to refine the control using the Rungewhere, mf and m0 are the final mass and the initial mass Kutta 4th order shooting method.13 The basic idea of the of the spacecraft-asteroid assembly, respectively. By us- Runge-Kutta 4th order shooting method is to convert oping the Pontryagin’s maximum principle, the optimal con- timal control problem to constrained parameter optimizatrol can be obtained. tion problem. The optimization parameters of this problem can be denoted as z = [u0 , u1 , · · · , uN ]T , where λv γ=− ui = [αi , βi , Ti ]T , whose initial guesses are already ob||λv || tained by using random generation method. The lower T is determined by the switching function: bound and upper bound of [αi , βi , Ti ]T are [0, − π2 , 0]T λ˙v = 5 63rd International Astronautical Congress, Naples, Italy. Copyright 2012 by the International Astronautical Federation. All rights reserved. and [2π, π2 , Tmax ]T , respectively. The constraints can be when executed in dual core Intel Xeon 3.2 GHz CPU. written as However, the case of 2012 FS35 Method 2 lower limit only needed 1:34.02 min of computational time. ∆ = [r(tf ), v(tf )]T − [rf , vf ]T (5) The outputs from global optimization (GO) (listed in where, [r(tf ), v(tf )]T is the final state obtained by nu- Table 4) problem are rough estimates for two reasons. On merical integration, and [rf , vf ]T is the target state. The the one hand, in some cases we did not get any converged performance index is the same as Eq. 3. result using directly the outputs from GO. Therefore, we Using the fixed step-size 4th order Runge-Kutta in- decided to change the time of flight between 6-33 Time tegral formula, the state equations can be propagated. Units (348.79 - 1918.37 days) with a time step of 3, which Meanwhile, the derivatives of the constraints and per- is approximate one orbit, and keep the converged result formance index with respect to optimization parameters, with the lowest time of flight. On the other hand, the pa∂J which are denoted as ∂∆ ∂z and ∂z , can be obtained. All rameters also change because the transfer of the asteroid is derivatives are expressed anallytically. In our research, considered from the first time that the engine opens until the local optimal solution of this parameters optimiza- it stops for the last time. Therefore, as the control profile tion problem can be found using sequential quadratic pro- we have obtained is different for each method, the time gramming (SQP) algorithm.14 of flight and departure time vary. Although, the objective function in this research is the maximization of the final mass. 5 Results In the trajectory optimization problem, the thrust of the engine is one input parameter. The criteria that we In the previous sections a methodology to capture asterhave followed is to start with the lowest one (5 N) and oids has been described. Starting from the estimation of search if for different time of flight it converged. Otherthe asteroids mass, finding a rough time of flight and dewise, we increased the thrust. parture time. Then, setting the procedure to generate a As we can see in Table 5 and 6, the thrust needed in suitable seed for the calculation of the optimal trajectory the lowest mass limit is in the range between 5 − 25 N. to capture an asteroid. This section is devoted to present However, in the upper limit mass cases the thrust needed the method of the computations itself and discuss the reis much bigger, in the range between 15 - 200 N. This sults found. one order of magnitude difference is obviously due to the The transfer opportunity search has been done one mass increment of the mass needed for the transfer. time for each asteroid, needing about 4 minutes per comTherefore, if the asteroids would have the lower limit putation. In section 4, a total of 156 trajectories has been of the mass estimated in Method 2, only using one VAScomputed, which is the combination of the 39 selected asMIR engine would be enough to capture them. Otherwise, teroids among the two methods and the two boundaries more than one engine would be required to move them. of each one. The computational time of the trajectory generation strongly depends on the time of flight of each The time when the thrust is open is a very interesting computation case. For instance, 2010 UY7 Method 2 up- parameter because it gives us an idea of the real cost to do per limit of mass estimation, 12:41.15 min time is needed the transfer. 6 63rd International Astronautical Congress, Naples, Italy. Copyright 2012 by the International Astronautical Federation. All rights reserved. Name DTgl [JD] 1991 VG 2000 LG6 2003 SW130 2003 WT153 2006 BV39 2006 JY26 2006 RH120 2007 EK 2007 UN12 2008 CM74 2008 GM2 2008 HU4 2008 JL24 2008 KT 2008 LD 2008 UA202 2008 UC20 2008 WO2 2009 BD 2009 WQ6 2457589.62 2458636.74 2456128.74 2456220.47 2457303.08 2456976.17 2460000.50 2458014.42 2458639.10 2457264.83 2458974.37 2457224.88 2459947.31 2456408.20 2459989.80 2459999.64 2460000.50 2456000.50 2459425.83 2457438.15 Global Optimization T OFgl [d] Name 976.84 1498.08 912.35 1569.58 1192.22 1500.00 1618.54 1220.78 883.32 1256.63 1453.31 1278.91 1628.75 1046.38 985.79 1443.84 1395.18 1799.58 1046.38 1569.58 2009 WW7 2009 WR52 2009 YR 2010 JW34 2010 RF12 2010 UY7 2010 UE51 2010 VL65 2010 VQ98 2011 AM37 2011 BQ50 2011 CA7 2011 CH22 2011 JV10 2011 MD 2011 UD21 2012 AQ 2012 EP10 2012 FS35 DTgl [JD] T OFgl [d] 2457444.51 2456024.82 2458046.12 2456006.76 2459962.67 2456970.35 2459989.98 2458982.79 2456000.50 2457412.01 2459410.93 2458886.04 2456895.09 2458967.43 2460000.50 2456000.50 2458872.57 2456000.75 2458396.83 1192.99 1154.47 1220.78 929.11 1499.15 1278.91 574.26 1245.14 1046.38 1307.15 1084.44 1064.23 1119.30 1351.09 1005.69 1226.28 1337.28 1453.31 1213.16 Table 4: List of global optimization section results ( departure time (DTgl ) and time of flight (T OFgl ) ) of the 39 pruned asteroids (H > 28 and 0.85 <semi-major axis< 1.15). Table 5: Outputs of the Method 1 in trajectory optimizalimit range. These 5 columns for each limit are the tion section of the 39 pruned asteroids (H > 28 Thrust (T), acceleration (a), departure time (DT), and 0.85 <semi-major axis< 1.15). A column of time of flight (TOF) and time when the thrust is on the asteroid stellar magnitude is written, specifying (tON ). in parenthesis in which group of Method 1 mass estimation belongs. A is in the group of H = 28, Table 6: List of results from Method 2 trajectory optimization in upper and lower limit of mass estimaB in the group of H = 28.5 , C in the group of tion, from the 39 pruned asteroids (H > 28 and H = 29, D in the group of H = 29.5 and E in 0.85 <semi-major axis< 1.15). The outputs specthe group of H = 30. Next to that, the following ified in the table are the thrust (T), the acceleration 5 columns belong to the minimum limit range of (a), the mass, the departure time, the time of flight Method 1 and the last 5 columns to the maximum and the time when the thrust is on (tON ). 7 Method 1: H 28.39 (B) 29.02 (C) 29.12 (C) 28.05 (A) 28.98 (C) 28.35 (B) 29.53 (D) 29.26 (C) 28.74 (C) 28.04 (A) 28.36 (B) 28.22 (A) 29.57 (D) 28.22 (A) 28.86 (C) 29.44 (D) 28.24 (A) 29.78 (E) 28.24 (A) 29.19 (C) 28.89 (C) 28.32 (B) 28.00 (A) 28.15 (A) 28.37 (B) 28.53 (B) 28.31 (B) 28.42 (B) 28.20 (A) 29.69 (D) 28.34 (B) 30.32 (E) 28.96 ( C) 29.71 (D) 28.07 (A) 28.48 (B) 30.70 (E) 29.17 (C) 30.29 (E) T [N] 15 15 15 25 15 15 5 15 15 15 15 15 5 15 15 5 25 5 15 15 15 15 15 15 15 15 15 15 15 5 15 5 15 5 15 15 5 15 5 Lower Limit Results a [m/s2 ] DT [JD] TOF [d] 8.81E-005 2457591.64 771.70 1.72E-004 2458638.80 1363.25 1.72E-004 2456128.90 884.98 5.35E-005 2456220.47 1491.10 1.72E-004 2457304.51 858.40 8.81E-005 2456980.30 1260.00 1.36E-004 2460000.50 1602.35 1.72E-004 2458016.94 854.55 1.72E-004 2458639.71 786.15 3.21E-005 2457266.12 1181.23 8.81E-005 2458977.62 1249.85 3.21E-005 2457224.88 1266.12 1.36E-004 2459955.43 1091.26 3.21E-005 2456408.92 1004.53 1.72E-004 2459992.17 739.34 1.36E-004 2460002.13 1169.51 5.35E-005 2460000.50 1395.18 1.36E-004 2456005.45 1511.64 3.21E-005 2459427.81 931.28 1.72E-004 2457438.15 755.72 1.72E-004 2457446.36 1025.97 8.81E-005 2456026.61 1050.57 3.21E-005 2458047.17 1147.53 3.21E-005 2456006.76 929.11 8.81E-005 2459967.57 989.44 8.81E-005 2456971.89 1176.60 8.81E-005 2459989.98 333.07 8.81E-005 2458983.86 1058.37 3.21E-005 2456000.68 1035.92 1.36E-004 2457413.81 928.08 8.81E-005 2459411.30 1062.75 1.36E-004 2458887.32 936.52 1.72E-004 2456898.36 917.83 1.36E-004 2458967.43 1067.36 3.21E-005 2460000.85 985.58 8.81E-005 2456003.88 993.29 1.36E-004 2458874.87 922.72 1.72E-004 2456005.25 1191.72 1.36E-004 2458399.55 788.56 tON [d] 126.99 209.73 392.31 674.92 274.21 285.00 129.48 305.20 97.16 1130.97 479.59 767.35 114.01 805.72 285.88 115.51 1088.24 287.93 481.34 619.69 357.90 565.69 793.51 724.71 359.80 345.31 126.34 510.51 523.19 274.50 260.26 372.48 257.44 283.73 683.87 171.68 240.71 130.80 254.76 T [N] 50 50 100 150 50 50 15 25 25 100 100 100 15 100 50 15 150 15 100 150 50 100 125 100 100 100 50 100 100 15 100 15 100 15 100 50 15 25 15 Upper Limit Results a [m/s2 ] DT [JD] TOF [d] 2.13E-005 2457590.97 849.85 5.04E-005 2458638.18 1311.47 1.01E-004 2456129.19 845.83 3.26E-005 2456220.56 517.96 5.04E-005 2457303.08 1353.32 2.13E-005 2456976.17 1500.00 3.21E-005 2460000.50 1602.35 2.52E-005 2458014.42 1918.37 2.52E-005 2458639.10 883.32 2.18E-005 2457266.48 1822.45 4.25E-005 2458975.77 1082.90 2.18E-005 2457224.88 1569.58 3.21E-005 2459949.55 1465.88 2.18E-005 2456408.51 1782.00 5.04E-005 2459989.80 1569.58 3.21E-005 2459999.89 1357.21 3.26E-005 2460000.50 655.73 5.10E-005 2456003.20 1365.53 2.18E-005 2459426.96 1188.24 1.51E-004 2457438.15 802.23 5.04E-005 2457445.47 1339.37 4.25E-005 2456026.37 1413.07 2.72E-005 2458046.12 1220.78 2.18E-005 2456006.76 1381.23 4.25E-005 2459962.67 1615.42 4.25E-005 2456970.80 2042.30 2.13E-005 2459989.98 871.99 4.25E-005 2458983.43 1207.79 2.18E-005 2456001.02 1469.87 3.21E-005 2457412.01 1307.15 4.25E-005 2459411.49 1051.90 5.10E-005 2458886.31 1412.62 1.01E-004 2456895.24 863.27 5.10E-005 2458967.43 1351.09 2.18E-005 2460000.74 1381.23 2.13E-005 2456000.50 1214.02 5.10E-005 2458873.49 1283.79 2.52E-005 2456000.75 1569.58 5.10E-005 2458396.83 1152.51 tON [d] 595.87 627.83 662.71 423.79 1074.29 1215.00 420.82 1841.64 582.99 1285.31 1024.68 1255.66 586.35 1170.00 1334.14 765.24 418.55 784.79 822.63 671.43 1102.19 1127.45 1062.08 948.72 726.94 793.50 828.39 996.11 674.94 1176.43 531.37 1020.22 531.91 756.61 1018.48 649.93 575.03 1130.09 655.11 63rd International Astronautical Congress, Naples, Italy. Copyright 2012 by the International Astronautical Federation. All rights reserved. 8 Name 1991 VG 2000 LG6 2003 SW130 2003 WT153 2006 BV39 2006 JY26 2006 RH120 2007 EK 2007 UN12 2008 CM74 2008 GM2 2008 HU4 2008 JL24 2008 KT 2008 LD 2008 UA202 2008 UC20 2008 WO2 2009 BD 2009 WQ6 2009 WW7 2009 WR52 2009 YR 2010 JW34 2010 RF12 2010 UY7 2010 UE51 2010 VL65 2010 VQ98 2011 AM37 2011 BQ50 2011 CA7 2011 CH22 2011 JV10 2011 MD 2011 UD21 2012 AQ 2012 EP10 2012 FS35 T [N] 5 5 5 10 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 0.1 5 5 Lower Limit Results a [m/s2 ] Mass [kg] DT 5.98E-005 83578 2457591.64 1.43E-004 35050 2458636.74 1.63E-004 30612 2456130.63 7.46E-005 134058 2456221.22 1.36E-004 36787 2457304.31 5.65E-005 88449 2456979.79 2.88E-004 17374 2460000.50 1.98E-004 25194 2458015.02 9.72E-005 51463 2458639.71 3.70E-005 134987 2457266.56 5.71E-005 87598 2458974.37 4.75E-005 105267 2457225.03 3.06E-004 16326 2459950.31 4.70E-005 106437 2456408.82 1.15E-004 43420 2459989.97 2.55E-004 19593 2459999.64 4.88E-005 102540 2460000.50 4.08E-004 12266 2456000.50 4.84E-005 103393 2459429.43 1.80E-004 27829 2457438.15 1.20E-004 41657 2457446.15 5.43E-005 92064 2456026.01 3.50E-005 142657 2458046.12 4.28E-005 116759 2456006.76 5.81E-005 86038 2459962.67 7.23E-005 69166 2456971.82 5.36E-005 93216 2459989.98 6.26E-005 79853 2458982.79 4.60E-005 108666 2456001.55 3.60E-004 13870 2457413.09 5.59E-005 89432 2459411.49 8.60E-004 5817 2458886.04 1.32E-004 37974 2456897.04 3.69E-004 13567 2458967.43 3.86E-005 129507 2460001.19 6.80E-005 73501 2456003.88 2.90E-005 3446 2458872.57 1.75E-004 28648 2456000.75 8.21E-004 6088 2458396.83 TOF [d] 771.70 685.96 766.37 776.07 882.24 1290.00 515.93 662.71 794.98 1156.10 871.99 854.55 523.19 1710.00 975.93 592.95 1569.58 455.76 1031.55 784.79 1014.04 1085.21 1162.65 1014.99 952.21 1086.50 459.41 1245.14 1123.12 627.83 1051.90 232.53 749.91 610.39 965.46 993.29 523.19 697.59 232.53 tON [d] 234.44 370.42 419.68 523.19 369.59 420.00 175.41 279.04 159.00 892.21 767.35 497.03 90.69 486.00 453.46 188.35 1224.27 176.72 313.38 540.63 477.19 969.76 732.47 429.02 742.93 427.27 246.93 722.18 317.40 111.61 368.71 223.23 374.95 140.39 543.07 220.73 292.99 174.40 137.19 T [N] 50 50 125 200 50 50 50 50 50 100 100 100 25 100 125 50 150 25 50 150 100 100 100 100 100 100 100 100 100 25 100 25 125 25 100 50 50 25 25 Upper Limit Results a [m/s2 ] Mass [kg] 1.89E-005 2642968 4.51E-005 1108388 1.29E-004 968037 4.72E-005 4239278 4.30E-005 1163300 1.79E-005 2796996 9.10E-005 549404 6.28E-005 796695 3.07E-005 1627391 2.34E-005 4268663 3.61E-005 2770077 3.00E-005 3328832 4.84E-005 516288 2.97E-005 3365828 9.10E-005 1373068 8.07E-005 619572 4.63E-005 3242589 6.45E-005 387875 1.53E-005 3269580 1.70E-004 880019 7.59E-005 1317322 3.43E-005 2911333 2.22E-005 4511198 2.71E-005 3692257 3.68E-005 2720770 4.57E-005 2187217 3.39E-005 2947759 3.96E-005 2525177 2.91E-005 3436307 5.70E-005 438624 3.54E-005 2828081 1.36E-004 183947 1.04E-004 1200859 5.83E-005 429034 2.44E-005 4095358 2.15E-005 2324299 4.59E-004 108966 2.76E-005 905924 1.30E-004 192527 DT 2457589.62 2458638.03 2456130.72 2456221.22 2457303.08 2456976.17 2460000.78 2458014.42 2458639.10 2457264.83 2458975.81 2457224.88 2459949.27 2456408.51 2459989.80 2460002.62 2460000.50 2456003.60 2459426.64 2457438.15 2457444.81 2456026.02 2458046.12 2456006.76 2459962.67 2456970.61 2459989.98 2458982.79 2456001.27 2457412.01 2459411.12 2458886.58 2456895.09 2458969.75 2460000.74 2456000.50 2458872.57 2456000.75 2458399.34 TOF [d] 928.00 1423.17 920.82 1515.84 1395.18 1691.65 1602.35 1046.38 874.48 1256.63 1311.47 1035.92 1433.30 1782.00 1241.71 1169.51 1743.97 1601.62 1506.79 793.51 1726.53 1674.21 1495.43 929.11 1569.58 1485.00 574.26 1395.18 1454.87 1267.93 1073.59 947.17 1220.78 1202.47 1381.23 1214.02 348.79 1495.44 824.95 tON [d] 722.86 764.02 481.34 866.20 1311.47 1639.33 161.85 983.60 441.66 879.64 1171.95 753.40 325.75 630.00 613.88 216.58 1255.66 593.86 1145.79 566.79 1430.06 1360.30 807.53 631.80 1083.01 1110.00 511.09 1269.61 419.96 575.15 694.04 393.77 964.42 702.57 878.96 637.67 136.03 1181.40 254.76 63rd International Astronautical Congress, Naples, Italy. Copyright 2012 by the International Astronautical Federation. All rights reserved. 9 Method 2: Name 1991 VG 2000 LG6 2003 SW130 2003 WT153 2006 BV39 2006 JY26 2006 RH120 2007 EK 2007 UN12 2008 CM74 2008 GM2 2008 HU4 2008 JL24 2008 KT 2008 LD 2008 UA202 2008 UC20 2008 WO2 2009 BD 2009 WQ6 2009 WW7 2009 WR52 2009 YR 2010 JW34 2010 RF12 2010 UY7 2010 UE51 2010 VL65 2010 VQ98 2011 AM37 2011 BQ50 2011 CA7 2011 CH22 2011 JV10 2011 MD 2011 UD21 2012 AQ 2012 EP10 2012 FS35 63rd International Astronautical Congress, Naples, Italy. Copyright 2012 by the International Astronautical Federation. All rights reserved. 6 Conclusions In this research, a procedure to capture an asteroid in the L2 libration point has been described. It has been divided in the transfer opportunity search and the trajectory generation sections, obtaining a low-thrust trajectory for each case with an on-off profile control. A list of the asteroids capture opportunities has been obtained indicating the results of the most important parameters such as the acceleration, the time of flight and the duration time of the thrust open. A converged solution has been found for all asteroids with all mass estimation combinations. Acknowledgements This work has been supported by the grants MTM201016425, MTM2009–06973 and 2009SGR859. We also acknowledge the use of the UPC Applied Math cluster system for research computing (see http://www.ma1.upc.edu/eixam/index.html). References [4] J.P. Sanchez. and C.R. McInnes On the Ballistic Capture of Asteroids for Resource Utilisation. International Astronautical Congress, 2011. [5] J.L. Hilton. Asteroid Masses and Densities. U.S. Naval Observatory. [6] Chesley, S. R., Chodas, P. W., Milani, A., Valsecchi, G. B. and Yeomans, D. K. Quantifying the Risk Posed by Potential Earth Impacts. Icarus, Vol. 159, 2002, pp.423-432., 2002. [7] Jet Propulsion Laboratory. Near Earth Observation Program http://neo.jpl.nasa.gov/glossary/h.html, [Accessed on 5th of September 2012]. [8] Fowler J.W. and Chillemi J.R. The IRAS Minor Planet Survey. 17,43, 1992. [9] J.A. Sims and S.N. Flanagan. Preliminary Design of Low-Thrust Interplanetary Missions. AAS/AIAA Astrodynamics Specialist Conference, 1999. [10] C. H. Yam, D. Izzo and F. Biscani Global Optimization of Low-Thrust Trajectories via Impulsive DeltaV Transcription. ESA - Advanced Concepts Team, 2009. [1] H. Baoyin, Y. Chen and J. Li Capturing Near Earth [11] C. H. Yam, D. Izzo and F. Biscani Towards a High Fidelity Direct Transcription Method for OptimisaObjects. Research in Astronomy and Astrophysics, tion of Low-Thrust Trajectories. ESA - Advanced 2011. Concepts Team, 2010. [2] J.P. Sanchez., D. Garcia Yarnoz and C.R. McInnes Near-Earth Asteroid Resource Accessibility and Fu- [12] K. Price and R. Storne Differential Evolution (DE)for Continuous Function Optimization. ture Capture Mission Opportunities. Global Space http://www1.icsi.berkeley.edu/ storn/code.html, Exploration Conference, 2012. [Accessed on 5th of September 2012]. [3] B.W. Longmier, J.P. Squire, L.D. Cassady, M.G. [13] Y. Gao Advances in low-thrust trajectory optimizaBallenger, M.D. Carter, C. Olsen, A. V. Ilin, T. W. tion and flight mechanics. PhD Thesis, 2003. Glover, G. E. McCaskill, F. R. Chang Dı́az, E.A. Bering and J. Del Valle. VASIMR VX-200 Perfor- [14] Philip E. GILL. Users Guide for SNOPT Version 7: mance Measurements and Helicon Throttle Tables Software for Large-Scale Nonlinear Programming. Using Argon and Krypton. 32nd International ElecDepartment of Mathematics. University of Califortric Propulsion Conference, 2011. nia, San Diego, La Jolla, CA 92093-0112, USA. 10