Romain Serra
email: serra.romain@gmail.com

Denis Arzelier
email: arzelier@laas.fr

Florent Bréhart
email: florent.brehard@laas.fr

Mioara Joldes
email: mioara.joldes@laas.fr

FUEL-OPTIMAL IMPULSIVE FIXED-TIME TRAJECTORIES IN THE LINEARIZED CIRCULAR RESTRICTED 3-BODY PROBLEM

The problem of fixed-time fuel-optimal trajectories with high-thrust propulsion in the vicinity of a Lagrange point is tackled via the linear version of the primer vector theory. More precisely, the proximity to a Lagrange point i.e. any equilibrium point -stable or not -in the circular restricted three-body problem allows for a linearization of the dynamics. Furthermore, it is assumed that the spacecraft has ungimbaled thrusters, leading to a formulation of the cost function with the 1-norm for space coordinates, even though a generalization exists for steerable thrust and the 2-norm. In this context, the primer vector theory gives necessary and sufficient optimality conditions for admissible solutions to two-value boundary problems. Similarly to the case of rendezvous in the restricted two-body problem, the in-plane and out-of-plane trajectories being uncoupled, they can be treated independently. As a matter of fact, the out-of-plane dynamics is simple enough for the optimal control problem to be solved analytically via this indirect approach. As for the in-plane dynamics, the primer vector solution of the so-called primal problem is derived by solving a hierarchy of linear programs, as proposed recently for the aforementioned rendezvous. The optimal thrusting strategy is then numerically obtained from the necessary and sufficient conditions. Finally, in-plane and out-of-plane control laws are combined to form the complete 3-D fuel-optimal solution. Results are compared to the direct approach that consists in working on a discrete set of times in order to perform optimization in finite dimension. Examples are provided near various Lagrange points in the Sun-Earth and Earth-Moon systems, hinting at the extensive span of possible applications of this technique in station-keeping as well as mission analysis, for instance when connecting manifolds to achieve escape or capture.

I. INTRODUCTION

Lagrange points, a.k.a. libration points, are interesting in binary systems e.g. Earth-Moon, since they are motionless in the corotating frame, and thus offer continuous visibility over one or two of the main bodies. A favoured type of orbit around Lagrange points is the so-called halo one, due to its periodicity. Depending on the binary system, such orbits can be dynamically unstable, creating a need for station-keeping. Even in a stable configuration, orbital perturbations such as radiation pressure can be enough to require trajectory corrections. Hence the minimization of fuel consumption for trajectories in the vicinity of libration points is a topic of high interest. This problem can be formulated as an optimal control problem [START_REF] Davis | Optimal transfers between unstable periodic orbits using invariant manifolds[END_REF][START_REF] Hiday-Johnston | Transfers between libration-point orbits in the elliptic restricted problem[END_REF] and tack-IAC-18-C1,IP,32,x42624 led by a direct or indirect approach. Generally speaking, the direct approach consists in discretizing the continuous control variables into a finite number of parameters that can be fed to a nonlinear programming solver [START_REF] Hull | Conversion of optimal control problems into parameter optimization problems[END_REF] . On the other hand, the indirect approach builds on the nature of the problem and aims at guaranteing the optimality of the solution by looking at the so-called co-state variables, a somewhat generalisation of the Lagrange multipliers [START_REF] Betts | Survey of Numerical Methods for Trajectory Optimization[END_REF] . For space trajectories with high-thrust propulsion and parametrized with Cartesian coordinates, the indirect approach is best known as the primer vector theory, due to the work of Lawden in the early 1960's [START_REF] Lawden | Optimal trajectories for space navigation[END_REF] . In short, it states optimality conditions for a trajectory in terms of a vector (primer vector) whose magnitude has to be always smaller or equal to one for optimality certification. The dates or locations where it reaches unit norm are the ones where the optimal maneuvers should occur. While these optimality conditions are only necessary with nonlinear dynamics, they become sufficient in a linear setting, as proven rigorously by Neustadt [START_REF] Neustadt | Optimization, a moment problem, and nonlinear programming[END_REF] . The general, non-linear theory has already been applied to high-thrust trajectories in the restricted 3-body problem, see for instance the work of Hiday-Johnston & Howell 8 or Davis et al [START_REF] Davis | Optimal transfers between unstable periodic orbits using invariant manifolds[END_REF] . Despite their more general framework, the optimality conditions in this formulation lack the sufficiency of the linearized dynamics that are yet available in the vicinity of a Lagrange point. Furthermore, the optimal solution to this control problem under the linear approximation, along with its co-state, could be a good initial guess for the nonlinear case to be tackled indirectly. In this work, a modern numerical approach, based on convergent discretization methods for the solution of semi-infinite convex problem (SICP) [START_REF] Reemtsen | Numerical techniques for semi-infinite programming: A survey[END_REF] and heuristically improved from the one used for the restricted 2-body problem in [START_REF] Arzelier | Lin-IAC-18-C1,IP,32,x42624 earized impulsive fixed-time fuel-optimal space rendezvous: A new numerical approach[END_REF] , is introduced for trajectories near any of the libration points in the circular restricted 3-body problem (CR3BP). In order to speed up the procedure, it is coupled with recent analytical results in the linear theory of the primer vector [START_REF] Serra | Analytical solutions for impulsive elliptic out-ofplane rendezvous problem via primer vector theory[END_REF] .

