Near-optimal capture in the planar circular restricted Pluto-Charon system

TL;DR

This paper investigates near-optimal capture points around Charon in the Pluto-Charon system using numerical methods within the PCRTBP framework, aiming to aid future space mission planning.

Contribution

It introduces a combined numerical approach to identify optimal capture points and orbits in the Pluto-Charon system, extending previous work on capture strategies.

Findings

Identified optimal capture points around Charon.

Classified initial conditions into bounded, escaping, and collisional.

Compared capture dynamics between Pluto-Charon and Earth-Moon systems.

Abstract

In this paper, the near-optimal capture problem is numerically investigated to locate the optimal capture points around Charon using the Pluto-Charon planar circular restricted three-body problem (PCRTBP). Capture orbits are divided into the pre- and post-maneuver portions. In the pre-maneuver portion, the gravitational capture conditions are discussed by backward numerical integration. In the post-maneuver portion, the smaller alignment index (SALI) is applied to numerically investigate the orbital characters of long-term capture orbits around Charon. The initial conditions corresponding to three types of motion: (i) bounded (regular or chaotic), (ii) escaping and (iii) collisional are classified and studied. Combining results of the pre- and post-maneuver portions, the near-optimal capture method is presented to find optimal capture points by overlay figures, and then the corresponding capture orbits are constructed in the PCRTBP. The different results between the Pluto-Charon system and the Earth-Moon system are compared and analyzed. Our results obtained in this paper, including the optimal capture points and the corresponding capture orbits, could be applied in future space mission design.

Figures

Tables

Table 1. Table 1: Data of prograde near-optimal capture orbits.

Fig. 17	$(x,y)$	$C_{max}$	$C$	$\Delta V$ (m/s)
(a)	$(0.92575062,-0.00083565)$	3.687	3.90	10.4064
(b)	$(0.79538850,0.00083565)$	3.697	3.71	1.2940
(c)	$(0.99093168,0.15961003)$	3.627	3.60	6.0325
(d)	$(0.82547207,-0.00083565)$	3.677	$C_2$	7.3500

Table 2. Table 2: Data of retrograde near-optimal capture orbits.

Fig. 18	$(x,y)$	$C_{max}$	$C$	$\Delta V$ (m/s)
(a)	$(0.90736622,-0.03091922)$	3.097	$C_2$	22.6293
(b)	$(0.89733836,0.07103064)$	2.947	$C_4$	3.5683
(c)	$(0.88062527,-0.09442897)$	2.867	2.90	2.5561
(d)	$(0.80207374,0.27158774)$	2.777	2.80	3.1574

Equations

$$U_x(x,y)$$

$$U_y(x,y)$$

$$U(x,y)=\frac{m_1}{r_1}+\frac{m_2}{r_2}+\frac{\mu \left(1-\mu \right)}{2}+\frac{1}{2}\left(x^2+y^2\right),$$

$$r_1$$

$$r_2$$

$$\dot{(\delta x)}$$

$$\dot{(\delta y)}$$

$$\dot{(\delta z)}$$

$$(\dot{\delta \dot{x}})$$

$$(\dot{\delta \dot{y}})$$

$$(\dot{\delta \dot{z}})$$

$$J(x,y,\dot{x},\dot{y})=2U(x,y)-\left({\dot{x}}^2+{\dot{y}}^2\right)=C,$$

$$\begin{array}{rl}& V_1=\sqrt{2U(x,y)-C},\\ & \dot{x}=-k{V}_{1}y/r_2,\\ & \dot{y}=k{V}_1(x+\mu -1)/r_2,\end{array}$$

$$V_{min}=\sqrt{2U(x,y)-C_{max}}+k{r}_2.$$

$$e_{min}=\frac{r_2{V}_{min}^2}{\mu }-1.$$

$$f(x,C)=2U(x,y=0)=C.$$

$$\begin{array}{rl}\Delta V& =|V_1-V_2|\\ & =|\sqrt{2U(x,y)-C_{max}}-\sqrt{2U(x,y)-C}|\end{array}$$

Full text

††journal: Planetary and Space Science

1 Introduction

