1. INTRODUCTION

Recently, authors such as (^{Gang et al. 2014}; ^{Shan et al. 2014}; ^{Zotes et al. 2012}) have given important contributions to the study of spatial trajectory optimization. The latter have considered a geometric method with results far from those presented in this article, by having a large flight time with a v similar to ours.

Thus, to minimize fuel consumption in the acceleration v [as in (^{Sharaf & Saad 2016})], the present article considers the methods proposed in ^{Leeghim (2013)}, which are variations of the Lambert problem, with neither direction of motion in a plane nor the time t of the transfer at the beginning given; both are calculated. Another difference of our approach is the obtainment of the wait time t₁ (^{Vallado 1997}), where the interceptor does not leave its initial orbit until the relative positions of the bodies are convenient. This is presented in § 2.

The problem becomes a constrained optimization system. Its direct solution leads to several systems of nonlinear equations (9 equations as in § 3, 11 equations as in § 4). The solution method, presented in § 5 (Newton’s method) needs an appropriate starting point (seed). This seed is usually found from an exhaustive model that takes a long time to locate a suitable value (grids). Therefore, an heuristic technique, the particle swarm optimization presented in § 5.1, is used to approximate the seed in a reasonable time. Such methods were applied to trips with real data in the solar system and are presented in § 6. The results of the first case (immediate transfer) were compared with information released from missions of several space agencies in ^{Villanueva (2013)} and here the two problems are compared with the data shown in ^{Kemble (2006)}, who deals with the Lambert problem directly from a geometric point of view.

2. PRELIMINARIES

The two-body problem (Keplerian) described in ^{Bate (1971)} was considered. The equation

models the dynamic of the system, with µ = GM . The position of the particle on an ellipse is r→(t) = [x₁(t) x ₂(t)], such that

where a is the semi-major axis, ǫ the eccentricity, E the eccentric anomaly. Using the fact that the norm of the angular momentum of the particle is constant we get

and therefore Kepler’s equation is obtained:

where M is the mean anomaly (^{Bate 1971}). To find a solution E for Kepler’s equation (2) at any time normally Newton’s method is used, but it does not converge or converges too slowly when ǫ ≈ 1 (^{Elipe et al. 2017}). From Danby (1962, p. 168) we use Taylor’s expansion in (2) to derive the equations of the dynamics related to the problem for any type of trajectory (classical curves). Then

With x = aE, we obtain

This equation accepts any value of 𝜖 and ȧ ≠ 0. Therefore, using the fact r = a(1 − ȧ 𝜖 cos E), (^{Bate 1971, p. 187}), and with the expression of Ė, an universal variable, x ∈ ℝ is defined as ẋ:= µr .

So, the following expressions are obtained for t and r (^{Bate 1971}) using equation (1), integrating the universal variable as µdt = rdx:

Since the initial position and velocity vectors r→₀ and v→₀ are linearly independent, and having r→, and r→ = r→ in the same plane, the vectors of position and velocity are expressed in terms of the initial vectors and the universal variable (^{Bate 1971}):

3. PROBLEM STATEMENT

The target travels in its own orbit, and the interceptor orbits around the Earth. To calculate the orbit to be taken by the interceptor to fly by the target the initial position and initial velocity of the interceptor r→₀, v→₀ and the target R→₀, V→₀ are required. We look for the minimum change of velocity Δv→₀, v→_transfer = v→₀ + Δv→₀, allowing overflight to be achieved. See Figure 1.

Fig. 1 Description of the first problem. 

Using equation (3), let A and X be the semi-major axis and the universal variable, respectively, associated with R→, and a and x the semi-major axis and the universal variable, respectively, associated with r→ (in the orbit transfer). The following functions modelling the movement of the two bodies in the solar system are defined:

We call η→_s: = (η₁ η ₂). Let R→ and r→be the positions of the target and the interceptor, respectively, shown in equation (4) and let the functional J = 1/2 Δv→₀ ^T Δv→₀. Then we have the following optimization problem:

To solve the problem, Lagrange multipliers are used in the following functional:

where λ→ ∈ ℝ² and ϕ→ ∈ ℝ ³.

After setting ∇H_s = 0 as ^{Leeghim (2013)}, we obtain the following system of nonlinear equations:

4. WAIT TIME

From ^{Vallado (1997, p. 318)} the wait time for planar and circular orbits is given by

where θ_i and θ are the initial and final angles (after t₁) between R→ and r→ , and w and W are the angular speeds of the interceptor and the target, respectively. In this paper we calculate t₁ (the wait time), which is one of the results of the problem to be minimized for any kind of transfer orbit (in space) and maintains the same functional to be minimized. See Figure 2.

Let η→_c = (η₁ η₂ η₃) such that:

Fig. 2 Description of the second problem. 

where η₁ is the equation of motion of the target, η₂ the equation of motion for the interceptor in its initial orbit during the wait time, and η₃ is the equation of motion of the interceptor after the wait time. Then we have a problem similar to the previous one:

Where J =1/2 Δv→₁ ^T Δv→₁. Therefore, the new functional is:

where λ→, ϕ→ ∈ ℝ³.

The system we obtain after setting ∇H_c = 0 as ^{Leeghim (2013)} is:

Fig. 3 Particle Swarm Optimization. 

5. SOLUTION METHOD

To find the roots of the systems (5) and (6) we use Newton’s multivariate method (^{Bate 1971}; ^{Leeghim 2013}). That is, we want to calculate f→ (x) = 0 using