Notations: The velocity increment at ν i will be denoted by ∆V (ν i ).

{b i } i=1,•••,N is a sequence of vari- ables b i , i = 1, • • • , N ,
and sgn(z) is the sign function of the variable z. The prime denotes differentiation with respect to the true anomaly ν. O p×m and 1 m denote respectively the null matrix of dimensions p × m and the identity matrix of dimension m. Let r ∈ N * and (p, q) ∈ (R × {∞}) 2 such that:

1 ≤ p ≤ ∞ and 1 p + 1 q = 1. Classically, C([ν 0 , ν f ], R r ) is the Ba- nach space of continuous functions f : [ν 0 , ν f ] → R r equipped with the norm f q = sup ν0≤ν≤ν f f (ν) q . Denote by L 1,p ([ν 0 , ν f ], R r
) the normed linear space of Lebesgue integrable functions from [ν 0 , ν f ] to R r with the norm given by: u

1,p = ν f ν0 u(ν) p dν. Let BV([ν 0 , ν f ], R r )
be the space of functions of bounded variation over the in-

terval [ν 0 , ν f ] with: g tv,p = sup Pκ κ i=1 g(ν i )-g(ν i-1 ) p ,
where the supremum is taken over all finite partitions

P κ = (ν i ) i=1,...,κ of [ν 0 , ν f ].
For a symmetric real matrix S ∈ R n×n , the notation S 0 (S 0) stands for the negative (positive) semi-definiteness of S. Finally, χ A is the indicator function of the set A.

II. LINEARIZED DYNAMICS OF THE CR3BP

The circular restricted 3-body problem deals with the motion of a negligible mass m (a spacecraft for instance) affected by Keplerian gravity field of two primary massive bodies m 2 < m 1 orbiting in a circular orbit around their barycentre B. This relative motion is classically studied in a non-inertial synodic frame (B, i, j, k) rotating with the two major point-masses m 1 and m 2 which origin is located at B. The x-axis is directed from the primary m 1 to the primary m 2 , the z-axis is parallel to the direction of the orbital angular velocity vector of the primary system with respect to an inertial frame and the y-axis completes the right-handed system as depicted on Figure 1. The mean motion ω of the primary system is constant and given by:

ω 2 R 3 = G(m 1 + m 2 ), [ 1
]
where G is the gravitational constant and R the distance between m 1 and m 2 .

R L 3 L 1 L 2 L 4 L 5 i j k x x 1 x 2 m 2 m 1 m ρ 2 ρ 1 B y
µ = m 2 m 1 + m 2 , [ 2]
IAC-18-C1,IP,32,x42624

as the mass ratio of the primary system (also defined as the non dimensional distance of m 2 to B), it is readily obtained that:

x 1 = - Rm 2 m 1 + m 2 = -Rµ, [3]
x 2 = Rm 1 m 1 + m 2 = R(1 -µ). [4]
Based on a Lagrangian formulation of the CR3TBP and changing the independent variable from the time t to the true anomaly ν of the primary system, the non-linear equations of motion of the third body m expressed in cartesian coordinates (x, y, z) in the rotating frame are given by [START_REF] Farquhar | The control and use of librationpoint satellites[END_REF][START_REF] Szebehely | Theory of orbit: The restricted problem of three Bodies[END_REF] :