Capture dynamics is a hot topic in the astronomy and astronautics. Using Kolmogorov-Arnold-Moser (KAM) tori, Astakhov et al. (2003) proposed the chaos-assisted capture (CAC) mechanism in the circular restricted three-body problem (CRTBP) to investigate irregular moons of Jupiter and Saturn, which identifies chaos as the essential feature responsible for initial temporary gravitational trapping within a planet’s Hill sphere. Their theory can explain why Saturn has more prograde irregular moons than Jupiter. In the Earth-Moon system, asteroid 2006 RH120 was observed to be captured temporarily from September 2006 to June 2007. Granvik et al. (2012) first proposed a concept of temporary capture (TC) to describe and investigate this phenomenon in the Earth-Moon system. Based on their computation, at any given time there should be at least one temporarily-captured irregular natural Earth satellite (NES) of 1-m diameter orbiting the Earth. The well-known gravitational capture or weak stability boundary (WSB) theory was first proposed by Belbruno (1987) to solve the lunar capture problem in Earth-to-Moon transfers. Then, this theory has been successfully applied to many lunar missions, such as JAXA’s Hiten rescue (Belbruno & Miller, 1990), ESA’s SMART-1 (Schoenmaekers et al., 2001), and NASA’s GRAIL (Roncoli & Fujii, 2016). The in-depth research indicates that the WSB can be extended to interplanetary transfers, such as to Mercury (Jehn et al., 2004) and Mars (Topputo & Belbruno, 2015).

In addition, many scholars focus on the capture problem of the near-Earth asteroids (NEAs), because the NEAs play important roles in the origin of the solar system and the establishment of early warning of asteroid impact. Baoyin et al. (2010) proposed a method to search for low-energy NEAs that may be temporarily captured by the Earth. NASA’s Asteroid Redirect Robotic Mission (ARRM) was proposed to capture an entire $4\sim 10$ m NEA and redirect it into an orbit around the Moon where astronauts in the Orion spacecraft could explore it (Strange et al., 2013). Yarnoz et al. (2013) found a family of so-called Easily Retrievable Objects (EROs), which can be transported from accessible heliocentric orbits into the Earth’s neighborhood at affordable costs. Bao et al. (2015) applied lunar flybys and Earth flybys to capture NEAs into bounded orbits of the Earth. Based on the CAC mechanism in (Astakhov et al., 2003), Verrier & McInnes (2014) studied the low-energy capture of asteroids onto KAM tori in the Sun-Earth system. Qi & de Ruiter (2018) proposed a concept called short-term capture (STC) to investigate temporary capture phenomena of NEAs in the EMS. Space and time conditions of STC were defined and studied. Numerical results indicated that there was a clear association between the distributions of the time probability of STC and its space conditions.

The impulsive maneuver, $\Delta V$, provided by the propulsion system of the spacecraft is the most practical maneuver method in the space mission to the present day. Hence, the optimal capture strategy usually can be evaluated by the magnitude of the impulsive maneuver. Qi & Xu (2014) proposed the near-optimal lunar capture theory in the Earth-Moon planar CRTBP to locate the optimal capture point and achieve the permanent lunar capture using the minimum $\Delta V$. This method combined the lunar gravitational capture and KAM tori around the Moon. The minimum capture velocity and the minimum capture eccentricity before the maneuver were defined and analysed. The permanent lunar capture after the maneuver was investigated by the distribution of KAM tori around the Moon. Furthermore, based on the near-optimal lunar capture, Qi & Xu (2016) constructed low-energy transfers from the Earth to permanent lunar orbits. The numerical results indicated that compared with traditional patched-conic transfers, those low-energy transfers can efficiently save fuel costs, which actually were very close to the theoretical minimum.

In this paper, we apply the near-optimal capture theory to investigate the capture problem of Charon in the Pluto-Charon planar CRTBP. There are two main differences between this paper and the previous research of the Earth-Moon system in (Qi & Xu, 2014). Firstly, the mass ratio $\mu$ of the Pluto-Charon system is significantly larger than that of the Earth-Moon system, so the gravity of Charon is more significant in its region than the case of the Moon. Consequently, we postulate that the near-optimal capture in the Pluto-Charon planar CRTBP can show different characters and results. Secondly, Qi & Xu (2014) applied KAM tori to analyze the long-term orbital characters after the impulsive maneuver $\Delta V$. Although this method is simple and visible, the position and coverage of KAM tori are based on the empirical observation rather than the quantitative analysis. In this paper, the smaller alignment index (SALI), a common index for determining the regular or chaotic character of an orbit, is applied to investigate and classify the long-term orbits around Charon more accurately instead of KAM tori.

In Zotos (2015) the orbital dynamics in the planar circular Pluto-Charon system were numerically investigated. More specifically, large sets of initial conditions of orbits were classified in an attempt to determine the dynamics of a test particle (e.g., a spacecraft, or a comet, or an asteroid) moving in the vicinity of Charon. It was found that the configuration $(x,y)$ space is covered by a highly complicated mixture of several types of basins, such as basins of escape, basins of collision and basins of bounded motion. In addition, by varying the value of the Jacobi constant (or in other words the energy level) he managed to monitor the parametric evolution of all types of orbits. In this paper, the quantitative method similar to (Zotos, 2015) will be proposed to analyze the long-term orbital dynamics in the Pluto-Charon system. Based on Qi & Xu (2014), the minimum $\Delta V$ satisfying the permanent capture can be obtained from two perspectives: decreasing the pre-maneuver velocity $V_{1}$ and increasing the post-maneuver velocity $V_{2}$. Under given constraints of the capture point and Jacobi constant, there exist a minimum $V_{1}$ and a maximum $V_{2}$. Combining two aspects, the optimal capture point with the minimum $\Delta V=V_{1}-V_{2}$ can be found.

