# Preliminary Capture Trajectory Design for Europa Tomography Probe

Authors: Lorenzo Federici; Alessandro Zavoli; Guido Colasurdo
Journal: International Journal of Aerospace Engineering
 (2018)

---

## Abstract

The objective of this work is the preliminary design of a low-ΔV transfer from an initial elliptical orbit around Jupiter into a final circular orbit around the moon Europa. This type of trajectory represents an excellent opportunity for a low-cost mission to Europa, accomplished through a small orbiter, as in the proposed Europa Tomography Probe mission, a European contribution to NASA’s Europa Multiple-Flyby Mission (or Europa Clipper). The mission strategy is based on the v∞ leveraging concept, and the use of resonant orbits to exploit multiple gravity-assist from the moon. Possible sequences of resonant orbits are selected with the help of the Tisserand graph. Suitable trajectories are provided by an optimization code based on the parallel running of several differential evolution algorithms. Different solutions are finally compared in terms of propellant consumption and flight time.


---

## Body

## 1. Introduction

The Jovian moon Europa is a celestial body of primary interest for astrophysicists. The likely existence of a global subsurface ocean, proved by measurements carried out during Galileo mission, makes Europa one of the most promising environments in the Solar System to sustain human habitability. The presence of an ocean may also imply that Europa hosts (or, at least, hosted) life [1]. The importance of the determination of the ice-water layer characteristics is clearly stated in NASA’s 2013–2022 Decadal Survey [2].Europa Clipper is the next mission planned by NASA with the aim of exploring Europa. Because of the extremely harsh Jovian environment in the proximity of Europa, the initial concept of an orbiter was abandoned in favour of a multi-flyby strategy, the same considered for Galileo mission. The present mission profile, with more than 40 flybys of Europa, allows for a paramount investigation of Europa surface and subsurface properties, but is not very favourable to the investigation of Europa’s deep interior structure.A scientific enhancement to Europa Clipper mission was investigated in [3]. There, a small probe deployed on a polar orbit around Europa, hosting just one scientific instrument (a magnetometer) and a transponder required for the Intersatellite Link (ISL) with the mother spacecraft, is proved to be capable of providing crucial information on the interior structure of the moon, such as depth, thickness, and conductivity of the subsurface ocean. Also, ISL could support the reconstruction of the mother spacecraft orbit, hence significantly improving the accuracy of the topographic reconstruction of Europa’s surface.Standing on these arguments, a scientific and engineering team at Sapienza University of Rome, in collaboration with the Imperial College of London, carried out a feasibility study for a probe that could be hosted by the main spacecraft during the interplanetary cruise and released in the Jovian system with the aim at entering into a low-altitude circular quasi-polar orbit around Europa [4]. The result is a small spacecraft named Europa Tomography Probe (or ETP), which could fit the provisional 250 kg allowance that NASA has assigned to a secondary flight element hosted by the main spacecraft.The feasibility study was carried out under the design philosophy of determining the minimum total mass and volume that allows for the scientific measurements considered in [3]. All subsystems have been described with some details, with the relevant exception of only two elements: (a) the transponder and (b) the trajectory which moves ETP into a polar orbit around Europa, which should interfere as less as possible with the mother spacecraft mission plan. This paper investigates the latter point, that is, the capture strategy for an Europa orbiter at the level of preliminary mission analysis.The problem of optimizing the capture trajectory of an orbiter (or a lander) directed towards a moon of an outer planet, such as Europa [5], Enceladus [6], or Titan [7], has been the subject of many investigations. A two-body patched-conic approximation is usually assumed for interplanetary missions and transfers in a multibody planetary system [8]. This dynamical model retains the most prominent features of the real system, while keeping the numerical difficulties low. Three-dimensionality and eccentricity of the planetary bodies can be easily taken into account, and several flybys of different moons can be dealt with.The same kind of missions have been also studied by using dynamical system techniques, which rely on circular restricted three-body problem (CRTBP). Low-energy trajectories are searched for, attempting the construction of “transfer tubes,” whose boundaries are typically given by invariant manifolds originating from invariant sets (such asL1 and L2 Lyapunov orbits). As an example, dynamical chains formed by linking heteroclinic connections and homoclinic orbits [9] are proposed for the analysis of fast resonance transitions between exterior and interior resonant orbits (in the Sun-Jupiter system) [10] or “loose” capture trajectories [11]. Similar concepts are exploited for Halo-to-Halo [12] or libration-to-libration [13] transfers between planetary moons in the Jovian system, adopting a “patched” CRTBP model.In the present problem, the probe comes from a high-energy condition and approaches the moon with a high hyperbolic excess velocity. The latter techniques are thus not efficient for attaining a solution, while a patched-conic approximation can be profitably adopted.Delta velocity gravity assist (ΔV-GA) or v∞ leveraging [14] has proved a powerful concept to improve the design of capture (or escape) trajectories. Large changes of the hyperbolic excess velocity at the encounter (v∞) are obtained by using small deep-space maneuvers. When this strategy is used in conjunction with a series of resonant gravity assists, a significant reduction of propellant requirement can be achieved, with respect to a direct insertion maneuver [15].In the present paper, the design of ETP capture trajectory usingv∞ leveraging is pursued by blending the patched-conic model and a modern global optimization procedure based on a differential evolution algorithm. The trajectory is modeled as a sequence of legs between two moon encounters; only one deep-space maneuver is permitted in each leg. This approach, proposed in [16] and hereafter referred as “MGA-1DSM,” permits a quite general parameterization of the whole trajectory, which is not limited to ΔV-GA maneuvers. A preliminary solution (the sequence of resonant orbits and intercepted bodies) is defined by using two simple tools: the suboptimal solution of the v∞ leveraging problem proposed by Sims and Longuski [17] and the Tisserand graph [18]. The former permits an easy design of a mission based only on v∞ leveraging maneuvers, by suggesting a viable sequence of resonant orbits. The latter is a powerful graphical aid for the design of the same class of missions, when multiple bodies are intercepted, and some deep-space maneuvers are conveniently replaced by gravity assists of other moons in the planetary system.The paper organization is here outlined. In Section2, the physical problem of interest is described, and the adopted dynamical model and relevant assumptions are stated. A mathematical formulation (MGA-1DSM), which parameterizes a generic interplanetary trajectory as a sequence of legs containing a gravity assist and one deep-space maneuver, is outlined, leading to the purposeful definition of a global optimization problem. Section 3 presents the multipopulation differential evolution algorithm that has been used to solve the optimization problem. Fundamental tools for preliminary mission design, that is, v∞ leveraging and Tisserand graph, are discussed in Section 4. A tentative solution is devised and used to prune global optimization search. Numerical results of this investigation are presented in Section 5. A conclusion section ends the paper.


## 2. Problem Statement and Mathematical Modeling

### 2.1. Problem Overview

According to the ongoing proposal [4], the probe is assumed to be released by the main spacecraft after a few Europa flybys have been completed. In particular, ETP starts its own transfer at the apocenter of a Jovian orbit of period four times the period of Europa (TEu), pericenter equal to the Europa semi-major axis, and coplanar with the Europa orbit; the orientation of the major axis is left free. A target circular quasi-polar orbit around Europa is desired, of assigned altitude hEOI=250 km over Europa’s surface. Four R-6D bipropellant thrusters form the primary propulsion system of the probe, which allows for a total thrust of 88 N with specific impulse Isp=294 s. This propulsion system will be used for deep-space maneuvers, Europa Orbit Insertion (EOI) maneuver, and orbit maintenance during the scientific part of the mission.A probe “net” massmu=146.6 kg, which does not account for the propellant mp and tank ms masses, was estimated in [4]. Assuming a structural coefficient ϵ=ms+mp/ms=6, which is a reasonable value for liquid propellant systems, a maximum value of velocity increment ΔVmax=1240.8 m/s can be obtained if the spacecraft wet mass m0=mu+mp+ms is constrained at 250 kg. This value of ΔVmax must cover the orbit maintenance (about 43.2 m/s for a 6-month mission) and capture cost. The goal is to reduce the ΔV required for the capture, so that a convenient safety margin is left.


### 2.2. Dynamical Model

A patched-conic model is assumed for the present analysis. Flybys are modeled as instantaneous changes in velocity. Subscripts “−” and “+” are used to distinguish between values immediately before or after the discontinuity, respectively. The radius of the sphere of influence (SOI) of the secondary bodies and the travel time inside these regions are assumed to be negligible. Powered flybys are neglected, as considered useless to reduce propellant consumption [19]. An impulsive-thrust model is adopted. This assumption well suits deep space maneuvers (DSMs), which are performed at a large distance from the main body and require usually a quite short time if compared to the orbital period, as chemical engines are here considered. This assumption is also used for the EOI maneuver, even though finite-thrust losses might be considered. Only one DSM is permitted between a flyby and the other. As a further assumption, Jovian moons move on Keplerian orbits (even though the proposed procedure is soon applicable to the general case that uses planetary ephemeris).Despite its simplicity, this model allows to capture the most prominent features of the mission, while keeping the analysis simple enough. In fact, under the hypothesis of impulsive thrust, the trajectory can be computed analytically, without involving the numerical integration of the complete equations of motion.


### 2.3. Trajectory Parameterization