x ′′ = - R 3 (1 -µ)(x + Rµ) ρ 1 3 - R 3 µ(x -R(1 -µ)) ρ 2 3 +2y ′ + x, [5]
y ′′ = - (1 -µ)R 3 y ρ 1 3 - µR 3 y ρ 2 3 -2x ′ + y, [6]
z ′′ = - (1 -µ)R 3 z ρ 1 3 - µR 3 z ρ 2 3 , [7]
where ρ 1 and ρ 2 are the distances from the third body to respectively m 1 and m 2 , (•) ′ is the derivative with respect to the true anomaly ν and ν = ω. A normalized form (R = 1 and ω = 1) of these equations may be found in a more recent reference [START_REF] Davis | Optimal transfers between unstable periodic orbits using invariant manifolds[END_REF] .

II.ii Linearized equations around the libration points

It is well known that the dynamical system represented by Equations [5]- [7] exhibits five equilibrium points, the so-called libration points, which are stationary points with respect to the rotating system. The negligible mass m placed at these points with zero initial velocity will remain there if there is no external perturbation. All five equilibrium points lie in the xy plane and two families of libration points are distinguished. The collinear points L 1 -L 3 for which y = 0 and z = 0, satisfy the following equation:

x = (1 -µ)R 3 (x + Rµ) 2 + µR 3 (x -R(1 -µ)) 2 .
[8]

The triangular or equilateral points L 4 , L 5 are defined by z = 0 and:

x = R 2 -µR, y = ± √ 3 2 R. [9]
A typical geometrical configuration of all libration points is given in Figure 1.

Defining the small perturbation of a nearby trajectory around an equilibrium point L i , where i ∈ {1, . . . , 5}, the linearized autonomous equations of the relative motion of m around the Libration point L i are then given in statespace form as: X(ν) ′ = AX(ν), [10] where X T = δx δy δz δx ′ δy ′ δz ′ and

A =         0 0 01 0 0 0 0 0 0 1 0 0 0 0 0 0 1 1 -U xx -U xy 0 0 2 0 -U xy 1 -U yy 0 -2 0 0 0 0 -U zz 0 0 0         , [11]
U zz = R 3 (1 -µ) r 3 1Li + R 3 µ r 3 2Li , [ 12
]
U xx = U zz - 3(1 -µ)R 3 (x Li -R(1 -µ)) 2 r 5 1Li - 3µR(x Li + Rµ) 2 r 5 2Li , [13]
U xy = - 3(1 -µ)R 3 (x Li -R(1 -µ))y Li r 5 1Li - 3µR 3 (x Li + Rµ)y Li r 5 2Li , [14]
U yy = U zz - 3(1 -µ)R 3 y 2 Li r 5 1Li - 3µR 3 y 2 Li r 5 2Li
, [15] and

r 1Li = (x Li + Rµ) 2 + y 2 Li , [16]
r 2Li = (x Li -R(1 -µ)) 2 + y 2 Li . [17]
The dynamical matrix A is a time-invariant matrix and its associated transition matrix Φ(ν, ν 0 ) for Equation [11] is given by exp(A(νν 0 )). Note that the xy inplane and out-of-plane motions are completely decoupled. Since U zz > 0, the out-of-plane motion is an harmonic oscillator with pulsation √ U zz while the in-plane trajectories may be bounded or not depending on µ and the initial conditions.

III. A PRIMER VECTOR APPROACH

The primer vector theory is the historical name given to an indirect approach to the design of fuel-optimal space trajectories. The seminal results dates back to the fifties and the beginning of the sixties with the works of Lawden [START_REF] Lawden | Impulsive transfer between elliptical orbits[END_REF][START_REF] Lawden | Optimal trajectories for space navigation[END_REF] which have been rigorously confirmed by Neustadt in linear [START_REF] Neustadt | Optimization, a moment problem, and nonlinear programming[END_REF] and nonlinear 18 settings right after. In this paper, since we are interested by designing minimumfuel transfer trajectories around the libration points, only its linear version is used and described. The main lines of the review of the theoretical background of the proposed results are based on the reference [START_REF] Arzelier | Lin-IAC-18-C1,IP,32,x42624 earized impulsive fixed-time fuel-optimal space rendezvous: A new numerical approach[END_REF] .