According to the above discussion, the structure of this paper is separated into five parts. Section 2 introduces the main properties of the dynamical system. Section 3 focuses on the gravitational capture conditions in pre-maneuver portion and obtain the minimum $V_{1}$. In Section 4, we investigate the long-term capture in post-maneuver potion and obtain the maximum $V_{2}$. After that, the near-optimal capture of Charon in the Pluto-Charon planar CRTBP is proposed and analyzed in Section 5. Finally, conclusions of our investigation are presented in Section 6.

2 Description of the dynamical model

According to the theory of the classical circular restricted three-body problem (Szebehely, 1967), the two primary bodies $P_{1}$ (Pluto) and $P_{2}$ (Charon) move on circular orbits (Buie et al., 2012), with the same angular velocity, around their common center of mass. The third body, also known as the test particle (e.g., a spacecraft, or a comet, or an asteroid), moves on the same plane under the resultant Newtonian gravitational field of the two main bodies. Taking into account that the mass of the test particle is extremely small, compared to the masses of the two primaries, we can reasonably assume that the third body does not perturb the motion of the primaries. The mass ratio $\mu=m_{2}/(m_{1}+m_{2})$, with $m_{1}>m_{2}$ defines the non-dimensional masses of the two primaries $m_{1}=1-\mu$ and $m_{2}=\mu$. In this notation, the centers of Pluto and Charon are located at $(-\mu,0)$ and $(1-\mu,0)$, respectively (see Fig. 1).

The planar motion of a test particle can be described by the following set of differential equations

[TABLE]

where

[TABLE]

is the time-independent effective potential function of the system, while

[TABLE]

are the distances to the respective primaries.

Similarly, the set of the variational equations is given by

[TABLE]

The set of the equations of motion (1) admits only the following integral of motion

[TABLE]

which is also known as the Jacobi integral. The numerical value $C$ of the Hamiltonian remains constant and it is called the Jacobi constant.

According to NASA Space Science Data Coordinated Archive (Update: 23 December 2016) the value of the mass ratio for the Pluto-Charon system is $\mu$ = 0.10851122058111It should be noted that in Zotos (2015), the value $\mu=0.099876695437731$ was adopted. However, over the last years the official value of the mass ratio regarding the Pluto-Charon planetary system has been revised., which is significantly larger than that of the Earth-Moon planetary system, 0.0121506683 (Qi & Xu, 2014). The values of the Jacobi constant at the five Lagrange points $(L_{i},i=1,...,5)$ are denoted by $C_{i}$, $i=1,...,5$ and they are equal to: $C_{1}$ = 3.717080, $C_{2}$ = 3.576079, $C_{3}$ = 3.204728, and $C_{4}=C_{5}$ = 3.

3 Gravitational capture

In this section, the gravitational capture in the pre-maneuver portion is investigated to decease the pre-maneuver velocity $V_{1}$. Similar to Qi & Xu (2014), we also assume that capture points are apsis points with respect to Charon, including periapsides and apoapsides. Since the velocity direction of the test particle at the capture point is perpendicular to the vector from the center of Charon to the test particle, the velocity direction has only two possibilities in the orbital plane: the prograde motion (anticlockwise motion with respect to Charon) or the retrograde motion (clockwise motion with respect to Charon).

The gravitational capture requires that the test particle starting from the capture point can leave the sphere of influence (SOI) of Charon under a given flight time by backward integration. The radius of the SOI of Charon is set as 10000 km (see Fig. 2), about double the distance of $L_{1}$ from the center of the smaller primary, which is similar to the lunar SOI proposed by Yamakawa et al. (1992). The flight time from the boundary of the SOI to the capture point is restricted to 15 time units (approximately 15.24 days). For the given position and velocity direction of the capture point, we postulate that a larger the velocity magnitude can make the escape from Charon (that is the Charon capture in the forward propagation) easier. Under the limitation of the flight time, we can use the following method to find the minimum $V_{1}$ for the given capture point.

The initial point of backward integration is an apsis point with respect to Charon chosen from a scattering region. The scattering region of initial points is an annulus region around Charon, whose inner and outer radii are equal to the radius of Charon and the distance of $L_{1}$ from Charon, i.e., 606.0 and 5848.3 km, respectively (see Fig. 2). Inside this annulus we define a uniform and dense grid of $1024\times 1024$ initial conditions. There are two kinds of direction of velocity: the prograde and retrograde motion, but its velocity magnitude in the Pluto-Charon rotating frame $V_{1}$ can be changed. 2. 2.