Let us assume that a sequenceM=Mj∣j=1,…,N+1 of N+1 body encounters, where Mj∈1..4 identifies the encountered body (1 = Io, 2 = Eu, 3 = Ga, and 4 = Ca), has been established. Europa is the first and last body in the series. The spacecraft trajectory can be modeled according to the multiple gravity assist-one deep space maneuver (MGA-1DSM) formulation [16]. The trajectory is broken down into a series of body-to-body legs. Each leg starts with a flyby and is made up of two ballistic arcs, joined by an impulsive maneuver.This general formulation for a multigravity assist trajectory can be adapted to the problem at hand by adding an initial leg, which moves the probe from the assigned initial conditions to the first encounter with Europa. The mission ends with a last approach to Europa’s surface, where an impulsive maneuver inserts the probe into the assigned polar orbit.


#### 2.3.1. Departure Leg

The departure leg, which connects ETP release position to the first encounter with Europa, is modeled as a Lambert arc, where release epocht0, flight angle Δθ, and flight time ΔT0 are design parameters to optimize.Leta^r,a^t,a^n be a Jovicentric radial-traversal-normal reference frame connected to a Europa position at epoch t1=t0+ΔT0, that is,
(1)a^r=r1r1a^n=r1×vm,1r1×vm,1a^t=a^n×a^r,where r1=rM1t1 and vm,1=vM1t1 indicates position and velocity vectors of Europa (M1=2) at time t1, respectively, which are provided by the ephemeris.The probe departure point is located on a circle of radiusr~0 which lies on the plane a^r−a^t; hence, it can be expressed as
(2)r0=r~a4:1cosΔθa^r−sinΔθa^t,while the velocity immediately before the release maneuver is
(3)v0−=v~a4:1sinΔθa^r+cosΔθa^t,where the values r~a4:1 and v~a4:1 are, respectively, the radius and velocity magnitude at the apocenter of a Jovian orbit with pericenter equal to the Europa semi-major axis, and period 4TEu. In this respect, the problem solution will eventually define the optimal orientation of the line of apses of the initial orbit with respect to Europa’s orbit.The velocity vectors after the release maneuverv0+ and immediately before the first flyby v1− can now be evaluated by solving the associated Lambert problem:
(4)Lambertr0,r1,ΔT0⟶v0+,v1−.The propulsive cost of the release maneuver is evaluated as(5)ΔV0=v0+−v0−.


#### 2.3.2. Intermediate Legs

Thekth leg can be parametrized by using four parameters rπ,k, βk, ΔTk, and ηk, which represent, respectively, the flyby radius, the flyby plane orientation, the (overall) leg flight time, and the fraction of the leg flight time at which DSM occurs. The trajectory associated to a generic intermediate leg is presented in Figure 1.Figure 1
Trajectory sketch for an intermediate leg.Beginning from timet1 at the first encounter, the epochs of the following encounters can be evaluated recursively as tk+1=tk+ΔTk, while the DSM epochs are tk+1/2=tk+ηkΔTk. Moon position rk=rMktk and velocity vm,k=vMktk are obtained from the body ephemeris. The position of the spacecraft at the flyby must be the same as the intercepted body; hence, rk−=rk+=rk.Let’s define a body-centered reference systeme^1,e^2,e^3, where e^1 is directed along the spacecraft incoming relative velocity, e^2 is orthogonal to the moon orbital plane, and e^3 closes the right-hand side term
(6)e^1=v∞k−v∞k−e^2=e^1×vm,ke^1×vm,ke^3=e^1×e^2,where v∞k=vk−−vm,k is the probe relative velocity before the flyby. The probe velocity vk+ after the flyby is evaluated as
(7)vk+=vm,k+v∞,kcosδke^1+cosβksinδke^2+sinβksinδke^3,where the rotation δk of the hyperbolic excess velocity in the flyby plane is
(8)δk=2asinμMkμMk+rπ,kv∞k2.Once the spacecraft state after the flyby is fully known, the position and velocity just before DSM, which arerk+1/2 and vk+1/2−, respectively, can be evaluated analytically by using propagation formulas for Keplerian orbits [20].Position of the spacecraft is known at both DSM maneuver and next flyby; the transfer time1−ηkΔTk is also known. Velocity vectors vk+1/2+ immediately after the DSM and vk+1− just before the subsequent flyby can be evaluated by solving the associated Lambert problem
(9)Lambertrk+1/2,rk+1,ΔTk⟶vk+1/2+,vk+1−Evaluation of the DSM propulsive cost is now straightforward(10)ΔVk=vk+1/2+−vk+1/2−


#### 2.3.3. Europa Orbit Insertion

The ETP injection into the target polar orbit around Europa is modeled as an impulsive burn applied at the hyperbola pericenter, which is purposely located on the target circular orbit.(11)vEOI−=v∞2+2μEurEOIvEOI+=μEurEOI.The cost of the EOI maneuver is evaluated as(12)ΔVN+1=vEOI+−vEOI−.


### 2.4. Optimization Problem

For an assigned sequence of encountered bodyMseq, the capture problem can be formally defined as follows:
(13)PMseq=minxΔVtotxs.t.xL<x<xuwhere ΔVtot=∑k=0N+1ΔVk is the overall cost of the capture trajectory, accounting for deep-space, release, and EOI maneuvers; x is a vector of design variables given by
(14)x=t0,ΔT0,Δθ∪rπ,j,βj,ΔTj,ηjj=1,..,N,and xL and xU are, respectively, lower and upper bounds of the design variables. The whole trajectory is thus parametrized by using 3+4N parameters.


## 2.1. Problem Overview

According to the ongoing proposal [4], the probe is assumed to be released by the main spacecraft after a few Europa flybys have been completed. In particular, ETP starts its own transfer at the apocenter of a Jovian orbit of period four times the period of Europa (TEu), pericenter equal to the Europa semi-major axis, and coplanar with the Europa orbit; the orientation of the major axis is left free. A target circular quasi-polar orbit around Europa is desired, of assigned altitude hEOI=250 km over Europa’s surface. Four R-6D bipropellant thrusters form the primary propulsion system of the probe, which allows for a total thrust of 88 N with specific impulse Isp=294 s. This propulsion system will be used for deep-space maneuvers, Europa Orbit Insertion (EOI) maneuver, and orbit maintenance during the scientific part of the mission.A probe “net” massmu=146.6 kg, which does not account for the propellant mp and tank ms masses, was estimated in [4]. Assuming a structural coefficient ϵ=ms+mp/ms=6, which is a reasonable value for liquid propellant systems, a maximum value of velocity increment ΔVmax=1240.8 m/s can be obtained if the spacecraft wet mass m0=mu+mp+ms is constrained at 250 kg. This value of ΔVmax must cover the orbit maintenance (about 43.2 m/s for a 6-month mission) and capture cost. The goal is to reduce the ΔV required for the capture, so that a convenient safety margin is left.


## 2.2. Dynamical Model

A patched-conic model is assumed for the present analysis. Flybys are modeled as instantaneous changes in velocity. Subscripts “−” and “+” are used to distinguish between values immediately before or after the discontinuity, respectively. The radius of the sphere of influence (SOI) of the secondary bodies and the travel time inside these regions are assumed to be negligible. Powered flybys are neglected, as considered useless to reduce propellant consumption [19]. An impulsive-thrust model is adopted. This assumption well suits deep space maneuvers (DSMs), which are performed at a large distance from the main body and require usually a quite short time if compared to the orbital period, as chemical engines are here considered. This assumption is also used for the EOI maneuver, even though finite-thrust losses might be considered. Only one DSM is permitted between a flyby and the other. As a further assumption, Jovian moons move on Keplerian orbits (even though the proposed procedure is soon applicable to the general case that uses planetary ephemeris).Despite its simplicity, this model allows to capture the most prominent features of the mission, while keeping the analysis simple enough. In fact, under the hypothesis of impulsive thrust, the trajectory can be computed analytically, without involving the numerical integration of the complete equations of motion.


## 2.3. Trajectory Parameterization

Let us assume that a sequenceM=Mj∣j=1,…,N+1 of N+1 body encounters, where Mj∈1..4 identifies the encountered body (1 = Io, 2 = Eu, 3 = Ga, and 4 = Ca), has been established. Europa is the first and last body in the series. The spacecraft trajectory can be modeled according to the multiple gravity assist-one deep space maneuver (MGA-1DSM) formulation [16]. The trajectory is broken down into a series of body-to-body legs. Each leg starts with a flyby and is made up of two ballistic arcs, joined by an impulsive maneuver.This general formulation for a multigravity assist trajectory can be adapted to the problem at hand by adding an initial leg, which moves the probe from the assigned initial conditions to the first encounter with Europa. The mission ends with a last approach to Europa’s surface, where an impulsive maneuver inserts the probe into the assigned polar orbit.


### 2.3.1. Departure Leg

The departure leg, which connects ETP release position to the first encounter with Europa, is modeled as a Lambert arc, where release epocht0, flight angle Δθ, and flight time ΔT0 are design parameters to optimize.Leta^r,a^t,a^n be a Jovicentric radial-traversal-normal reference frame connected to a Europa position at epoch t1=t0+ΔT0, that is,
(1)a^r=r1r1a^n=r1×vm,1r1×vm,1a^t=a^n×a^r,where r1=rM1t1 and vm,1=vM1t1 indicates position and velocity vectors of Europa (M1=2) at time t1, respectively, which are provided by the ephemeris.The probe departure point is located on a circle of radiusr~0 which lies on the plane a^r−a^t; hence, it can be expressed as
(2)r0=r~a4:1cosΔθa^r−sinΔθa^t,while the velocity immediately before the release maneuver is
(3)v0−=v~a4:1sinΔθa^r+cosΔθa^t,where the values r~a4:1 and v~a4:1 are, respectively, the radius and velocity magnitude at the apocenter of a Jovian orbit with pericenter equal to the Europa semi-major axis, and period 4TEu. In this respect, the problem solution will eventually define the optimal orientation of the line of apses of the initial orbit with respect to Europa’s orbit.The velocity vectors after the release maneuverv0+ and immediately before the first flyby v1− can now be evaluated by solving the associated Lambert problem:
(4)Lambertr0,r1,ΔT0⟶v0+,v1−.The propulsive cost of the release maneuver is evaluated as(5)ΔV0=v0+−v0−.