III.i Problem statement and its moment counterpart

Designing minimum-fuel transfer trajectories around the libration points in the CR3BP setting naturally boils down to solving the following linear optimal control problem:

Problem 3.1. (Optimal control problem) Find ū ∈ L 1,p ([ν 0 , ν f ], R 3 ) solution of the optimal con- trol problem: inf u u 1,p = inf u ν f ν0 u(ν) p dν s.t. X ′ (ν) = AX(ν) + Bu(ν), X(ν 0 ) = X 0 , X(ν f ) = X f ∈ R n , ∀ ν ∈ [ν 0 , ν f ], ν 0 , ν f fixed, [18]
where B T = O (n-r)×r 1 r . The state vector dimension and the number of inputs in [10] are denoted n and r, respectively with n = 2, r = 1 for the out-of-plane case and n = 4, r = 2 for the in-plane case. The value of the integer p actually depends on the geometric configuration of the thrusters. For a steerable thruster, p = 2, while for 6 ungimbaled identical thrusters, p = 1.

Following the approach from 17 , Problem 3.1 is now transformed into an equivalent problem of moment by integrating equation [18]. Let us define a fundamental matrix ϕ(ν) for Equation [10]

(i.e. Φ(ν, ν 0 ) = ϕ(ν)ϕ -1 (ν 0 )) and the matrix Y (ν) = ϕ -1 (ν)B = y 1 (ν) • • • y n (ν) T ∈ R n×r , then: c = ϕ -1 (ν f )X(ν f ) -ϕ -1 (ν 0 )X 0 = ν f ν0 ϕ -1 (σ)B(σ)u(σ)dσ = ν f ν0 Y (σ)u(σ)dσ. [19]
It is important to notice for the remainder of the analysis that for the specific matrices Y (ν) encountered in the studied problem,

y 1 (ν) • • • y n (ν) are linearly independent elements of C([ν 0 , ν f ], R r
). This will be assumed in the rest of the paper. It follows from [19] that Problem 3.1 can be equivalently written as:

Problem 3.2. (Minimum norm moment problem) Find ū(t) ∈ L 1,p ([ν 0 , ν f ], R r ) solution of the minimum norm moment problem: inf u u 1,p = inf u ν f ν0 u(ν) p dν s.t. ν f ν0 Y (σ)u(σ)dσ = c, ν 0 , ν f fixed.
[20]

It is well-known that Problem 3.2 may not reach its optimal solution due to concentration effects (see the reference [START_REF] Görner | Numerical techniques in relaxed optimization problems[END_REF] ) . It is then necessary to resort to a relaxation scheme by embedding the space

L 1,p ([ν 0 , ν f ], R r ) in the dual space C * ([ν 0 , ν f ], R r ) of the Banach space C([ν 0 , ν f ], R r ).

III.ii A relaxed problem and its SICP formulation

Generalized solutions of the original Problem 3.2 may be obtained as the solutions of the following relaxed moment problem. [START_REF] Neustadt | Optimization, a moment problem, and nonlinear programming[END_REF] that the infimum of Problem 3.3 is reached and that it is equal to the infimum of Problem 3.2, denoted by η in what follows. The association between the space BV([ν 0 , ν f ], R r ) and the dual

P κ = {ν 0 = ν 1 < ν 2 , • • • , < ν κ = ν f } is any finite partition of [ν 0 , ν f ]. It is shown in
C * ([ν 0 , ν f ], R r ) of the space C([ν 0 , ν f ], R r
) defined in the Riesz Representation Theorem and the use of a duality principle based on the extension form of the Hahn-Banach theorem [START_REF] Luenberger | Optimization by Vector Space Methods[END_REF] are the basic elements for the next result, originally given in [START_REF] Neustadt | Optimization, a moment problem, and nonlinear programming[END_REF] in its complete form and partially in [START_REF] Krasovskii | On the theory of optimum regulation[END_REF] for particular optimization problems. Here, we follow the lines developed in the textbook of [16, Chapter 5] .

Theorem 1. Let y i (•) ∈ C ([ν 0 , ν f ], R r ), ∀ i = 1, • • • , n, q = p
p-1 the conjugate index and suppose that