According to Eq. (5), $V_{1}=\sqrt{2U(x,y)-C}$, so the minimum $V_{1}$ can translate to the maximum $C$. For the given capture point $(x,y)$ and Jacobi constant $C$, its velocity in the Pluto-Charon rotating frame, $(\dot{x},\dot{y})$, can be obtained from

[TABLE]

where $k$ denotes the direction of the velocity: $k=1$ for the prograde motion with respect to Charon; $k=-1$ for the retrograde motion with respect to Charon. $r_{2}$ is the distance of the initial point to the center of Charon, and can be obtained from Eq. (3). In this way, the initial state of the capture point $(x,y,\dot{x},\dot{y})$ in the rotating frame can be derived. 3. 3.

We choose the initial $C$ as $C_{1}$. Then, using the method in step 2, we can obtain the corresponding initial state $(x,y,\dot{x},\dot{y})$. By backward integration, a capture orbit can be obtained from the initial state. If the orbit is bounded in the SOI within 15 time units or collides with Charon, we postulate that the given $C$ cannot satisfy the gravitational capture. Then we test whether the backward escape can be achieved for a smaller $C$. In this way, $C$ is decreasing until the backward escape can be satisfied. The corresponding maximum $C$ is defined as $C_{max}$.

It should be noted that the above numerical method is similar to that in (Qi & Xu, 2014), but in step 3, we use $C$ as the variable rather than $V_{1}$ in (Qi & Xu, 2014). The advantage of $C$ is that its initial value or upper boundary is specific, i.e. $C_{1}$.

In Fig. 3 we display the distributions of $C_{max}$ in the Pluto-Charon rotating frame, where left and right panels correspond to prograde and retrograde capture points, respectively. From Fig. 3, we discover that $C_{max}$s of the prograde motion are larger than those of the retrograde motion in most areas. The lower boundary of $C_{max}$ of the retrograde motion is significantly smaller than that of prograde motion. Compared with distributions of $C_{max}$ of Earth-Moon PCRTBP in (Qi & Xu, 2014), structures of $C_{max}$ for retrograde motion are very similar. For the prograde motion, the distribution of $C_{max}$ in the Earth-Moon PCRTBP (Qi & Xu, 2014) is more regular and consecutive, and there exist two branches with smaller $C_{max}$ in two sides of the Moon; however, in the Pluto-Charon PCRTBP, highly-fractal diffused structures of $C_{max}$ are more complicated in the Charon realm.

In the calculation of $C_{max}$, the integration time from the initial point to the boundary of SOI of Charon is denoted by $t_{esc}$. Fig. 4 shows the distributions of $t_{esc}$ in the Pluto-Charon rotating frame. Similarly, left and right panels correspond to prograde and retrograde capture points, respectively. By looking at the figure, it can be seen that in most region, orbits starting from prograde capture points are easier than those from the retrograde points to escape the Charon’s SOI by backward integration. In addition, we can find some similarities between structures of $t_{esc}$ and $C_{max}$ in Fig. 3.

According to Eq. (6), the minimum $V_{1}$ we seek can be obtained by $C_{max}$. It should be noted that this minimum velocity is a magnitude in the Pluto-Charon rotating frame rather than in the Charon-centered inertial frame. We define the corresponding minimum velocity in the Charon-centered inertial frame as $V_{min}$. Since capture points are apsides with respect to Charon, we can obtain the following relationship,

[TABLE]

Fig. 5 shows the distributions of $V_{min}$ in the Pluto-Charon dimensionless rotating frame. The left panel corresponds to prograde capture points, while the right one corresponds to retrograde capture points.

From Fig. 5, we observe that, generally, if the capture point is farther from Charon, the minimum capture velocity $V_{min}$ is smaller. This result is accordant with the fact: the attraction from Charon at a capture point farther to Charon is smaller, which results in a smaller escape velocity. As we can see from the figure, for both prograde and retrograde motions, minimum capture velocity $V_{min}$s are mainly distributed in the range of 0.28 $\sim$ 2.57. However, in Fig. 3 $C_{max}$ in the prograde orbits is obviously larger than that in the retrograde orbits in most areas. We can explain this phenomenon by Eq. (7): for the same $V_{min}$ at a given capture point, $C_{max}$ of prograde motion is larger than that of retrograde motion. $C_{max}$ is derived from the state of the rotating frame, so the influence of the Pluto-Charon rotation can influence $C_{max}$.