### 2.3.2. Intermediate Legs

Thekth leg can be parametrized by using four parameters rπ,k, βk, ΔTk, and ηk, which represent, respectively, the flyby radius, the flyby plane orientation, the (overall) leg flight time, and the fraction of the leg flight time at which DSM occurs. The trajectory associated to a generic intermediate leg is presented in Figure 1.Figure 1
Trajectory sketch for an intermediate leg.Beginning from timet1 at the first encounter, the epochs of the following encounters can be evaluated recursively as tk+1=tk+ΔTk, while the DSM epochs are tk+1/2=tk+ηkΔTk. Moon position rk=rMktk and velocity vm,k=vMktk are obtained from the body ephemeris. The position of the spacecraft at the flyby must be the same as the intercepted body; hence, rk−=rk+=rk.Let’s define a body-centered reference systeme^1,e^2,e^3, where e^1 is directed along the spacecraft incoming relative velocity, e^2 is orthogonal to the moon orbital plane, and e^3 closes the right-hand side term
(6)e^1=v∞k−v∞k−e^2=e^1×vm,ke^1×vm,ke^3=e^1×e^2,where v∞k=vk−−vm,k is the probe relative velocity before the flyby. The probe velocity vk+ after the flyby is evaluated as
(7)vk+=vm,k+v∞,kcosδke^1+cosβksinδke^2+sinβksinδke^3,where the rotation δk of the hyperbolic excess velocity in the flyby plane is
(8)δk=2asinμMkμMk+rπ,kv∞k2.Once the spacecraft state after the flyby is fully known, the position and velocity just before DSM, which arerk+1/2 and vk+1/2−, respectively, can be evaluated analytically by using propagation formulas for Keplerian orbits [20].Position of the spacecraft is known at both DSM maneuver and next flyby; the transfer time1−ηkΔTk is also known. Velocity vectors vk+1/2+ immediately after the DSM and vk+1− just before the subsequent flyby can be evaluated by solving the associated Lambert problem
(9)Lambertrk+1/2,rk+1,ΔTk⟶vk+1/2+,vk+1−Evaluation of the DSM propulsive cost is now straightforward(10)ΔVk=vk+1/2+−vk+1/2−


### 2.3.3. Europa Orbit Insertion

The ETP injection into the target polar orbit around Europa is modeled as an impulsive burn applied at the hyperbola pericenter, which is purposely located on the target circular orbit.(11)vEOI−=v∞2+2μEurEOIvEOI+=μEurEOI.The cost of the EOI maneuver is evaluated as(12)ΔVN+1=vEOI+−vEOI−.


## 2.3.1. Departure Leg

The departure leg, which connects ETP release position to the first encounter with Europa, is modeled as a Lambert arc, where release epocht0, flight angle Δθ, and flight time ΔT0 are design parameters to optimize.Leta^r,a^t,a^n be a Jovicentric radial-traversal-normal reference frame connected to a Europa position at epoch t1=t0+ΔT0, that is,
(1)a^r=r1r1a^n=r1×vm,1r1×vm,1a^t=a^n×a^r,where r1=rM1t1 and vm,1=vM1t1 indicates position and velocity vectors of Europa (M1=2) at time t1, respectively, which are provided by the ephemeris.The probe departure point is located on a circle of radiusr~0 which lies on the plane a^r−a^t; hence, it can be expressed as
(2)r0=r~a4:1cosΔθa^r−sinΔθa^t,while the velocity immediately before the release maneuver is
(3)v0−=v~a4:1sinΔθa^r+cosΔθa^t,where the values r~a4:1 and v~a4:1 are, respectively, the radius and velocity magnitude at the apocenter of a Jovian orbit with pericenter equal to the Europa semi-major axis, and period 4TEu. In this respect, the problem solution will eventually define the optimal orientation of the line of apses of the initial orbit with respect to Europa’s orbit.The velocity vectors after the release maneuverv0+ and immediately before the first flyby v1− can now be evaluated by solving the associated Lambert problem:
(4)Lambertr0,r1,ΔT0⟶v0+,v1−.The propulsive cost of the release maneuver is evaluated as(5)ΔV0=v0+−v0−.


## 2.3.2. Intermediate Legs

Thekth leg can be parametrized by using four parameters rπ,k, βk, ΔTk, and ηk, which represent, respectively, the flyby radius, the flyby plane orientation, the (overall) leg flight time, and the fraction of the leg flight time at which DSM occurs. The trajectory associated to a generic intermediate leg is presented in Figure 1.Figure 1
Trajectory sketch for an intermediate leg.Beginning from timet1 at the first encounter, the epochs of the following encounters can be evaluated recursively as tk+1=tk+ΔTk, while the DSM epochs are tk+1/2=tk+ηkΔTk. Moon position rk=rMktk and velocity vm,k=vMktk are obtained from the body ephemeris. The position of the spacecraft at the flyby must be the same as the intercepted body; hence, rk−=rk+=rk.Let’s define a body-centered reference systeme^1,e^2,e^3, where e^1 is directed along the spacecraft incoming relative velocity, e^2 is orthogonal to the moon orbital plane, and e^3 closes the right-hand side term
(6)e^1=v∞k−v∞k−e^2=e^1×vm,ke^1×vm,ke^3=e^1×e^2,where v∞k=vk−−vm,k is the probe relative velocity before the flyby. The probe velocity vk+ after the flyby is evaluated as
(7)vk+=vm,k+v∞,kcosδke^1+cosβksinδke^2+sinβksinδke^3,where the rotation δk of the hyperbolic excess velocity in the flyby plane is
(8)δk=2asinμMkμMk+rπ,kv∞k2.Once the spacecraft state after the flyby is fully known, the position and velocity just before DSM, which arerk+1/2 and vk+1/2−, respectively, can be evaluated analytically by using propagation formulas for Keplerian orbits [20].Position of the spacecraft is known at both DSM maneuver and next flyby; the transfer time1−ηkΔTk is also known. Velocity vectors vk+1/2+ immediately after the DSM and vk+1− just before the subsequent flyby can be evaluated by solving the associated Lambert problem
(9)Lambertrk+1/2,rk+1,ΔTk⟶vk+1/2+,vk+1−Evaluation of the DSM propulsive cost is now straightforward(10)ΔVk=vk+1/2+−vk+1/2−


## 2.3.3. Europa Orbit Insertion

The ETP injection into the target polar orbit around Europa is modeled as an impulsive burn applied at the hyperbola pericenter, which is purposely located on the target circular orbit.(11)vEOI−=v∞2+2μEurEOIvEOI+=μEurEOI.The cost of the EOI maneuver is evaluated as(12)ΔVN+1=vEOI+−vEOI−.


## 2.4. Optimization Problem

For an assigned sequence of encountered bodyMseq, the capture problem can be formally defined as follows:
(13)PMseq=minxΔVtotxs.t.xL<x<xuwhere ΔVtot=∑k=0N+1ΔVk is the overall cost of the capture trajectory, accounting for deep-space, release, and EOI maneuvers; x is a vector of design variables given by
(14)x=t0,ΔT0,Δθ∪rπ,j,βj,ΔTj,ηjj=1,..,N,and xL and xU are, respectively, lower and upper bounds of the design variables. The whole trajectory is thus parametrized by using 3+4N parameters.


## 3. Optimization Algorithm

The global optimization problem stated in Section2 presents several features that make it hard to solve with local optimization approaches (e.g., Newton-like methods) as (i) multiple local optima exist, (ii) the solution may not be defined for some “unfortunate” set of optimization variables, and (iii) the gradient of the objective function is often not available or does not exist at all. In these cases, stochastic algorithms are usually preferred.In the present application, an optimization algorithm based on differential evolution (DE) has been employed. DE is a population-based algorithm, featuring simple and efficient heuristics for global optimization problems defined over a continuous space [21]. Its good performance on several benchmarks and real-world problems drew the attention of many researchers all over the world, who further improved the effectiveness of the algorithm, by devising many variants. A review of the state of art can be found in [22]. The implementation adopted in this study collects several of these ideas.


### 3.1. Standard Differential Evolution

A brief description of the standard DE algorithm is here provided. Let us consider the minimization problem(15)minxfxs.t.xl≤x≤xU.A population ofNP candidate solutions Pop=xi,i=1..NP is randomly created, and for each individual (or agent) xi∈ℝND, the corresponding fitness fxi is evaluated. Then, a new population Popnew is constructed by repeating, for each vector xi belonging to Pop, a sequence of mutation/crossover/selection steps, defined as follows.


#### 3.1.1. Mutation