D = g ∈ BV([ν 0 , ν f ], R r ) : ν f ν0 Y (ν)dg(ν) = c , [ 22
] is a non empty set then η = min g∈D g tv,p = max

Y T (ν)λ q ≤1 c T λ. [23]
IAC-18-C1,IP,32,x42624

In addition, let ḡ and λ be optimal solutions of [23],

λ = Arg[ max Y T (ν)λ q ≤1 c T λ] and let ȳ(ν) = n i=1 λi y i (ν) = Y T (ν) λ ∈ R r .
Then the optimal ḡ is such that:

ν f ν0 λT Y (ν)dḡ(ν) = sup ν0≤ν≤ν f ȳ(•) q ḡ tv,p [24]
The two problems defined in Eq. [23] may be considered as dual through the equality of the optimal values of their respective objectives and the relation between their solutions thanks to the condition in Eq. [24]. This results in a significant simplification: The infinite-dimensional optimization Problem 3.3 has been converted to a search of an optimal vector λ in a finite-dimensional vector space submitted to a continuum of constraints, yielding a semiinfinite convex problem (SICP):

Problem 3.4. (SICP problem) Find λ ∈ R n solution of μ = min λ∈R n -c T λ Y T (ν)λ q ≤ 1.
[25]

Note that μ = -η. Once its solution is obtained, the relation [24] between the function ȳ(•) element of the Banach space C([ν 0 , ν f ], R r ) and the optimal vector λ is particularly important to get back to the optimal bounded variation solution of the relaxed Problem 3.3 as shown in the following result from [START_REF] Neustadt | Optimization, a moment problem, and nonlinear programming[END_REF] . Theorem 2. Let y i (•) ∈ C ([ν 0 , ν f ], R r ), i = 1, . . . , n and λ ∈ R n be an optimal solution of Problem (25). Define the sets

Γ s = {ν ∈ [ν 0 , ν f ] : |ȳ s (ν)|= 1} and Γ = ν ∈ [ν 0 , ν f ], ȳ(ν) q = max ν0≤ν≤ν f ȳ(ν) q = 1 . Note that Γ = ∪ s Γ s for p = 1.
There is an optimal solution ḡ(•) ∈ BV ([ν 0 , ν f ], R r ) of the relaxed Problem 3.3, which is a step function with at most n points of discontinuity νj ∈ Γ, j = 1, • • • , N ≤ n. Its jumps are given by: ḡs (ν j ) -ḡs (ν - j ) = α νj sgn(ȳ s (ν j ))χ Γj , α νj > 0, when p = 1, or ḡs (ν j ) -ḡs (ν - j ) = α νj |ȳ s (ν j )| q-1 sgn(ȳ s (ν j )), when 1 < p < ∞,

[26] for s = 1, • • • , r and α νj solutions of the linear system:

N j=1 β i (ν j )α νj = c i , i = 1, • • • , n [27]
where β i (ν j ) are given by:

β i (ν j ) = r s=1 y i,s (ν j )sgn(ȳ s (ν j )), when p = 1, or β i (ν j ) = r s=1 y i,s (ν j )|ȳ s (ν j )| q-1 sgn(ȳ s (ν j )), when 1 < p < ∞, [28] for all j = 1, • • • , N .
This theorem states important results that have been known for a while in the aerospace community but whose value has not been completely exploited to derive efficient numerical algorithms for impulsive maneuvers design. First, it says that the optimal controlled trajectory for the minimum-fuel linearized optimal control problem is purely impulsive and that the number of impulses is upper-limited by n which is the dimension of the fixed final conditions of the optimal control problem.

Remark 1.

It is also shown in [START_REF] Neustadt | Optimization, a moment problem, and nonlinear programming[END_REF] that a sequence of functions u ǫ (•) ∈ L 1,p ([ν 0 , ν f ], R r ) converges to a linear combination of δ(•) functions corresponding to the function ḡ(•) with equal norms. Let ∆V (ν j ) = ḡ(ν j ) -ḡ(ν - j ), then roughly speaking, this may be described by:

ūǫ (ν) → N j=1
∆V (ν j )δ(ν jν), ǫ → 0.