In the two-body model, we can use the orbital eccentricity $e$ to determine whether the orbit is elliptic or hyperbolic. In the vicinity of a primary body of the PCRTBP, since the gravity of that celestial body is dominant, the concept of the eccentricity $e$ can still be valid to study the orbital character of the test particle. However, it should be noted that in the PCRTBP, $e$ is time-varying. When $e$ is alternated from larger than 1 to smaller than 1, we postulate that gravitational capture happens. The minimum capture eccentricity $e_{min}$ is defined as the eccentricity at the capture point corresponding to the minimum $V_{min}$, and can be obtained from

[TABLE]

Fig. 6 illustrates distributions of the minimum capture eccentricity $e_{min}$ in the Pluto-Charon rotating frame. Based on the definition of the gravitational capture, capture points with $e_{min}>1$ are infeasible for the gravitational capture. We should avoid these regions when we choose the capture points. As we can see from the figure, in the right figure, all retrograde capture points are available, and the structures are similar to those of $C_{max}$. Qi & Xu (2014) found that structures of $e_{min}$ of the Earth-Moon PCRTBP for prograde motion have two branches in two sides of the Moon with $e_{max}>1$. However, as we can see from the left figure, the infeasible prograde capture points for gravitational capture are only distributed in a small region on the left edge.

4 Long-term capture

In this section, we focus on the maximum post-maneuver velocity $V_{2}$ satisfying the long-term capture around Charon. According to (Qi & Xu, 2014), a regular orbit in KAM tori around Charon can be bounded in the Charon realm permanently in the Pluto-Charon PCRTBP, otherwise, an orbit located in chaotic areas escapes the Charon realm. Therefore, we require that a test particle after $\Delta V$ should be inserted into a regular orbit around Charon. In (Qi & Xu, 2014), KAM tori were applied to find regular bounded regions. Although this method is simple and visible, the position and coverage of KAM tori rely on the empirical observation rather than the quantitative analysis. The SALI is a common index for quantitatively determining the regular or chaotic character of an orbit (Skokos, 2001). The analytical definition of SALI, as well as the particular threshold values for distinguishing between ordered and chaotic motion are given in Zotos & Jung (2018). The classifications of the long-term orbits around Charon can be obtained by the value of the SALI quantitatively, including (i) bounded (regular or chaotic), (ii) escaping and (iii) collision (Zotos, 2015). Essentially, other diagnostic indexes for the orbital chaos, such as the the fast Lyapunov indicator (FLI) (Astakhov & Farrelly, 2004), the correlation dimension (Belbruno et al., 2008), and the frequency analysis method (Laskar, 1990), can also be applied to distinguish stable regular orbits from chaotic ones. In this section, we use the SALI to investigate bounded basins of the long-term orbits around Charon.

In Fig. 7 we present the orbital structure of the configuration $(x,y)$ plane for the prograde motion. It is seen that when $C>C_{1}$ all the energetically allowed area of motion is covered either by collision or by regular initial conditions of orbits. As long as $C<C_{1}$ the channel around $L_{1}$ opens and therefore orbits can escape to the Pluto realm. The stable areas of regular initial conditions are destroyed abruptly. It should be noted that when $C=C_{2}$ two new stability islands emerge at both sides of Charon. We observe that as the value of the Jacobi constant decreases the amount of both regular and collision orbits decrease. On the other hand, the amount of escaping orbits grows rapidly and for $C\leq C_{4}$ they dominate the entire configuration plane. In the Earth-Moon PCRTBP (Qi & Xu, 2014), KAM tori rapidly disappear when $C$ decreases to $C_{2}$, but apparently, regular orbits in the Pluto-Charon system can still survive when $C=C_{2}$. We postulate that this difference is derived from the larger $\mu$ of Pluto-Charon system, because orbits rotating around the smaller primary with a larger $\mu$ have a greater resistance against the perturbation from the larger primary in the PCRTBP.

When $C<C_{1}$, the abrupt emergence of chaotic regions in Fig. 7 looks quite similar to the abrupt bifurcation to fully developed chaotic scattering in (Bleher et al., 1990). However, we postulate that the chaotic phenomenon in Fig. 7 essentially is different from the chaotic scattering. Firstly, the fully developed chaotic scattering denotes a situation in which all periodic orbits are unstable and there are no KAM surfaces (i.e., the dynamics is hyperbolic) (Bleher et al., 1989), but in Fig. 7, obviously, there exist KAM tori in the bounded regular region (denoted by green) when $C$ is smaller than $C_{1}$ (Qi & Xu, 2014). In addition, in the CRTBP we consider, not all periodic orbits are unstable, for example the distant retrograde orbits around the smaller primary are stable (Henon, 1969; Lam & Whiffen, 2005). Secondly, Koon et al. (2000) pointed out that the chaotic phenomenon in the region around the smaller primary of the CRTBP is ascribed to the existence of transversal homoclinic and heteroclinic points. We can explain these chaotic phenomena by the symbolic dynamics or the horseshoe-like dynamics (Koon et al., 2000).