A mutated vectorvi is created as a linear combination of a few population members. More precisely, several mutation rules were proposed in order to attain either a better exploration of the search space or a faster convergence (exploitation). In the present implementation, the following four strategies, among those available in literature, are adopted:
(16)vi=xr1+Fxr2−xr3,vi=xGbest+Fxr1−xr2,vi=xi+FxGbest−xi+Fxr1−xr2,vi=xGbest+Fxr1−xr2+Fxr3−xr4,where F is a parameter controlling the mutation scale, xGbest is the best individual (or agent) of the current generation, and r1,…,r4 represent randomly chosen, nonrepeated, indexes belonging to 1..NP.Each strategy has weaknesses and strengths: Strategies based on mutation of the best individual (strategies 2 and 4) typically show a faster converge rate toward an (often local) minimum, whereas strategies based on randomly chosen individuals (strategies 1 and 3) better explore the whole search space.


#### 3.1.2. Crossover

A trial vectorui is obtained by a crossover between the target vector xi and the mutated vector vi
(17)ui,j=vi,jifpj≤CRorj=jrandxi,jelsej=1,…,ND,where pj is a random number between 0 and 1, CR is an algorithm parameter (typical values about 0.5), and jrand is a randomly chosen index in the range 1..ND.


#### 3.1.3. Selection

Target and trial vectors are compared. The best one is retained and inserted in the new population(18)xinew=uiiffui<fxixielse.This process is repeated iteratively, creating at each generation a new population which replaces the previous one. The procedure ends after a fixed number of generations (NG). The best attained solution is deemed the optimal problem solution.


### 3.2. Self-Adaptation of Control Parameters

A common practical issue for many stochastic algorithms concerns the selection of suitable values for the control parameters. A fine tuning is often required in order to make the algorithm suitable for complex, real-world problems. DE is privileged in this respect, as the number of its parameters is very low. Apart from the the population sizeNP, the performance of the DE algorithm depends on an appropriate selection of the scale factor F, which controls the mutation phase, and crossover probability CR, which controls the crossover phase.In order to avoid the manual tuning of DE control parameters, a self-adaptation scheme [23] has been implemented. The values of F and CR are encoded into the individuals, which enter the optimization procedure, and randomly initialized within the intervals Fmin,Fmax and CRmin,CRmax, respectively. Better values of these encoded control parameters will presumably lead to better individuals who, in turn, are more prone to survive and produce offspring, thus propagating these “superior” control parameters.In order to maintain a certain diversity of the control parameters among the population, at the end of each generation, each individual undergoes a random uniform mutation of his control parameters which happens with probabilitypτ, that is,
(19)FG+1=UFmin,Fmax,ifp1≤pτFG,otherwise,CRG+1=UCRmin,CRmax,ifp2≤pτCRG,otherwise,where p1 and p2 are random numbers in 0,1, the apexes G and G+1 refer to the current and next generation, respectively, and Ua,b indicates to a randomly sampled number in the range a,b. In the present implementation, the following values are used: Fmin=0.1, Fmax=1, CRmin=0.5, CRmax=1, and pτ=0.1.


### 3.3. Balancing Local and Global Search

A key point in designing a global optimization algorithm concerns the delicate balance between two opposite needs: “global” exploration and “local” exploitation. Exploration refers to the algorithm capability of probing wide portions of the search space, with the hope of finding promising solutions that are yet to refine, while exploitation is the ability of improving a previous solution by probing a small region of the search space around it. A proper balance is required for the success of the algorithm: favouring local search reduces computational time, at the risk of being trapped in local optima, whereas favouring the global search requires longer computation time as a wide portion of the search space has to be tested. In order to achieve a good balance between exploration and exploitation, the proposed algorithm encompasses the creation of different populations, or tribes, each located on an “island.”Each tribe evolves independently from the others and features one specific mutation strategy among the four proposed variants. As shown in Figure2, the different islands are arranged in a radial configuration, so that tribes on inner islands feature the less exploiting (i.e., most exploring) mutation strategy, whereas tribes on the outer islands feature the most exploiting mutation strategy. The order of the strategies from the outer to the inner islands is then 2, 4, 3, and 1. Every 100 generations, a migration is performed: Each tribe passes its best three agents to the “following one” (if it does exist, according to the direction of migration), in which these agents replace the three worst agents. Outward and inward migration tides alternate at each migration event. The proposed scheme allows an easy parallelization. A minimum of 4 islands are required to make it effective, but it scales well on any architecture with 4x cores.Figure 2
Migration scheme: forward (red) and backward (blue), for the 16-island case.
(a)
Outward tide
(b)
Inward tide


### 3.4. Termination Criterion

The termination criterion is mainly based on the generation number, that is, on the available computational budget. This parameter strongly depends on the complexity of the analyzed problem. The maximum number of generation is chosen in such a way that the outcomes of independent runs of the code bring, in almost all cases, to similar results. In this respect, maintaining a certain diversity of the population in terms of its distribution on the search space is mandatory, to avoid a premature convergence on a local optimum.A partial-restart mechanism, hereafter named “Epidemic,” is adopted to handle this issue. More precisely, population diversity is evaluated for each tribe at the end of each generation, by using as metric the Euclidean distance between pairs of solutions. If the diversity score falls under a certain threshold, a large part (90%) of the population is randomly reinitialized over the entire search space. The maximum number of epidemic events that may occur in any run is fixed (2, in this application) in order not to compromise the overall efficiency of the search.The reported solutions have been obtained by exploiting an 8-tribe optimization engine, with 512 agents per tribe and a maximum number of generations equal to 10000. In order to mitigate potential issues due to the stochastic nature of the algorithm, that is, an unfortunate premature convergence on a suboptimal solution, several runs can be performed starting from a different initial population, increasing the confidence on the attained result. In the present application, each optimization is repeated 25 times, and the best found solution is assumed as putative optimal solution.


## 3.1. Standard Differential Evolution

A brief description of the standard DE algorithm is here provided. Let us consider the minimization problem(15)minxfxs.t.xl≤x≤xU.A population ofNP candidate solutions Pop=xi,i=1..NP is randomly created, and for each individual (or agent) xi∈ℝND, the corresponding fitness fxi is evaluated. Then, a new population Popnew is constructed by repeating, for each vector xi belonging to Pop, a sequence of mutation/crossover/selection steps, defined as follows.


### 3.1.1. Mutation

A mutated vectorvi is created as a linear combination of a few population members. More precisely, several mutation rules were proposed in order to attain either a better exploration of the search space or a faster convergence (exploitation). In the present implementation, the following four strategies, among those available in literature, are adopted:
(16)vi=xr1+Fxr2−xr3,vi=xGbest+Fxr1−xr2,vi=xi+FxGbest−xi+Fxr1−xr2,vi=xGbest+Fxr1−xr2+Fxr3−xr4,where F is a parameter controlling the mutation scale, xGbest is the best individual (or agent) of the current generation, and r1,…,r4 represent randomly chosen, nonrepeated, indexes belonging to 1..NP.Each strategy has weaknesses and strengths: Strategies based on mutation of the best individual (strategies 2 and 4) typically show a faster converge rate toward an (often local) minimum, whereas strategies based on randomly chosen individuals (strategies 1 and 3) better explore the whole search space.


### 3.1.2. Crossover

A trial vectorui is obtained by a crossover between the target vector xi and the mutated vector vi
(17)ui,j=vi,jifpj≤CRorj=jrandxi,jelsej=1,…,ND,where pj is a random number between 0 and 1, CR is an algorithm parameter (typical values about 0.5), and jrand is a randomly chosen index in the range 1..ND.


### 3.1.3. Selection

Target and trial vectors are compared. The best one is retained and inserted in the new population(18)xinew=uiiffui<fxixielse.This process is repeated iteratively, creating at each generation a new population which replaces the previous one. The procedure ends after a fixed number of generations (NG). The best attained solution is deemed the optimal problem solution.


## 3.1.1. Mutation

A mutated vectorvi is created as a linear combination of a few population members. More precisely, several mutation rules were proposed in order to attain either a better exploration of the search space or a faster convergence (exploitation). In the present implementation, the following four strategies, among those available in literature, are adopted:
(16)vi=xr1+Fxr2−xr3,vi=xGbest+Fxr1−xr2,vi=xi+FxGbest−xi+Fxr1−xr2,vi=xGbest+Fxr1−xr2+Fxr3−xr4,where F is a parameter controlling the mutation scale, xGbest is the best individual (or agent) of the current generation, and r1,…,r4 represent randomly chosen, nonrepeated, indexes belonging to 1..NP.Each strategy has weaknesses and strengths: Strategies based on mutation of the best individual (strategies 2 and 4) typically show a faster converge rate toward an (often local) minimum, whereas strategies based on randomly chosen individuals (strategies 1 and 3) better explore the whole search space.


## 3.1.2. Crossover

A trial vectorui is obtained by a crossover between the target vector xi and the mutated vector vi
(17)ui,j=vi,jifpj≤CRorj=jrandxi,jelsej=1,…,ND,where pj is a random number between 0 and 1, CR is an algorithm parameter (typical values about 0.5), and jrand is a randomly chosen index in the range 1..ND.


## 3.1.3. Selection

Target and trial vectors are compared. The best one is retained and inserted in the new population(18)xinew=uiiffui<fxixielse.This process is repeated iteratively, creating at each generation a new population which replaces the previous one. The procedure ends after a fixed number of generations (NG). The best attained solution is deemed the optimal problem solution.


## 3.2. Self-Adaptation of Control Parameters