[29] Indeed, the initial optimal control problem amounts to find the sequences of optimal impulse locations {ν i } i=1,•••,N and optimal impulse vectors {∆V (ν i )} i=1,•••,N verifying the boundary equation:

c = N i=1 Y (ν i )∆V (ν i ).
[30]

The vector p v (ν) = Y T (ν)λ involved in [25] is nothing but the primer vector initially defined in the seminal work of Lawden [START_REF] Lawden | Optimal trajectories for space navigation[END_REF] . In this reference, the primer vector p v (ν) is defined as the velocity adjoint vector arising from applying the Pontryagin Maximum Principle to optimal trajectory problems or Lagrangian duality as in [START_REF] Carter | Linearized impulsive rendezvous problem[END_REF] where the vector λ is the optimal Lagrange multiplier. For an optimal impulsive trajectory, the primer vector p v (ν) must satisfy the well-known Lawden's necessary and sufficient optimality conditions recalled in [START_REF] Carter | Linearized impulsive rendezvous problem[END_REF] . In summary, finding the optimal solution to the primal problem is an indirect way of solving the original trajectory problem. As a result, efficient algorithms to solve the former are extremely valuable. The next section proposes a new procedure based on a discretization algorithm for the solution of the semi-infinite programming Problem 3.4.

IAC-18-C1,IP,32,x42624

IV. A MODERN NUMERICAL APPROACH

The resolution of Problem 3.4 when dealing with trajectories around Lagrange points is analogous to the solution of the circular linearized rendezvous 2-body problem [START_REF] Arzelier | Lin-IAC-18-C1,IP,32,x42624 earized impulsive fixed-time fuel-optimal space rendezvous: A new numerical approach[END_REF] . Actually, the out-of-plane dynamics being a harmonic oscillator as well, the optimal solution for n = 1 is analytical and identical to the one given by Prussing [START_REF] Prussing | Illustration of the primer vector in time-fixed orbit transfer[END_REF] . As for n ≥ 2, it can be solved numerically as in Arzelier et al [START_REF] Arzelier | Lin-IAC-18-C1,IP,32,x42624 earized impulsive fixed-time fuel-optimal space rendezvous: A new numerical approach[END_REF] . Basically, the primer vector is obtained iteratively by solving a sequence of convex problems in R 2n .

An improved heuristic version is proposed in this paper. Moreover, by relying on properties of the 1 and ∞ norms, concatenating the projected in-plane and out-of-plane optimal solution yields the optimal solution for n = 3 when p = 1. Taking advantage of this property makes for a faster resolution under the complete dynamics. This section details the main results at the basis of the algorithm used to compute the minimum-fuel sequence of impulsive maneuvers.

The general, numerical approach proposed in 1 and based on discretization methods or exchange methods reminds the numerical solution of the Chebyshev approximation problem [START_REF] Reemtsen | Numerical techniques for semi-infinite programming: A survey[END_REF] . It intends as solving a series of convex problems 3.4 for which the feasible set described by infinitely many constraints is replaced by a discrete approximation involving a finite set of constraints defined on a given grid Θ i of locations in the interval [ν 0 , ν f ]. For p = 1, the discretized version of 3.4 boils down to a linear program:

Problem 4.1. Discretized SICP (p = 1) max λ c T • λ s.t. -1 ≤ |Y T • λ(ν)| j ≤ 1, ∀ ν ∈ Θ i , j = 1, . . . , n,
while for p = 2, it falls into the scope of Semi-Definite Programming :

Problem 4.2. Discretized SICP (p = 2) max λ c T • λ s.t. 1 Y T • λ(ν) Y • λ T (ν) 1 n ≥ 0 ∀ν ∈ Θ i .
The key point in the method is that, by properly building a sequence of sets Θ i , the solutions λ * (i) of the discretized problems 4.1 or 4.2 will converge to the solution of 3.4. Its efficiency strongly relies on the fact that Problems 4.1 and 4.2 can be solved efficiently by dedicated methods like simplex or interior-point methods [START_REF] Ben-Tal | Lectures on Modern Convex Optimization[END_REF] .

It has been shown in [START_REF] Arzelier | Lin-IAC-18-C1,IP,32,x42624 earized impulsive fixed-time fuel-optimal space rendezvous: A new numerical approach[END_REF] that the sequence of discretized sets