It should be noted that trapped chaotic motion is very limited. This means that the corresponding initial conditions appear only as lonely and isolated points, mainly in the vicinity of the stability islands of regular motion. Indeed, in panel (b) of Fig. 7 one can distinguish several near the boundaries of the bounded basins of regular motion the presence of several initial conditions which correspond to trapped chaotic motion.

In the same vein, the orbital structure of the configuration plane for the retrograde motion is depicted in Fig. 8. It is evident that the orbital structure corresponding to the retrograde motion has many differences with respect to that of the prograde motion. The most notable ones are: (i) bounded regular motion is possible for much lower values of the Jacobi constant, or alternatively, for higher values of the total orbital energy; (ii) the overall orbital structure seems less complicated, without the presence of regions with highly fractal basin boundaries. Those results are similar to those of the Earth-Moon PCRTBP in(Qi & Xu, 2014).

The distributions of the corresponding escape and collision time of the orbits of both the prograde and retrograde motion are given in Figs. 9 and 10, respectively. As expected, the lowest values of both escape and collision time were measured inside the corresponding basins of escape and collision, while the highest times were identified mainly in the fractal basin boundaries, where the degree of the unpredictability of the system maximizes.

The color-coded diagrams of Figs. 7 and 8 reveal the orbital structure of the configuration $(x,y)$ plane for a fixed value of the Jacobi constant (or equivalently, for a fixed orbital energy), for both prograde and retrograde motion. In an attempt to obtain the same information for a continuous interval of $C$ values, we define a plane of representation in which the $x$-coordinate is the abscissa, while the ordinate is given by the Jacobi constant. This means that all orbits start with initial position $(x_{0},0)$, i.e, on the $x$-axis, with initial velocities $(0,\dot{y}_{0})$, where the $y$-component of the velocity is derived from the third of Eqs. (6).

In Fig. 11 we present the basin diagram on the $(x,C)$ plane for $C\in[2.5,5]$, for both the prograde and retrograde motions. The black continuous line indicates the corresponding Zero Velocity Curve (ZVC), and it is defined as

[TABLE]

We see that for the prograde motion bounded regular orbits exist only when $C>C_{2}$, while for the retrograde motion the same type of orbits is still possible even below $C=C_{4}$. However, for extremely high and low values of $C$ we encounter the same orbital structure for both types of motion (prograde and retrograde). In particular, for extremely low values of the Jacobi constant escaping orbits completely dominate, while for relatively high values of $C$ bounded and collision types of orbits share the energetically allowed phase space. In Fig. 12 we present the distribution of times employed by the test particle in the escape and collision orbits.

Based on the previous discussion, for both prograde and retrograde motions, bounded orbits generally erode or disappear with the decrease of $C$ (see Figs 7 and 8). Hence, for a given capture point $(x,y)$ and its direction of motion, there exists a minimum $C$ satisfying long-term capture in the Charon region. According to Eq. (5), $V_{2}$ corresponding to the minimum $C$ can be obtained, which is the maximum $V_{2}$ we seek for.

5 Near-optimal capture of Charon

In the previous two sections, we have investigated how to obtain the minimum $V_{1}$ and the maximum $V_{2}$, successively. In this section, we synthesize those results to analyse the minimum $\Delta V$. According to Qi & Xu (2014), the key to find the minimum $\Delta V$ is to locate the near-optimal capture point. Combining results obtained in previous sections, near-optimal capture points can be found by the following method:

For the given velocity direction, $\Delta V$ corresponding to the capture point $(x,y)$ can be calculated by

[TABLE]

where $C_{max}$ is obtained from Figs. 3 or 4, and $C$ is the Jacobi constant of the post-maneuver portion. It should be noted that if the given capture point $(x,y)$ is located in the energetically forbidden regions for the given $C$, Eq. (10) is invalid and the corresponding $\Delta V$ does not exist. 2. 2.

Using Eq. (10), the distribution of $\Delta V$ in the scattering region for given velocity direction and $C$ can be obtained. For example, Figs. 13 and 14 display the distributions of $\Delta V$ in the Pluto-Charon rotating frame corresponding to the cases presented in Figs. 7 and 8, respectively. In addition, for the same velocity direction and $C$, the corresponding basin diagrams can be obtained using the method in Section 4 (see Figs. 7 and 8). Then, we put the bounded basins onto the distribution of $\Delta V$, and get an overlay figure of them. For example, Figs. 15 and 16 display overlays for different $C$s and motions. The background of the overlay is the distribution of $\Delta V$, and cyan regions denote bounded basins. 3. 3.