A common practical issue for many stochastic algorithms concerns the selection of suitable values for the control parameters. A fine tuning is often required in order to make the algorithm suitable for complex, real-world problems. DE is privileged in this respect, as the number of its parameters is very low. Apart from the the population sizeNP, the performance of the DE algorithm depends on an appropriate selection of the scale factor F, which controls the mutation phase, and crossover probability CR, which controls the crossover phase.In order to avoid the manual tuning of DE control parameters, a self-adaptation scheme [23] has been implemented. The values of F and CR are encoded into the individuals, which enter the optimization procedure, and randomly initialized within the intervals Fmin,Fmax and CRmin,CRmax, respectively. Better values of these encoded control parameters will presumably lead to better individuals who, in turn, are more prone to survive and produce offspring, thus propagating these “superior” control parameters.In order to maintain a certain diversity of the control parameters among the population, at the end of each generation, each individual undergoes a random uniform mutation of his control parameters which happens with probabilitypτ, that is,
(19)FG+1=UFmin,Fmax,ifp1≤pτFG,otherwise,CRG+1=UCRmin,CRmax,ifp2≤pτCRG,otherwise,where p1 and p2 are random numbers in 0,1, the apexes G and G+1 refer to the current and next generation, respectively, and Ua,b indicates to a randomly sampled number in the range a,b. In the present implementation, the following values are used: Fmin=0.1, Fmax=1, CRmin=0.5, CRmax=1, and pτ=0.1.


## 3.3. Balancing Local and Global Search

A key point in designing a global optimization algorithm concerns the delicate balance between two opposite needs: “global” exploration and “local” exploitation. Exploration refers to the algorithm capability of probing wide portions of the search space, with the hope of finding promising solutions that are yet to refine, while exploitation is the ability of improving a previous solution by probing a small region of the search space around it. A proper balance is required for the success of the algorithm: favouring local search reduces computational time, at the risk of being trapped in local optima, whereas favouring the global search requires longer computation time as a wide portion of the search space has to be tested. In order to achieve a good balance between exploration and exploitation, the proposed algorithm encompasses the creation of different populations, or tribes, each located on an “island.”Each tribe evolves independently from the others and features one specific mutation strategy among the four proposed variants. As shown in Figure2, the different islands are arranged in a radial configuration, so that tribes on inner islands feature the less exploiting (i.e., most exploring) mutation strategy, whereas tribes on the outer islands feature the most exploiting mutation strategy. The order of the strategies from the outer to the inner islands is then 2, 4, 3, and 1. Every 100 generations, a migration is performed: Each tribe passes its best three agents to the “following one” (if it does exist, according to the direction of migration), in which these agents replace the three worst agents. Outward and inward migration tides alternate at each migration event. The proposed scheme allows an easy parallelization. A minimum of 4 islands are required to make it effective, but it scales well on any architecture with 4x cores.Figure 2
Migration scheme: forward (red) and backward (blue), for the 16-island case.
(a)
Outward tide
(b)
Inward tide


## 3.4. Termination Criterion

The termination criterion is mainly based on the generation number, that is, on the available computational budget. This parameter strongly depends on the complexity of the analyzed problem. The maximum number of generation is chosen in such a way that the outcomes of independent runs of the code bring, in almost all cases, to similar results. In this respect, maintaining a certain diversity of the population in terms of its distribution on the search space is mandatory, to avoid a premature convergence on a local optimum.A partial-restart mechanism, hereafter named “Epidemic,” is adopted to handle this issue. More precisely, population diversity is evaluated for each tribe at the end of each generation, by using as metric the Euclidean distance between pairs of solutions. If the diversity score falls under a certain threshold, a large part (90%) of the population is randomly reinitialized over the entire search space. The maximum number of epidemic events that may occur in any run is fixed (2, in this application) in order not to compromise the overall efficiency of the search.The reported solutions have been obtained by exploiting an 8-tribe optimization engine, with 512 agents per tribe and a maximum number of generations equal to 10000. In order to mitigate potential issues due to the stochastic nature of the algorithm, that is, an unfortunate premature convergence on a suboptimal solution, several runs can be performed starting from a different initial population, increasing the confidence on the attained result. In the present application, each optimization is repeated 25 times, and the best found solution is assumed as putative optimal solution.


## 4. Preliminary Analysis

A preliminary analysis is carried out in order to define the main features of the tentative trajectory that undergoes the DE optimization process. For the sake of simplicity, the orbits of the Galilean moons are assumed to be circular and coplanar. The Tisserand graph and the numerical suboptimal solution of thev∞ leveraging or ΔV-Europa gravity assist (ΔV-EGA) maneuver are effective tools in defining the preliminary solution. These tools are here described with reference to the dual problem of escaping from Europa, because it is more intuitive and several papers are found in literature, dealing with v∞ leveraging to move away from the Earth.


### 4.1. Tisserand Graph

The Tisserand graph is a powerful instrument for the preliminary design of trajectories of a spacecraft that exploits multiple gravity assists in a multibody planetary system. Several versions have been discussed by different authors. The two-dimensional apoapsis-periapsis (orra−rp) plot is adopted here.Any point on the graph (ra>rp>0) identifies a Keplerian elliptical orbit of the spacecraft around the central body (Jupiter). For a given secondary body M, moving on a circular orbit of radius rM, the region AM=ra,rp∣rp<rM<ra encompasses all spacecraft orbits that intersect the moon orbit. For each point ra,rp in this region, one can easily compute the values of magnitude and direction of the hyperbolic excess velocity v∞ in case of encounter. The direction of v∞ with respect to the moon velocity vm is defined by the pump angle α (see Figure 3).Figure 3
Velocity before (−) and after (+) the flyby.Curves of constant hyperbolic excess velocity are typically superimposed on this plot. Each curve collects all Jovicentric orbits approaching the moon with the samev∞ magnitude and pump angle ranging from 0 degrees (rp=rM, rightmost point) to 180 degrees (ra=rM, leftmost point).Lines of constant orbital period (i.e., constant energy) can also be plotted. This feature is mainly used to representK : L resonant orbits, where K and L indicate the number of revolutions completed by moon and spacecraft, respectively, between two successive encounters.Starting from any point on the graph, a shift to another orbit along a curve of constantv∞ can be obtained at no cost, as result of a gravity assist. The displacement along the curve depends on the rotation of the hyperbolic excess velocity δ, which is given by 8.Sequences of flybys can be performed, provided that the spacecraft is able to re-encounter the moon. The role of resonant orbits is now apparent, as hopping between resonant orbits guarantees the re-encounter in a known time.For each resonant orbit, one can draw a pair of maximum deviation curves, assuming that a minimum-altitude flyby is enforced. These curves encompass all Jovicentric orbits achievable from that resonant orbit, by means of a gravity assist. Maximum deviation curves show, for any initialv∞ value, which resonances are available for the next jump.Points that do not lay on the samev∞ curve can be connected only by means of some propulsive maneuver. As an example, vertical displacements on the graph can be obtained as a result of a tangential impulse at the apocenter.Figure4 shows a Tisserand plot for the Jovian planetary system. Contour lines of constant v∞ value are displayed for Europa and Ganymede. Resonant orbits corresponding to resonance 6 : 5, 3 : 2, 7 : 4, 2 : 1, and 5 : 2 with Europa and associated maximum-deviation curves are shown as an example.Figure 4
Tisserand graph for the Jovian planetary system. Contour lines of constantv∞, at steps of 2% of Vm. Sample resonant orbits with Europa and corresponding maximum-deviation curves for a minimum flyby altitude hπ,min=30km are also depicted.


### 4.2. Suboptimal Solution forV Leveraging Transfers

The maneuver known asv∞ leveraging, illustrated in Figure 5, provides an efficient way to move between points on the Tisserand graph that do not lay on the same v∞ curve, as a large change of v∞ is obtained by means of a small ΔV. The spacecraft leaves the secondary body (Europa) with a given hyperbolic excess velocity v∞, entering a near-resonant orbit around the central body (Jupiter). A deep-space maneuver is performed near the apoapsis in order to modify v∞ (increase, if one is considering an escape trajectory, or decrease, if one instead considers a capture) at the next encounter of the probe with the secondary body. Here, a flyby permits to rotate the hyperbolic excess of velocity, allowing the probe to enter a Jovicentric orbit with a larger apocenter. A generic ΔV-EGA maneuver can be labeled as K : LM±, where K and L indicate the number of revolutions that moon and probe complete, respectively, before the next encounter; M≤L denotes the spacecraft orbit where DSM occours; the sign ± is used to distinguish between the possible encounter locations: just before (−) or after (+) traversing the line of apsides.Figure 5
ΔV-EGA trajectory.A suboptimal solution forΔV-EGA maneuvers is easily achieved [14] under the assumption that (i) the spacecraft initial orbit is tangent to the moon’s orbit and (ii) a tangential ΔV is applied at the apocenter of the nearly-resonant K : L orbit. The maneuver is slightly improved if the aforementioned assumptions are removed [24]. Dealing with high-v∞ flybys, the replacement of the patched-conic model with the RTBP model does not modify the mission and the numerical performance [25].


## 4.1. Tisserand Graph