Θ i+1 = Θ i {argmax ν∈[ν0,ν f ] Y T • λ * (i) (ν) q }
verifies the necessary assumptions of 1 and converges to the optimal solution of 3.4. In this paper, following incentives from exchange methods [START_REF] Reemtsen | Numerical techniques for semi-infinite programming: A survey[END_REF] and ideas from [START_REF] Arzelier | A new mixed iterative algorithm to solve the fuel-optimal linear impulsive rendezvous problem[END_REF] , a heuristic improvement is made by also removing all dates whose corresponding constraint is not saturated. In other words, [START_REF] Koenig | Real-time algorithm for globally optimal impulsive control of linear timevariant systems[END_REF] , this modification leads to a reduction of the number of constraints involved in Problems 4.1 or 4.2. This is summarized in the algorithm below written in pseudo-code.

Θ i+1 = {ν ∈ Θ i , Y T • λ(ν) q = 1} {argmax ν∈[ν0,ν f ] Y T • λ * (i) (ν) q }. Also used in
Input: interval Θ = [ν 0 , ν f ], matrix Y (ν), initial condition c, accuracy ε Output: µ (i) and λ (i) numerical solution of Pb. 3.4 Init: i ← 0; Convergence ← F alse; Θ 0 ← {θ 0 ; θ 1 } ⊂ Θ s.t. θ 0 -θ 1 = kπ; while Convergence=False do Find λ (i)
solution of discretized problem:

µ (i) = inf λ∈R n -c T λ s.t. Y T (θ k )λ q ≤ 1 for all θ k ∈ Θ i Convergence = max θ∈Θ Y (θ) T λ (i) q -1 > ε if Convergence=False then i ← i + 1; Θ i ← Θ i-1 ∪ arg max θ∈Θ Y T (θ)λ (i) q ; Θ i ← Θ i -θ k | Y (θ k ) T λ (i) q < 1 -ε ; end end return µ (i) , λ (i) .
Algorithm 1: Numerical procedure for solving Problem 3.4

Once Algorithm 1 returns the near-optimal value of λ, which gives the near-optimal primer vector, Theorem 2 is used to reconstruct the instantaneous velocity increments.

V. RESULTS AND DISCUSSION

The algorithm described in the previous Section has been implemented in Python, making use of the scipy.linprog routine for Linear Programming and cvx solver for Semi-Definite Programming. The practical results presented thereafter concern the L 2 libration point of the Earth-Moon system which is of particular interest as the tentative location of a fuel-depot to be used by interplanetary spacecraft. The method is nonetheless applicable to any Lagrange point in any binary system. A set of initial and final values for the transfer (in-and out-IAC-18-C1,IP,32,x42624 of-plane), reported in Tables 1, was generated randomly, but could be chosen for specific purposes such as stationkeeping or transfer between Halo orbits. Comparisons are twofold: between 1 and 2-norm costs on the one hand and between indirect and direct approaches on the other hand. For an a priori fixed number of impulsive maneuvers at given locations, a convex problem (LP or SDP) is formulated and solved numerically by discretizing the possible locations of burns. Its solution is therefore suboptimal, depending strongly upon the number of impulses [START_REF] Larsson | Fuel efficient relative orbit control stategies for formation flying and rendezvous within prisma[END_REF][START_REF] Louembet | Robust rendezvous planning under maneuver execution errors[END_REF] . Introducing slack variables, 1-and 2-norm minimization respectively boils down to a linear and semi-definite program. Due to the greater computational cost of the latter, the grid is chosen to be sparser here e.g. 1000 points against 10000. On the other hand, the indirect approach uses the thinner version when checking the magnitude of the primer vector. Figure 3 depicts the primer vector history for p = 1 whereas Figure 2 is for the fuel-consumption, the relative positions and velocities obtained via the indirect approach. The corresponding trajectory is showed in Figure 4. The same plots for p = 2 are respectively on Figures 6, 5 and 7. Tables 2 and3 contain the control laws respectively for p = 1 and 2. The optimality of the solutions is certified by the norm of the primer vector being always less or equal to one, the latter occurring exclusively at times of burns. control has more burns i.e. four rather than two. This can be explained by the fact that it is more efficient with p = 2 to minimize the number of maneuvers. Here one can see that somehow the two middle impulses of p = 1 are merged into a single one for p = 2. This also translates into an optimal solution for the 2-norm with no zero component, unlike for the 1-norm. Moreover, the mid-course maneuvers exhibited by both minimal control laws is a clear setback for the analytic yet generally sub-optimal 2-impulse method that consists in designing an initial maneuver to achieve transfer in position at the final location where the second and last impulse matches the desired velocity. The fuel-optimal solution obtained for p = 2 in particular shows the advantage of an initial coasting arc. Objective functions and timings are reported for direct and indirect approaches in Table 4 for the 1-norm and Table 5 for the 2-norm. For the two cases, both approaches return an identical fuel-consumption. However, due to need for a fine locations grid for the possible burns, the direct method requires in general a larger computing effort.