From the overlay figures, points in bounded basins denote feasible points in terms of the long-term capture after $\Delta V$. Among these feasible capture points, those located in the minimum $\Delta V$ regions are chosen as near-optimal capture points. For example, blue points in Figs. 15 and 16 are near-optimal capture points we chose. The corresponding $\Delta V$ is the minimum $\Delta V$ we seek for.

Figs. 13 and 14 display the distributions of $\Delta V$ corresponding to Figs. 7 and 8, respectively. Gary regions in those figures are the energetically forbidden regions in the Charon realm for larger $C$s. Compared Fig. 13 with Fig. 3 (a), the structure of the distribution of $\Delta V$ for the prograde motion actually is quite similar to that of $C_{max}$ for prograde motion. It can be explained by Eq. (10). Similarly, the structure of the distribution of $\Delta V$ in Fig. 14 is similar to that of $C_{max}$ for retrograde motion in Fig. 3 (b).

Figs. 15 and 16 show the overlays of $\Delta V$ and the bounded basins of different $C$ for the prograde and retrograde motions, respectively. Gray and cyan regions in these figures denote energetically forbidden regions and bounded basins, respectively. Dark blue points are near-optimal capture points we chose. Based on the analysis of Fig. 7, only for a larger $C$ can bounded basins of the prograde motion exist. Hence, near-optimal prograde capture points will disappear with the decrease of $C$. But near-optimal retrograde capture points can survive in smaller $C$s. By looking at these figures, it can be seen that since distributions of $\Delta V$ and bounded basins of prograde motion are both quite different from those of retrograde motion, the distributions of the corresponding near-optimal capture points have significant differences between two kinds of motions. For example, in Fig. 16, near-optimal retrograde capture points are mostly located in the inner boundary of bounded basins near Charon, and similar to the distribution of bounded basins, near-optimal capture points gradually move away from Charon with the decrease of $C$. However, situations of near-optimal prograde capture points are more complicated according to Fig. 15, where regions with highly fractal basin boundaries exist, especially when $C=3.71$. In Fig. 15 (a), near-optimal capture points are distributed in the inner boundary of bounded basin nearest Charon. But when $C$ decreases to 3.71, i.e., necks around $L_{1}$ just open, near-optimal capture points are distributed in a highly fractal pattern around Charon, which is similar to the structures near Charon in Figs. 3 and 13. When $C$ decreases to 3.60, corresponding bounded basins are distributed in two separated stable islands, we can obtain two near-optimal capture regions in two separated bounded basins, respectively. But when $C$ decreases to $C_{2}$, near-optimal capture points can only be distributed in a small region of the left island. Numerical computation indicates that minimum $\Delta V$s corresponding to near-optimal capture points in Figs. 15 and 16 are all smaller than 23 m/s.

Using the near-optimal capture points as the initial state, the corresponding pre- and post-maneuver capture orbits can be constructed. In Figs. 17 and 18, we illustrate the capture orbits corresponding to the near-optimal capture points in Figs. 15 and 16, respectively. The figures are depicted in the dimensionless Pluto-Charon rotating frame. Red dashed lines denote pre-maneuver portions of capture orbits, while green solid lines denote post-maneuver portions of capture orbits (or long-term capture orbits). Red star points are the optimal capture points (or impulsive points) chosen from the near-optimal capture regions in Figs. 15 and 16 correspondingly. By looking at the figures, we observe that both pre- and post-maneuver portions of capture orbits have significant differences between prograde and retrograde motions. In Fig. 17, since the corresponding $C$s of prograde motion are larger than $C_{2}$, the sizes of post-maneuver orbits are smaller than those of retrograde motion in Fig. 18. In addition, all pre-maneuver portions of prograde near-optimal capture points pass through the gateway of $L_{1}$. However, in Fig. 18, most pre-maneuver portions of retrograde near-optimal capture points pass through the gateway of $L_{2}$. As we can see Fig. 18, multiple Charon flybys occur during the pre-maneuver portions. Data of prograde and retrograde near-optimal capture orbits in Figs. 17 and 18 are listed in Tables 1 and 2. From these tables, $\Delta V$s in Fig. 17 are smaller than 11 m/s, while $\Delta V$s in Fig. 18 are smaller than 23 m/s. Of particular note is that based on the data of Fig. 17 (c) and (d), $C_{max}$s of pre-maneuver portions are larger than the corresponding $C$s of post-maneuver portions, so the impulsive maneuvers at the near-optimal capture points increase the speed instead of brake.