The Tisserand graph is a powerful instrument for the preliminary design of trajectories of a spacecraft that exploits multiple gravity assists in a multibody planetary system. Several versions have been discussed by different authors. The two-dimensional apoapsis-periapsis (orra−rp) plot is adopted here.Any point on the graph (ra>rp>0) identifies a Keplerian elliptical orbit of the spacecraft around the central body (Jupiter). For a given secondary body M, moving on a circular orbit of radius rM, the region AM=ra,rp∣rp<rM<ra encompasses all spacecraft orbits that intersect the moon orbit. For each point ra,rp in this region, one can easily compute the values of magnitude and direction of the hyperbolic excess velocity v∞ in case of encounter. The direction of v∞ with respect to the moon velocity vm is defined by the pump angle α (see Figure 3).Figure 3
Velocity before (−) and after (+) the flyby.Curves of constant hyperbolic excess velocity are typically superimposed on this plot. Each curve collects all Jovicentric orbits approaching the moon with the samev∞ magnitude and pump angle ranging from 0 degrees (rp=rM, rightmost point) to 180 degrees (ra=rM, leftmost point).Lines of constant orbital period (i.e., constant energy) can also be plotted. This feature is mainly used to representK : L resonant orbits, where K and L indicate the number of revolutions completed by moon and spacecraft, respectively, between two successive encounters.Starting from any point on the graph, a shift to another orbit along a curve of constantv∞ can be obtained at no cost, as result of a gravity assist. The displacement along the curve depends on the rotation of the hyperbolic excess velocity δ, which is given by 8.Sequences of flybys can be performed, provided that the spacecraft is able to re-encounter the moon. The role of resonant orbits is now apparent, as hopping between resonant orbits guarantees the re-encounter in a known time.For each resonant orbit, one can draw a pair of maximum deviation curves, assuming that a minimum-altitude flyby is enforced. These curves encompass all Jovicentric orbits achievable from that resonant orbit, by means of a gravity assist. Maximum deviation curves show, for any initialv∞ value, which resonances are available for the next jump.Points that do not lay on the samev∞ curve can be connected only by means of some propulsive maneuver. As an example, vertical displacements on the graph can be obtained as a result of a tangential impulse at the apocenter.Figure4 shows a Tisserand plot for the Jovian planetary system. Contour lines of constant v∞ value are displayed for Europa and Ganymede. Resonant orbits corresponding to resonance 6 : 5, 3 : 2, 7 : 4, 2 : 1, and 5 : 2 with Europa and associated maximum-deviation curves are shown as an example.Figure 4
Tisserand graph for the Jovian planetary system. Contour lines of constantv∞, at steps of 2% of Vm. Sample resonant orbits with Europa and corresponding maximum-deviation curves for a minimum flyby altitude hπ,min=30km are also depicted.


## 4.2. Suboptimal Solution forV Leveraging Transfers

The maneuver known asv∞ leveraging, illustrated in Figure 5, provides an efficient way to move between points on the Tisserand graph that do not lay on the same v∞ curve, as a large change of v∞ is obtained by means of a small ΔV. The spacecraft leaves the secondary body (Europa) with a given hyperbolic excess velocity v∞, entering a near-resonant orbit around the central body (Jupiter). A deep-space maneuver is performed near the apoapsis in order to modify v∞ (increase, if one is considering an escape trajectory, or decrease, if one instead considers a capture) at the next encounter of the probe with the secondary body. Here, a flyby permits to rotate the hyperbolic excess of velocity, allowing the probe to enter a Jovicentric orbit with a larger apocenter. A generic ΔV-EGA maneuver can be labeled as K : LM±, where K and L indicate the number of revolutions that moon and probe complete, respectively, before the next encounter; M≤L denotes the spacecraft orbit where DSM occours; the sign ± is used to distinguish between the possible encounter locations: just before (−) or after (+) traversing the line of apsides.Figure 5
ΔV-EGA trajectory.A suboptimal solution forΔV-EGA maneuvers is easily achieved [14] under the assumption that (i) the spacecraft initial orbit is tangent to the moon’s orbit and (ii) a tangential ΔV is applied at the apocenter of the nearly-resonant K : L orbit. The maneuver is slightly improved if the aforementioned assumptions are removed [24]. Dealing with high-v∞ flybys, the replacement of the patched-conic model with the RTBP model does not modify the mission and the numerical performance [25].


## 5. Numerical Results

### 5.1. Tentative Solution

The dual problem of escaping from Europa is considered in this section. A sequence ofΔV-EGA maneuvers with increasing K : L resonances is used to progressively augment the probe apocenter, by targeting, after each flyby, a larger near-resonant orbit. Eventually, the spacecraft enters a hyperbolic trajectory and escapes Jupiter SOI [15]. In the present case, the probe trajectory ends at the apocenter of a 4 : 1 resonant orbit.The tentative solution depends on the performance index that is minimized. In a minimum-propellant problem, the totalΔV is minimized by aligning vector v∞, after each flyby, to Europa orbital velocity. This corresponds to maximizing the energy after the flyby (max-E strategy), assuming assigned v∞. Low increments of v∞ correspond to very small deep-space impulses. The optimization problem is ill-defined without any kind of time constraints, as the optimal solution would require an infinite number of ΔV-EGA maneuvers, each with infinitesimal v∞ increment. Overall flight time is approximately equal to TEu∑K. As a consequence, high-K resonances should be avoided.AssumingKmax=16, the suboptimal minimum-propellant solution presents ΔV = 1.04318 km/s with flight time of 627.8 days. However, by skipping some resonances, while keeping max-E strategy, flight time can be considerably reduced with a minimum increment of propellant expenditure. Europa, due to its small mass, can provide only small rotations of the hyperbolic excess velocity, and long resonances, such as 9 : 5 and 11 : 5, must be included in order to exploit max-E strategy. Flight time remains incompatible with ETP mission requirements.High-K resonances can be skipped by privileging the energy increment provided by each gravity assist, instead of the energy after the flyby. The maximum increment of energy (max-ΔE strategy) implies the maximum change of the semi-major axis, that is, period. In an unconstrained gravity assist, max-ΔE strategy requires the alignment of the hyperbola axis to the velocity of the secondary body (α+=π−α−). Moreover, a v∞ magnitude exists that maximizes ΔE. The present problem is, at a large extent, different. The resonant orbit after the previous gravity assist is assigned, and the pump angle α− before the flyby is not free but a function of v∞. Moreover, the new resonance is “a priori” selected in order to contain the trip time (low K is preferred). Therefore, ΔE is assigned and the mission designer just selects a hyperbolic excess velocity in the range that permits the desired resonance jump. The best v∞ value depends on the entire mission and will be a result of the optimization process. When the max-ΔE strategy (here improperly so called) is adopted, v∞ and pump angle after the gravity assist are large, and a further flyby, without any DSM (and leveraging), could align the spacecraft velocity to the moon orbital velocity.The mission can be further improved by taking advantage from Europa’s eccentricity and by removing constraints on position and direction of the deep-space impulse. This task is carried out by means of the proposed optimization algorithm, based on the general MGA-1DSM formulation described in Section2. A reasonable tentative sequence of resonances is used to prune the global search. A careful exam of the Tisserand graph, in the light of the previous concepts, permits the exclusion of resonances with K>7. Lower and upper bounds of the optimization variables are conveniently adjusted to improve convergence. In particular, for a generic K : L(M) ΔV-EGA maneuver, one has
(20)rπ∈rπ,min,rπ,max,β∈−π,π,ΔT∈K−0.1TEu,K+0.1TEu,η∈M−1L+1e−5,ML.Flyby parametersrπ and β do not depend on the prescribed resonance.A capture trajectory is first attained under the hypothesis of a circular orbit for Europa (solution “A”). Next, this assumption is removed, and a slightly different trajectory is attained (solution “B”). Eventually, a third trajectory, which also exploits a flyby of Ganymede, is proposed (solution “C”).


### 5.2. Solution “A”—Europa Circular Orbit

Solutions according to max-E strategy are impractical as the overall flight time is incompatible with ETP mission requirements. The Tisserand graph actually suggests that a reduction of the flight time can be achieved performing maneuvers based on max-ΔE strategy. Solution “A” is obtained by assuming that Europa moves on a circular orbit, in order to match hypothesis and indications of the Tisserand graph. A suitable sequence of resonant orbits is assumed, and an optimal solution is provided by the DE algorithm. Features of the mission are summarized in Table 1. Each row presents initial time, v∞ magnitude at flyby, resonance, and deep-space impulse of each leg. The last row refers to the impulsive injection into the target orbit around Europa.Table 1
Solution “A”—mission features.
Event
Time [days]
v∞[m/s]
Resonance
DSM [m/s]
0
Departure
0.00
4 : 1
1.931
A
Eu flyby
7.07
3683.4
7 : 2
13.09
B
Eu flyby
32.13
3519.7
3 : 1
15.16
C
Eu flyby
42.74
3351.2
5 : 2
21.00
D
Eu flyby
60.46
3150.8
2 : 1
97.19
E
Eu flyby
67.29
2353.0
7 : 4
40.55
F
Eu flyby
91.98
2016.7
3 : 2
62.49
G
Eu flyby
102.39
1560.1
4 : 3
56.11
H
Eu flyby
116.33
1176.8
6:5
60.56
I
EOI
137.32
803.2
715.18
Totals
137.32
1083.27CaptureΔV is reduced substantially with respect to a direct insertion (2781.9 m/s), allowing to reduce ETP total mass (including propellant and structures) to 239.58 kg. The Jovicentric trajectory is plotted in Figure 6. Deep-space maneuvers are performed substantially at the apogee of the resonant orbits. Flyby locations on Europa orbit move clockwise, indicating that all maneuvers belong to class (−). Only the first flyby does not follow this rule, probably due to the peculiar features of the transition from the initial leg to the first intermediate leg.Figure 6
Jovicentric trajectory of solution A.The Tisserand graph is presented in Figure7. For resonances lower than 2 : 1, after the flyby the spacecraft enters a Jovicentric orbit with perigee equal to Europa radius. This matches the suboptimal solution proposed by Longuski, where the hyperbolic excess of velocity is aligned with the moon velocity. The alignment condition is not verified for resonance higher than 2 : 1, when max-ΔE maneuvers are carried out in order to minimize the number of ΔV-EGA maneuvers, skipping the most time-consuming legs.Figure 7
Tisserand graph of solutionA. Axes are normalized with respect to Europa semi-major axis.