VI. CONCLUSIONS

In this paper, a modern approach to the linear primer vector theory has been proposed for fuel-optimal trajectories near any Lagrange point, extending the previous work of the authors with the linearized restricted elliptical 2body problem. It combines recent results that are both analytical and numerical to provide the user with a fast, efficient way of computing the optimal burns. Output of a sequence of linear or semi-definite programs (depending on the geometrical configuration of the thrusters), the minimal consumption can be evaluated without the actual need of calculating the maneuvers. The latter, which cannot be achieved by a direct approach, turns this method into a powerful tool for mission analysis. Still, the convergence of the used algorithm remains an open problem to be studied in coming developments. Future work includes also considering the elliptical case i.e. the non-circular generalization of the Lagrange points. Perspectives are to initialize an indirect solving approach to the control problem under the non-linear dynamics with the solution of the linearized version.

Fig. 1 :

 1 Fig. 1: Libration points geometry in the corotating frame

Problem 3 . 3 .

 33 (Relaxed moment problem) Determine ḡ ∈ BV([ν 0 , ν f ], R r ) solution of the following problem: i )g(ν i-1 ) p , s.t.

ν f ν0 Y

 ν0 (ν)dg(ν) = c.[21]

Fig. 3 :Fig. 4 :

 34 Fig. 3: Optimal primer vector in the 1-norm case

Fig. 6 :Fig. 7 :

 67 Fig. 6: Optimal primer vector in the 2-norm case

Table 1 :

 1 Initial and final conditions

	ν 0 (rad) δx 0 (m)	δy 0 (m)	δz 0 (m)
	3.322	6449.40	65117.03 22814.91
		δ ẋ0 (m/s) δ ẏ0 (m/s) δ ż0 (m/s)
		-0.0312	0.0392	0.2114
	ν f (rad) δx f (m)	δy f (m)	δz f (m)
	4.737	59066.09 67728.64 84015.47
		δ ẋf (m/s) δ ẏf (m/s) δ żf (m/s)
		-0.1087	0.1616	-0.1730

Table 2 :

 2 Location and components of optimal impulses for p = 1, N * = 4

	ν (rad)	3.322	3.987	4.030	4.737
	∆V x (m/s)	0.0126 0	0	-0.5540
	∆V y (m/s)	0	0	0.5530 0.3617
	∆V z (m/s)	0	0.1570 0	0
	∆V 1 (m/s) 0.0126 0.1570 0.5530 0.9158

Table 3 :

 3 Location and components of optimal impulses for p = 2, N * = 2

	ν (rad)	3.928	4.737
	∆V x (m/s)	-0.0181 -0.5677
	∆V y (m/s)	0.5173	0.4595
	∆V z (m/s)	0.1541	0.0165
	∆V 2 (m/s) 0.5401	0.6850
	The minimal 1-and 2-norm trajectories obtained via
	the indirect approach are visually very similar. The main
	difference comes from the fact that the minimal 1-norm

Table 4 :

 4 Comparison between direct and indirect approaches for p = 1

	Approach Cost (m/s) Computation time (s)
	Indirect	1.6384	0.8447
	Direct	1.6384	1.3974

Table 5 :

 5 Comparison between direct and indirect approaches for p = 2

	Approach Cost (m/s) Computation time (s)
	Indirect	1.2251	1.3731
	Direct	1.2251	160.5454

th International Astronautical Congress, Bremen, Germany. Copyright c 2018 by the authors. Published by the International Astronautical Federation with permission.