The above results we obtained, including optimal capture points and corresponding capture orbits, could be applied in future space mission design (like the New Horizons interplanetary space probe). Although in real space missions there are several types of plane perturbations mainly due to the presence of other giant celestial bodies, such as Jupiter and Saturn, all these perturbations are extremely weak in the vicinity of Charon (which was the scattering region in our work). Therefore, we argue that the numerical results of the Pluto-Charon CRTBP obtained in this paper are structurally stable against out of plane perturbations.

6 Conclusions

In this paper, the near-optimal capture problem in the Pluto-Charon system was investigated and discussed. In an attempt to achieve long-term capture around Charon with the minimum $\Delta V$, optimal capture points were found and located by the near-optimal capture theory. According to the definition of $\Delta V$, the capture orbit was divided into the pre- and post-maneuver portions. In the pre-maneuver portion, the maximum escape Jacobi constant $C_{max}$ was obtained and discussed by backward numerical integration. In the post-maneuver portion, the SALI was applied to numerically investigate the orbital characters of long-term capture orbits around Charon. The initial conditions corresponding to three types of motion: (1) bounded, (2) escaping and (3) collisional were classified and distinguished. Combining the results of pre- and post-maneuver portions, optimal capture points were located by overlay figures. Furthermore, corresponding capture orbits were constructed in the PCRTBP. In addition, the near-optimal capture results in the Pluto-Charon system of this paper were compared with those in the Earth-Moon system. We postulated that the significant difference of mass ratio $\mu$ resulted in the different results between two planet-moon systems.

The main results of our numerical research can be summarized as follows:

Numerical computation indicated that the $C_{max}$s of the prograde motion are larger than those of the retrograde motion in most areas. Compared with distributions of $C_{max}$ of the Earth-Moon PCRTBP in (Qi & Xu, 2014), the structures of $C_{max}$ for prograde motion in the Pluto-Charon PCRTBP are more complicated in the Charon realm, where highly-fractal diffused structures of $C_{max}$ exist. Based on the results of $t_{esc}$, in most region, orbits starting from prograde capture points are easier than those from retrograde points to escape the Charon’s SOI by backward integration. 2. 2.

Based on the definition of gravitational capture and distributions of $e_{min}$, we found that all retrograde capture points are available for gravitational capture, and its structure of $e_{min}$ is similar to that of $C_{max}$. Different from the two branches in two sides of the Moon for prograde motion in the Earth-Moon PCRTBP (Qi & Xu, 2014), the prograde infeasible capture points for gravitational capture in the Pluto-Charon PCRTBP are only distributed in a small region on the left edge. 3. 3.

According to the orbital structure of long-term orbits around Charon, we found that the orbital structure corresponding to the retrograde motion has obvious differences with respect to that of the prograde motion: (1) bounded regular motion is possible for much lower values of the Jacobi constant, or alternatively, for higher values of the total orbital energy; (2) the overall orbital structure seems less complicated, without the presence of regions with highly fractal basin boundaries. In addition, the lowest values of both escape and collision time were measured inside the corresponding basins of escape and collision, while the highest times were identified mainly in the fractal basin boundaries, where the degree of the unpredictability of the system maximizes. 4. 4.

The distributions of near-optimal capture points indicated that near-optimal prograde capture points disappear rapidly with the decrease of $C$, but near-optimal retrograde capture points can survive in smaller $C$s. Near-optimal retrograde capture points are mostly located in the inner boundary of bounded basins near Charon, and gradually move away from Charon with the decrease of $C$. However, situations of prograde motion are more complicated, where highly fractal patterns around Charon and separated stable islands exist. 5. 5.

Numerical results indicated that all pre-maneuver portions of prograde near-optimal capture points pass through the gateway of $L_{1}$. However, most pre-maneuver portions of retrograde near-optimal capture points pass through the gateway of $L_{2}$. Data of both prograde and retrograde near-optimal capture orbits showed that, $\Delta V$s are smaller than 23 m/s.

For the numerical integration (both backward and forward in time) of the equations of motion (1), as well as of the variational equations (4) we used a double precision Bulirsch-Stoer FORTRAN 77 algorithm (Press et al., 1992), with a fixed time step equal to $10^{-3}$. In all our calculations the numerical error, related to the conservation of the Jacobi integral (5) was generally smaller than $10^{-12}$, while in most of the cases it was smaller than $10^{-14}$. All the graphical illustration has been created using the latest version 11.2 of the software Mathematica*®* (Wolfram, 2003).

The obtained results in this paper, including optimal capture points and corresponding capture orbits, could be applied in future space mission design.

Acknowledgments

Yi Qi wishes to acknowledge the support of the Natural Sciences and Engineering Research Council of Canada through a Discovery Accelerator Supplement under Grant RGPAS-493042-2016. The authors would like to express their warmest thanks to the two anonymous reviewers for the careful reading of the manuscript as well as for all the valuable suggestions and comments which allowed us to improve both the quality and the clarity of the paper.