### 5.3. Solution “B”—Europa Elliptic Orbit

The hypothesis on circularity of Europa’s orbit is removed, and a new solution for the same mission scheme is searched for. Table2 describes the attained solution. Figures 8 and 9 present the spacecraft trajectory and the corresponding Tisserand graph, respectively. Minor, yet interesting differences, can be observed. The transition between resonances 2 : 1 and 7 : 4 is now obtained by means of a max-ΔE maneuver. The final flybys and EOI maneuver occur in close proximity of Europa pericenter, in order to benefit from a lower v∞. The propellant consumption is coherent with solution A, and the mission takes slight advantage from the eccentricity of Europa orbit.Table 2
Solution “B”—mission features.
Event
Time [days]
v∞[m/s]
Resonance
DSM [m/s]
0
Departure
1.40
4 : 1
0.01
A
Eu flyby
8.56
3702.4
7 : 2
14.65
B
Eu flyby
33.36
3518.4
3 : 1
16.32
C
Eu flyby
43.96
3335.6
5 : 2
21.58
D
Eu flyby
61.67
3127.5
2:1
79.93
E
Eu flyby
68.60
2478.3
7:4
57.09
F
Eu flyby
93.26
2012.6
3:2
63.14
G
Eu flyby
103.67
1548.9
4 : 3
55.76
H
Eu flyby
117.62
1166.5
6:5
59.85
I
EOI
137.21
797.0
712.76
Totals
137.88
1081.10Figure 8
Jovicentric trajectory of solution B.Figure 9
Tisserand graph of solution B. Axes are normalized with respect to the Europa semi-major axis.


### 5.4. Solution “C”—Ganymede Flyby

Previous solutions show that ETP trajectory crosses Callisto and Ganymede orbits during its journey. An improved solution is searched for, aiming at replacing some deep-space maneuvers with Ganymede or Callisto flybys.In particular, the semi-major axis of Ganymede (1.595REu) almost matches the apoapsis radius of the nominal 3 : 2 orbit (1.621REu). The 3 : 2(2) arc from Europa to Europa is thus substituted by two legs: the first from Europa to Ganymede and the second from Ganymede to Europa. Also, a 5 : 3 resonant orbit is used instead of 7 : 4, reducing the flight time by 2TEu.Results for this mission are summarized in Table3. The Jovicentric trajectory and the Tisserand graph are presented in Figures 10 and 11, respectively. The flybys progressively rotate the initial hyperbolic excess velocity of the spacecraft, moving the spacecraft through a series of resonant orbits, until the apocenter is close to the Ganymede orbit, where a gravity assist from this moon makes the spacecraft orbit almost tangent to Europa orbit. No significant DSM is performed during the legs preceding the Ganymede flyby; thus, the Europa position at the encounters A–F is the same.Table 3
Solution “C”—mission features.
Event
Time [days]
v∞[m/s]
Resonance
DSM [m/s]
0
Departure
3.03
4 : 1
0.64
A
Eu flyby
10.05
3769.2
7 : 2
0.47
B
Eu flyby
34.92
3772.7
3 : 1
0.01
C
Eu flyby
45.58
3772.7
5 : 2
0.63
D
Eu flyby
63.34
3767.9
2 : 1
0.05
E
Eu flyby
70.45
3767.8
5:3
0.06
F
Eu flyby
88.22
3767.7
0.11
G
Ga flyby
95.65
2666.3
0.01
H
Eu flyby
99.68
1784.5
4 : 3
90.24
I
Eu flyby
113.61
1213.2
6:5
67.19
J
EOI
134.59
806.9
716.64
Totals
131.56
882.04Figure 10
Jovicentric trajectory of solution C.Figure 11
Tisserand graph of solution C. Axes are normalized with respect to the Europa semi-major axis.The Tisserand graph for this solution is presented in Figure11, confirming that one Ganymede flyby is sufficient to change v∞ in a so large amount that no deep space maneuver is necessary before two final Europa gravity assists permit the achievement of the 6 : 5 resonant orbit which ends with the EOI maneuver. A closer inspection shows small changes of v∞ in the first part of the mission, despite no impulsive maneuver is performed. This is explained by recalling that v∞ curves on the Tisserand plot are drawn assuming circular and coplanar orbits of the moons. In particular, Ganymede and Europa orbits are not coplanar (their inclinations are 0.2 and 0.4 degree, resp.); therefore, unlike previous solutions, this ETP trajectory is not coplanar with the Europa orbit. Spacecraft inclination and velocity vary at each flyby, resulting into a deceptive v∞ change on the planar circular Tisserand graph.The DE algorithm was able to find a solution with negligible deep-space impulses between flybys A and H. This kind of optimization method suffers the very low variation of the performance index in the proximity of the optimal solution and cannot be very accurate. The quality of the solution can be improved by replacing the general MGA-1DSM formulation with another one that considers ballistic legs between flybys, when required by theoretical considerations.By reducing the captureΔV to 882.04 m/s, ETP mass becomes 218.61 kg, saving approximately 21 kg with respect to previous solutions. The exploitation of a gravity assist by Ganymede also permits a 6-day reduction of the mission time length. A major drawback is the requirement of a suitable initial phasing between Europa and Ganymede, or a precise timing for ETP release, which repeats every 7.058 days, that is, the Ganymede-Europa synodic period. Even though this issue might be thought irrelevant, as nowadays missions withstand many complex trajectory constraints, this is quite a delicate point for a probe that aims at minimal impact on the mission of the mother spacecraft. In this respect, a similar solution should be searched for, only after the baseline mission has been defined. The ETP release condition will be accordingly constrained.


## 5.1. Tentative Solution

The dual problem of escaping from Europa is considered in this section. A sequence ofΔV-EGA maneuvers with increasing K : L resonances is used to progressively augment the probe apocenter, by targeting, after each flyby, a larger near-resonant orbit. Eventually, the spacecraft enters a hyperbolic trajectory and escapes Jupiter SOI [15]. In the present case, the probe trajectory ends at the apocenter of a 4 : 1 resonant orbit.The tentative solution depends on the performance index that is minimized. In a minimum-propellant problem, the totalΔV is minimized by aligning vector v∞, after each flyby, to Europa orbital velocity. This corresponds to maximizing the energy after the flyby (max-E strategy), assuming assigned v∞. Low increments of v∞ correspond to very small deep-space impulses. The optimization problem is ill-defined without any kind of time constraints, as the optimal solution would require an infinite number of ΔV-EGA maneuvers, each with infinitesimal v∞ increment. Overall flight time is approximately equal to TEu∑K. As a consequence, high-K resonances should be avoided.AssumingKmax=16, the suboptimal minimum-propellant solution presents ΔV = 1.04318 km/s with flight time of 627.8 days. However, by skipping some resonances, while keeping max-E strategy, flight time can be considerably reduced with a minimum increment of propellant expenditure. Europa, due to its small mass, can provide only small rotations of the hyperbolic excess velocity, and long resonances, such as 9 : 5 and 11 : 5, must be included in order to exploit max-E strategy. Flight time remains incompatible with ETP mission requirements.High-K resonances can be skipped by privileging the energy increment provided by each gravity assist, instead of the energy after the flyby. The maximum increment of energy (max-ΔE strategy) implies the maximum change of the semi-major axis, that is, period. In an unconstrained gravity assist, max-ΔE strategy requires the alignment of the hyperbola axis to the velocity of the secondary body (α+=π−α−). Moreover, a v∞ magnitude exists that maximizes ΔE. The present problem is, at a large extent, different. The resonant orbit after the previous gravity assist is assigned, and the pump angle α− before the flyby is not free but a function of v∞. Moreover, the new resonance is “a priori” selected in order to contain the trip time (low K is preferred). Therefore, ΔE is assigned and the mission designer just selects a hyperbolic excess velocity in the range that permits the desired resonance jump. The best v∞ value depends on the entire mission and will be a result of the optimization process. When the max-ΔE strategy (here improperly so called) is adopted, v∞ and pump angle after the gravity assist are large, and a further flyby, without any DSM (and leveraging), could align the spacecraft velocity to the moon orbital velocity.The mission can be further improved by taking advantage from Europa’s eccentricity and by removing constraints on position and direction of the deep-space impulse. This task is carried out by means of the proposed optimization algorithm, based on the general MGA-1DSM formulation described in Section2. A reasonable tentative sequence of resonances is used to prune the global search. A careful exam of the Tisserand graph, in the light of the previous concepts, permits the exclusion of resonances with K>7. Lower and upper bounds of the optimization variables are conveniently adjusted to improve convergence. In particular, for a generic K : L(M) ΔV-EGA maneuver, one has
(20)rπ∈rπ,min,rπ,max,β∈−π,π,ΔT∈K−0.1TEu,K+0.1TEu,η∈M−1L+1e−5,ML.Flyby parametersrπ and β do not depend on the prescribed resonance.A capture trajectory is first attained under the hypothesis of a circular orbit for Europa (solution “A”). Next, this assumption is removed, and a slightly different trajectory is attained (solution “B”). Eventually, a third trajectory, which also exploits a flyby of Ganymede, is proposed (solution “C”).