where Jf→is the Jacobian matrix of f→. The variables of the first problem are y→_s = (X, x, t, Δv→₀ ⃗{ϕ}) ∈ ℝ⁹ and the second one has y→_c = (X, x_1, x, t₁, t, Δv→₁ Δϕ→) ∈ ℝ¹¹.

For each case presented in § 6, it was necessary to find a seed (an initial value for Newton’s method). However, this seed must be very close to the solution for the method to converge (Local Convergence Theorem). Given the importance of calculating the seed, we used the following heuristic method to ap-proximate the solutions, which gave us very good results.

6. RESULTS

The methods were applied to the following problem, proposed in ^{Vallado (1997, p. 352)}, which deals with orbit transfers near Earth (Low Earth Orbit, LEO).

The Hubble telescope will be released from the space shuttle, which is in a circular orbit at 590 km from the Earth’s surface. The relative ejection speed (viewed from the shuttle) is [−0 1 −0 4 −0 2]^T m/s. After 4 minutes, the Hubble needs to meet the shuttle. The change of RSW to IJK coordinates is shown in ^{Vallado (1997, p. 367)}.

After applying the algorithm, the optimal orbit transfer without wait time has the parameters of the first row of Table 2 (for LEO), and has universal variables X = 136 206 and x = 136 207. With a wait time (shown by an asterisk in Table 2), X = 651 694, x₁ = 265 476 and x = 386 215 are obtained. The orbits are shown in Figure 4.

Fig 4. Transfers between Hubble (blue) - Shuttle (black). The color figure can be viewed online. 

The Hohmann transfer gives the optimal change of velocity in planar and circular orbits. Leeghim’s method deals with more general orbits: elliptic and hyperbolic; this method gives a smaller flight time.

It is also seen in Table 2 that for LEO orbits the time of flight of the transfer with wait time decreases 32 1% with respect to the one of Hohmann. In addition, Hohmann’s change in velocity is 33 3% higher than that with wait time, the latter evidencing significant fuel savings.

The following examples are the orbit transfers from Earth to the other planets of the solar system with a launch date on January 1^st, 2014.

Figure 5 shows the current positions of Earth and Mercury for the assumed launch date, and the subsequent flight times. The resulting wait time is 3 days; with this, the change of velocity is reduced by one third and the flight time is reduced by 3 days with respect to the calculation without wait time.

Fig. 5 Transfers between Earth (blue) - Mercury (black). The color figure can be viewed online. 

These results are better than the values obtained in ^{Kemble (2006, p. 50)}, where different launch dates for the optimization were assumed.

Computations were done without taking into account technical constraints. For example, the fastest space probe that NASA has launched is New Horizons with a |Δv→| of 16 26 km/s relative to Earth, so the first |Δv→| obtained without wait time is not feasible for such a mission, even though the flight time is the shortest.

Sometimes the direction of the Earth’s velocity [−29 61 − 6 41 0]^T for January 1^st, 2014) and the initial positions of the planets are not suitable for orbit transfers. This means that, in some cases, a very large change of velocity (which translates to a large fuel consumption) is required, as presented in the first row of Table 2 (for Venus). See Figure 6.

Fig. 6 Transfers between Earth (blue) - Venus (black). The color figure can be viewed online. 

Adding the wait time (167 days) results in a much lower change of velocity (12% of the first one, which requires less fuel), and the flight time increases by 57%, improving the results presented in ^{Kemble (2006, p. 51)} for comparable parameters.

The transfer of Earth-Mars orbit was calculated with wait time and without wait time. The difference lies in the fact that, with a wait time, the change of velocity was better and the flight time increased by a week. The results obtained are better than those of ^{Kemble (2006, p. 53)}. See Figure 7.

Fig. 7 Transfers between Earth (blue) - Mars (black). The color figure can be viewed online. 

The system proposed in § 5.1 can have different solutions (local maximums of the index of performance J). The algorithm used to solve it finds many of them; here the optimal are shown.

For the example Earth-Jupiter, Figure 8, the three changes in velocity without a wait time, with wait time and as found in ^{Kemble (2006, p. 55)} are very similar, but the time of flight with wait time is smaller than the one without, and it is very close to the values obtained in ^{Kemble (2006, p. 55)}. In this case, the eccentricity of the orbit transfers found is greater than the one of the previous cases (e = 0 75).

Fig. 8 Transfers between Earth (blue) - Jupiter (black). The color figure can be viewed online. 

The solution without a wait time for the system Earth-Saturn, Figure 9, is on a parabola (with eccentricity equal to 1, an uncommon result) and |Δ v→| is extremely large. The transfer with wait time has an eccentricity of 0.82, giving a better change of velocity than the one presented in ^{Kemble (2006, p. 56)} and the flight time is 3 years shorter.

Fig. 9 Transfers between Earth (blue) - Saturn (black). The color figure can be viewed online. 

7. CONCLUSIONS

The two techniques presented here make it possible to optimize the fuel consumption. However, its use is purely theoretical, leaving aside the technical constraints, but considering only constraints of the trajectory to make a more general model. For Newton’s multivariate method, the search of the starting point was performed with the help of the PSO with excellent results. It perfectly bounds the search region of the starting point. The PSO is an heuristic method that needs no derivatives and has a rapid convergence. As future work we propose to use this method it in other situations and to perform the con-vergence analysis (local and semi-local) of the methods.

The authors are grateful with the young researcher fellowship of the Sergio Arboleda University from which this paper is one of its results.