## 5.2. Solution “A”—Europa Circular Orbit

Solutions according to max-E strategy are impractical as the overall flight time is incompatible with ETP mission requirements. The Tisserand graph actually suggests that a reduction of the flight time can be achieved performing maneuvers based on max-ΔE strategy. Solution “A” is obtained by assuming that Europa moves on a circular orbit, in order to match hypothesis and indications of the Tisserand graph. A suitable sequence of resonant orbits is assumed, and an optimal solution is provided by the DE algorithm. Features of the mission are summarized in Table 1. Each row presents initial time, v∞ magnitude at flyby, resonance, and deep-space impulse of each leg. The last row refers to the impulsive injection into the target orbit around Europa.Table 1
Solution “A”—mission features.
Event
Time [days]
v∞[m/s]
Resonance
DSM [m/s]
0
Departure
0.00
4 : 1
1.931
A
Eu flyby
7.07
3683.4
7 : 2
13.09
B
Eu flyby
32.13
3519.7
3 : 1
15.16
C
Eu flyby
42.74
3351.2
5 : 2
21.00
D
Eu flyby
60.46
3150.8
2 : 1
97.19
E
Eu flyby
67.29
2353.0
7 : 4
40.55
F
Eu flyby
91.98
2016.7
3 : 2
62.49
G
Eu flyby
102.39
1560.1
4 : 3
56.11
H
Eu flyby
116.33
1176.8
6:5
60.56
I
EOI
137.32
803.2
715.18
Totals
137.32
1083.27CaptureΔV is reduced substantially with respect to a direct insertion (2781.9 m/s), allowing to reduce ETP total mass (including propellant and structures) to 239.58 kg. The Jovicentric trajectory is plotted in Figure 6. Deep-space maneuvers are performed substantially at the apogee of the resonant orbits. Flyby locations on Europa orbit move clockwise, indicating that all maneuvers belong to class (−). Only the first flyby does not follow this rule, probably due to the peculiar features of the transition from the initial leg to the first intermediate leg.Figure 6
Jovicentric trajectory of solution A.The Tisserand graph is presented in Figure7. For resonances lower than 2 : 1, after the flyby the spacecraft enters a Jovicentric orbit with perigee equal to Europa radius. This matches the suboptimal solution proposed by Longuski, where the hyperbolic excess of velocity is aligned with the moon velocity. The alignment condition is not verified for resonance higher than 2 : 1, when max-ΔE maneuvers are carried out in order to minimize the number of ΔV-EGA maneuvers, skipping the most time-consuming legs.Figure 7
Tisserand graph of solutionA. Axes are normalized with respect to Europa semi-major axis.


## 5.3. Solution “B”—Europa Elliptic Orbit

The hypothesis on circularity of Europa’s orbit is removed, and a new solution for the same mission scheme is searched for. Table2 describes the attained solution. Figures 8 and 9 present the spacecraft trajectory and the corresponding Tisserand graph, respectively. Minor, yet interesting differences, can be observed. The transition between resonances 2 : 1 and 7 : 4 is now obtained by means of a max-ΔE maneuver. The final flybys and EOI maneuver occur in close proximity of Europa pericenter, in order to benefit from a lower v∞. The propellant consumption is coherent with solution A, and the mission takes slight advantage from the eccentricity of Europa orbit.Table 2
Solution “B”—mission features.
Event
Time [days]
v∞[m/s]
Resonance
DSM [m/s]
0
Departure
1.40
4 : 1
0.01
A
Eu flyby
8.56
3702.4
7 : 2
14.65
B
Eu flyby
33.36
3518.4
3 : 1
16.32
C
Eu flyby
43.96
3335.6
5 : 2
21.58
D
Eu flyby
61.67
3127.5
2:1
79.93
E
Eu flyby
68.60
2478.3
7:4
57.09
F
Eu flyby
93.26
2012.6
3:2
63.14
G
Eu flyby
103.67
1548.9
4 : 3
55.76
H
Eu flyby
117.62
1166.5
6:5
59.85
I
EOI
137.21
797.0
712.76
Totals
137.88
1081.10Figure 8
Jovicentric trajectory of solution B.Figure 9
Tisserand graph of solution B. Axes are normalized with respect to the Europa semi-major axis.


## 5.4. Solution “C”—Ganymede Flyby

Previous solutions show that ETP trajectory crosses Callisto and Ganymede orbits during its journey. An improved solution is searched for, aiming at replacing some deep-space maneuvers with Ganymede or Callisto flybys.In particular, the semi-major axis of Ganymede (1.595REu) almost matches the apoapsis radius of the nominal 3 : 2 orbit (1.621REu). The 3 : 2(2) arc from Europa to Europa is thus substituted by two legs: the first from Europa to Ganymede and the second from Ganymede to Europa. Also, a 5 : 3 resonant orbit is used instead of 7 : 4, reducing the flight time by 2TEu.Results for this mission are summarized in Table3. The Jovicentric trajectory and the Tisserand graph are presented in Figures 10 and 11, respectively. The flybys progressively rotate the initial hyperbolic excess velocity of the spacecraft, moving the spacecraft through a series of resonant orbits, until the apocenter is close to the Ganymede orbit, where a gravity assist from this moon makes the spacecraft orbit almost tangent to Europa orbit. No significant DSM is performed during the legs preceding the Ganymede flyby; thus, the Europa position at the encounters A–F is the same.Table 3
Solution “C”—mission features.
Event
Time [days]
v∞[m/s]
Resonance
DSM [m/s]
0
Departure
3.03
4 : 1
0.64
A
Eu flyby
10.05
3769.2
7 : 2
0.47
B
Eu flyby
34.92
3772.7
3 : 1
0.01
C
Eu flyby
45.58
3772.7
5 : 2
0.63
D
Eu flyby
63.34
3767.9
2 : 1
0.05
E
Eu flyby
70.45
3767.8
5:3
0.06
F
Eu flyby
88.22
3767.7
0.11
G
Ga flyby
95.65
2666.3
0.01
H
Eu flyby
99.68
1784.5
4 : 3
90.24
I
Eu flyby
113.61
1213.2
6:5
67.19
J
EOI
134.59
806.9
716.64
Totals
131.56
882.04Figure 10
Jovicentric trajectory of solution C.Figure 11
Tisserand graph of solution C. Axes are normalized with respect to the Europa semi-major axis.The Tisserand graph for this solution is presented in Figure11, confirming that one Ganymede flyby is sufficient to change v∞ in a so large amount that no deep space maneuver is necessary before two final Europa gravity assists permit the achievement of the 6 : 5 resonant orbit which ends with the EOI maneuver. A closer inspection shows small changes of v∞ in the first part of the mission, despite no impulsive maneuver is performed. This is explained by recalling that v∞ curves on the Tisserand plot are drawn assuming circular and coplanar orbits of the moons. In particular, Ganymede and Europa orbits are not coplanar (their inclinations are 0.2 and 0.4 degree, resp.); therefore, unlike previous solutions, this ETP trajectory is not coplanar with the Europa orbit. Spacecraft inclination and velocity vary at each flyby, resulting into a deceptive v∞ change on the planar circular Tisserand graph.The DE algorithm was able to find a solution with negligible deep-space impulses between flybys A and H. This kind of optimization method suffers the very low variation of the performance index in the proximity of the optimal solution and cannot be very accurate. The quality of the solution can be improved by replacing the general MGA-1DSM formulation with another one that considers ballistic legs between flybys, when required by theoretical considerations.By reducing the captureΔV to 882.04 m/s, ETP mass becomes 218.61 kg, saving approximately 21 kg with respect to previous solutions. The exploitation of a gravity assist by Ganymede also permits a 6-day reduction of the mission time length. A major drawback is the requirement of a suitable initial phasing between Europa and Ganymede, or a precise timing for ETP release, which repeats every 7.058 days, that is, the Ganymede-Europa synodic period. Even though this issue might be thought irrelevant, as nowadays missions withstand many complex trajectory constraints, this is quite a delicate point for a probe that aims at minimal impact on the mission of the mother spacecraft. In this respect, a similar solution should be searched for, only after the baseline mission has been defined. The ETP release condition will be accordingly constrained.


## 6. Conclusions

In the search for a low-ΔV capture of a light-weight spacecraft into a polar orbit around Europa, v∞ leveraging represents an excellent technique, which reduces the propellant consumption, at the cost of increased flight time. Moreover, several deep-space maneuvers along Europa-resonant orbits can be replaced by a single Ganymede gravity assist; total ΔV is further reduced with a small saving of the flight time. The Tisserand graph has proved to be an effective tool for the preliminary analysis. The same graph is also useful to introduce in the mission the gravity assist from other Galilean moons.The search for the optimal solution is made complex by the need that the actual mission must take eccentricity and inclination of the secondary bodies into account. The exploitation of an effective heuristic optimization code allowed a wide exploration of the solution space, and an optimal solution is found in a relatively easy and fast manner. The performance of this Europa capture trajectory confirms the feasibility of a low-cost exploration mission, carried out by a small orbiter, compatible with the specifics for the secondary payload of Europa Clipper mission.



---

*Source: 6890173-2018-08-09.xml*