A study of temporary captures and collisions in the Circular Restricted Three-Body Problem with normalizations of the Levi-Civita Hamiltonian

TL;DR

This paper investigates the dynamics near Lagrange points in the Circular Restricted Three-Body Problem, focusing on how regularization techniques extend the analysis of tube manifolds and their role in temporary captures of celestial bodies.

Contribution

It introduces a novel analysis of Hamiltonian normalizations using Levi-Civita regularization, overcoming amplitude limits of Cartesian variables and revealing structural transitions in tube manifolds.

Findings

Regularizations extend the normalization amplitude beyond Cartesian limits.

Structural transitions in tube manifolds occur at higher energies.

Orbits can collide with the secondary body before full circulation.

Abstract

The dynamics near the Lagrange equilibria $L_1$ and $L_2$ of the Circular Restricted Three-body Problem has gained attention in the last decades due to its relevance in some topics such as the temporary captures of comets and asteroids and the design of trajectories for space missions. In this paper we investigate the temporary captures using the tube manifolds of the horizontal Lyapunov orbits originating at $L_1$ and $L_2$ of the CR3BP at energy values which have not been considered so far. After showing that the radius of convergence of any Hamiltonian normalization at $L_1$ or $L_2$ computed with the Cartesian variables is limited in amplitude by $\|1-\mu-x_{L_1}\|$ ($\mu$ denoting the reduced mass of the problem), we investigate if regularizations allow us to overcome this limit. In particular, we consider the Hamiltonian describing the planar three-body problem in the Levi-Civita regularization and we compute its normalization for the Sun-Jupiter reduced mass for an interval of energy which overcomes the limit of Cartesian normalizations. As a result, for the largest values of the energy that we consider, we notice a transition in the structure of the tubes manifolds emanating from the Lyapunov orbit, which can contain orbits that collide with the secondary body before performing one full circulation around it. We discuss the relevance of this transition for temporary captures.

Figures

Tables

Table 1. Table 1: Norm of the remainder R N superscript 𝑅 𝑁 R^{N} evaluated on the Lyapunov orbit for E a = − 1.516 subscript 𝐸 𝑎 1.516 E_{a}=-1.516 (left columns) and for E b = − 1.509 subscript 𝐸 𝑏 1.509 E_{b}=-1.509 (right columns), according to Eq. ( 32 ), produced with expansions limited to order 𝒩 p = 27 subscript 𝒩 𝑝 27 {\cal N}_{p}=27 and normalization order N = 3 , … , 24 𝑁 3 … 24 N=3,...,24 .

Red orbit Fig. 4		Blue orbit Fig. 4
$N$	$R^N$	X $N$	$R^N$
$3$	$1.61835\times {10}^{[}-5]$	X $3$	$7.28121\times {10}^{[}-4]$
$4$	$2.35741\times {10}^{[}-6]$	X $4$	$3.53145\times {10}^{[}-4]$
$5$	$4.90253\times {10}^{[}-7]$	X $5$	$2.61175\times {10}^{[}-4]$
$6$	$1.00808\times {10}^{[}-7]$	X $6$	$1.98361\times {10}^{[}-4]$
$7$	$2.34993\times {10}^{[}-8]$	X $7$	$1.45230\times {10}^{[}-4]$
$8$	$4.72191\times {10}^{[}-9]$	X $8$	$1.04454\times {10}^{[}-4]$
$9$	$2.13184\times {10}^{[}-9]$	X $9$	$8.94943\times {10}^{[}-5]$
$10$	$4.13041\times {10}^{[}-10]$	X $10$	$6.06359\times {10}^{[}-5]$
$11$	$1.51257\times {10}^{[}-10]$	X $11$	$4.88443\times {10}^{[}-5]$
$12$	$3.59044\times {10}^{[}-11]$	X $12$	$3.31315\times {10}^{[}-5]$
$13$	$8.67395\times {10}^{[}-12]$	X $13$	$2.28439\times {10}^{[}-5]$
$14$	$1.81004\times {10}^{[}-12]$	X $14$	$1.53451\times {10}^{[}-5]$
$15$	$5.38462\times {10}^{[}-13]$	X $15$	$1.14096\times {10}^{[}-5]$
$16$	$9.09849\times {10}^{[}-14]$	X $16$	$7.54539\times {10}^{[}-6]$
$17$	$2.32918\times {10}^{[}-14]$	X $17$	$5.71248\times {10}^{[}-6]$
$18$	$5.61067\times {10}^{[}-15]$	X $18$	$4.31952\times {10}^{[}-6]$
$19$	$2.70410\times {10}^{[}-15]$	X $19$	$3.64058\times {10}^{[}-6]$
$20$	$3.35487\times {10}^{[}-16]$	X $20$	$2.42248\times {10}^{[}-6]$
$21$	$1.67395\times {10}^{[}-16]$	X $21$	$2.22411\times {10}^{[}-6]$
$22$	$4.91488\times {10}^{[}-17]$	X $22$	$1.63076\times {10}^{[}-6]$
$23$	$1.75187\times {10}^{[}-17]$	X $23$	$1.31428\times {10}^{[}-6]$
$24$	$5.80105\times {10}^{[}-18]$	X $24$	$9.80338\times {10}^{[}-7]$

Table 2. Table 2: Normal form decomposition of K 3 subscript 𝐾 3 {K_{3}} for the example L ​ L 1 ​ ( E a ) 𝐿 subscript 𝐿 1 subscript 𝐸 𝑎 LL_{1}(E_{a})

monomial				strategy
$k_1$	$l_1$	$k_2$	$l_2$	(a)	(b)	(c)	${\alpha }_{k_1,l_1,k_2,l_2}^{\left(3\right)}$
0	0	0	3	x	x	x	$-0.239846$
0	1	0	2	x	x	x	$0.1371133+0.0334760i$
0	3	0	0	x	x	x	$0.144918$
1	0	0	2	x	x	x	$0.0334760+0.1371133i$
1	2	0	0	x	x	x	$-0.164176i$
2	1	0	0	x	x	x	$0.164176$
3	0	0	0	x	x	x	$-0.144918i$
0	0	1	2	x	x	x	$0.642141$
0	1	1	1	x	x	x	$-0.284350$
1	0	1	1	x	x	x	$-0.284350i$
0	0	2	1	x	x	x	$-0.642141$
0	1	2	0	x	x	x	$0.1371133-0.0334760i$
1	0	2	0	x	x	x	$-0.0334760+0.1371133i$
0	0	3	0	x	x	x	$0.239846$
Total of terms				6	10	14

Table 3. Table 3: Normal form decomposition of K 4 subscript 𝐾 4 {K_{4}} for the example L ​ L 1 ​ ( E a ) 𝐿 subscript 𝐿 1 subscript 𝐸 𝑎 LL_{1}(E_{a})

monomial				strategy
$k_1$	$l_1$	$k_2$	$l_2$	(a)	(b)	(c)	${\alpha }_{k_1,l_1,k_2,l_2}^{\left(4\right)}$
0	0	0	4	x	x	x	$-0.478035$	†
0	1	0	3	x	x	x	$0.732257+1.135126i$	†
0	2	0	2	x	x	x	$-0.342835-1.071743i$	†
0	4	0	0	x	x	x	$0.0583902$	†
1	0	0	3	x	x	x	$1.135126+0.732257i$	†
1	1	0	2	x	x	x	$-2.41906i$	†
1	3	0	0	x	x	x	$-1.81534i$	†
2	0	0	2	x	x	x	$0.342835-1.071743i$	†
2	2	0	0	x	x	x	$-3.37708$	†
3	1	0	0	x	x	x	$1.81534i$	†
4	0	0	0	x	x	x	$0.0583902$	†
0	0	1	3	x	x	x	$2.23041$	†
0	1	1	2	x	x	x	$-2.86948-1.59924i$	†
0	2	1	1	x	x	x	$2.31725$	⋆
1	0	1	2	x	x	x	$-1.59924-2.86948i$	†
1	1	1	1	x	x	x	$2.55115i$	⋆
2	0	1	1	x	x	x	$-2.31725$	⋆
0	0	2	2	x	x	x	$3.75798$	∙
0	1	2	1	x	x	x	$-2.86948-1.59924i$	†
0	2	2	0	x	x	x	$-0.342835+1.071743i$	†
1	0	2	1	x	x	x	$-1.59924+2.86948i$	†
1	1	2	0	x	x	x	small $-2.41906i$	†
2	0	2	0	x	x	x	$0.342835+1.071743i$	†
0	0	3	1	x	x	x	$2.23041$	†
0	1	3	0	x	x	x	$-0.732257+1.135126i$	†
1	0	3	0	x	x	x	$1.135126-0.732257i$	†
0	0	4	0	x	x	x	$-0.478035$	†
Total of terms				9	23	27

Table 4. Table 4: Normalization information for the orbits presented in Fig.  4 .

$E$	$N_p$	$N$	term ${𝒦}_E^{\left(N\right)}$	terms ${ℛ}_E^{\left(N\right)}$	Max ampl
$-1.518$	13	10	$1658$	$2056$	$2.38\times {10}^{[}-2]$
$-1.516$	13	10	$1658$	$2056$	$3.83\times {10}^{[}-2]$
$-1.512$	19	16	$6967$	$5335$	$6.53\times {10}^{[}-2]$
$-1.509$	21	18	$10200$	$6924$	$8.54\times {10}^{[}-2]$
$-1.507$	25	22	$19869$	$10990$	$9.73\times {10}^{[}-2]$
$-1.505$	33	30	$57843$	$23346$	$1.08\times {10}^{[}-1]$
$-1.504$	33	30	$57843$	$23346$	$1.17\times {10}^{[}-1]$

Equations

$$T_P=\frac{a_P}{a}+2\sqrt{\frac{a}{a_P}(1-e^2)}\cos i$$

$$h=\frac{p_x^2+p_y^2+p_z^2}{2}+p_{x}y-p_{y}x-\frac{\mu }{\sqrt{(x-1+\mu {)}^2+y^2+z^2}}$$

$$-\frac{1-\mu }{\sqrt{(x+\mu {)}^2+y^2+z^2}}\text{leavevmode}.$$

$$T_P=-2h+\mathcal{O}(\mu )\text{leavevmode}.$$

$$(x,y,p_x,p_y)=C_N(Q_1,Q_2,P_1,P_2)$$

$$h=$$

$$-\frac{1-\mu }{\sqrt{(x+\mu {)}^2+y^2}}$$

$$H=E_i+\sum _{j=2}^N{K}_j(Q_1,Q_2,P_1,P_2)+\sum _{j\ge N+1}{\mathcal{R}}_j(Q_1,Q_2,P_1,P_2)\text{leavevmode}.$$

$$K_2=\lambda Q_1{P}_1+\omega \frac{P_2^2+Q_2^2}{2}\text{leavevmode},$$

$${LL}_{i}(E)=\Big{\{}Q,P:\ \ Q_{1},P_{1}=0,\ \ E_{i}+\sum_{j=2}^{N}{K}_{j}(0,Q_{2},0,P_{2})=E\Big{\}}\leavevmode\nobreak\ ,$$

$$\displaystyle{W}^{s}_{i}(E)=\Big{\{}Q,P:\ \ Q_{1}=0,\ \ E_{i}+\sum_{j=2}^{N}{K}_{j}(0,Q_{2},P_{1},P_{2})=E\Big{\}},$$

$$\displaystyle{W}^{u}_{i}(E)=\Big{\{}Q,P:\ \ P_{1}=0,\ \ E_{i}+\sum_{j=2}^{N}{\cal K}_{j}(Q_{1},Q_{2},0,P_{2})=E\Big{\}},$$

$$X=x-x_2\text{leavevmode},Y=y\text{leavevmode},P_X=p_x\text{leavevmode},P_Y=p_y-x_2\text{leavevmode}\text{leavevmode},$$

$$X=u_1^2-u_2^2\text{leavevmode},Y=2u_1{u}_2\text{leavevmode},$$

$$\text{PX=2∣u∣2U1u1−U2u2}\text{leavevmode},\text{PY=2∣u∣2U1u2+U2u1}\text{leavevmode},$$

$$dt=|u{|}^{2}d\tau$$

$${\mathcal{K}}_E(u,U)=$$

$$\displaystyle=\left|{u}\right|^{2}\left[h\Big{(}X(u)+x_{2},Y(u),P_{X}(u,U),P_{Y}(u,U)+x_{2}\Big{)}-E\right]$$

$$={\left|u\right|}^2\left[h\left(\text{u12−u22+x2},\text{2u1u2},\text{2∣u∣2U1u1−U2u2},\text{2∣u∣2+x2U1u2+U2u1}\right)-E\right]$$

$$=\frac{1}{8}{\left(U_1+2{\left|u\right|}^2{u}_2\right)}^2+\frac{1}{8}{\left(U_2-2{\left|u\right|}^2{u}_1\right)}^2$$

$$-\frac{1}{2}{\left|u\right|}^6-\mu -|u{|}^2\left(E+\frac{(1-\mu {)}^2}{2}\right)$$

$$-(1-\mu )|u{|}^2\left[\frac{1}{\sqrt{1+2(u_1^2-u_2^2)+|u{|}^4}}+u_1^2-u_2^2\right]\text{leavevmode},$$

$$\frac{d{u}_j}{d\tau }=\frac{\partial }{\partial U_j}{\mathcal{K}}_E(u,U),\frac{d{U}_j}{d\tau }=-\frac{\partial }{\partial u_j}{\mathcal{K}}_E(u,U)$$

$$t={\int }_0^{\tau }{\left|u(s)\right|}^{2}ds,$$

$$(u_{\ast },U_{\ast }):=(0,u_{2_{\ast }}^{},-2{\left(u_{2_{\ast }}^{}\right)}^3,0)\text{leavevmode},$$

$$F($$

$$+(1-\mu )\left(4(u_{2_{\ast }}^{}{)}^2-\frac{2(u_{2_{\ast }}^{}{)}^2}{(1-(u_{2_{\ast }}^{}{)}^2{)}^2}-\frac{2}{1-(u_{2_{\ast }}^{}{)}^2}\right)\text{leavevmode}.$$

$$E_{\ast }=-\frac{3}{2}\left(1-\frac{4}{3}\mu +\frac{1}{3}{\mu }^2\right).$$

$$(u_1,u_2,U_1,U_2)=A(\xi ,q,\eta ,p)+(u_{\ast ,1},u_{\ast ,2},U_{\ast ,1},U_{\ast ,2})$$

$${\mathcal{K}}_E^{(2)}(\xi ,q,\eta ,p)=K_0+K_2(\xi ,q,\eta ,p)+\sum _{j\ge 3}{k}_j(\xi ,q,\eta ,p)\text{leavevmode},$$

$$K_2=\lambda \xi \eta +i\omega qp\text{leavevmode},$$

$${\alpha }_{k_1,l_1,k_2,l_2}^{(2)}{q}^{k_1}{p}^{l_1}{\xi }^{k_2}{\eta }^{l_2}.$$

$${\mathcal{K}}_E^{(N)}$$

$$=K_0+\lambda \xi \eta +i\omega qp+\sum _{j=3}^N{K}_j(\xi ,q,\eta ,p)$$

$$+\sum _{j\ge N+1}{\mathcal{R}}_j^{(N)}(\xi ,q,\eta ,p),$$

$${\alpha }_{k_1,l_1,k_2,l_2}^{(N)}{q}^{k_1}{p}^{l_1}{\xi }^{k_2}{\eta }^{l_2},$$

Full text

\JournalInfo

Accepted Manuscript - International Journal of Non-Linear Mechanics \Archive \PaperTitleA study of temporary captures and collisions in the Circular Restricted Three-Body Problem with normalizations of the Levi-Civita Hamiltonian \PaperShortTitleCaptures and collisions in the CRTBP with normalizations of the Levi-Civita Hamiltonian \AuthorsRocío I. Paez1, Massimiliano Guzzo1* \Keywordscelestial mechanics, captures, CR3BP \AbstractThe dynamics near the Lagrange equilibria $L_{1}$ and $L_{2}$ of the Circular Restricted Three-body Problem has gained attention in the last decades due to its relevance in some topics such as the temporary captures of comets and asteroids and the design of trajectories for space missions. In this paper we investigate the temporary captures using the tube manifolds of the horizontal Lyapunov orbits originating at $L_{1}$ and $L_{2}$ of the CR3BP at energy values which have not been considered so far. After showing that the radius of convergence of any Hamiltonian normalization at $L_{1}$ or $L_{2}$ computed with the Cartesian variables is limited in amplitude by $\left|1-\mu-x_{L_{1}}\right|$ ($\mu$ denoting the reduced mass of the problem), we investigate if regularizations allow us to overcome this limit. In particular, we consider the Hamiltonian describing the planar three-body problem in the Levi-Civita regularization and we compute its normalization for the Sun-Jupiter reduced mass for an interval of energy which overcomes the limit of Cartesian normalizations. As a result, for the largest values of the energy that we consider, we notice a transition in the structure of the tubes manifolds emanating from the Lyapunov orbit, which can contain orbits that collide with the secondary body before performing one full circulation around it. We discuss the relevance of this transition for temporary captures.

1 Introduction

In the last decades the close encounters of a small body with a planet have been investigated especially in connection with the dynamics of comets, of near-earth asteroids, and space mission design (see [1], [2], [3], [4], [5], [6], [7], [8], [9], [10], [11], [12], [13], [14], [15], [16], [17], [18], [19], [20], [21], [22], and references therein). A classical tool to classify a close encounter is given by the Tisserand parameter respect to a given planet,

[TABLE]

($a,e$ denote the semi-major axis and the eccentricity of the small body; $i$ denotes the inclination of the small body with respect to the orbit of the planet; $a_{P}$ denotes the semi-major axis of the planet) whose variation before and after each encounter is small. The Tisserand parameter gives a qualitative idea of the encounter: for increasing values of $T_{P}$ in a small interval around $T_{P}=3$, we have the transition between the ’fast’ encounters, occurring with an orbit of the small body which is hyperbolic in the planetocentric reference frame, and

©$<$2020$>$. This manuscript version is made available under the

CC-BY-NC-ND 4.0 license http://creativecommons.org/licenses/by-nc-nd/4.0/

the ’slow’ close encounters which can lead to a temporary capture of the small body. Öpik theory ([23], recently revisited in [24], [25]), provided good results in the study of the fast close encounters of comets with Jupiter and of near-Earth asteroids with the Earth ([26],[8]). The slow close encounters are instead better studied in the framework of the dynamics generating at the Lagrangian points $L_{1},L_{2}$ of the Circular Restricted Three–Body Problem, defined by the Hamiltonian

[TABLE]

[TABLE]

This Hamiltonian is written in the barycentric rotating reference frame and with the usual units of measure for this problem: the masses of the primaries $P_{1}$ and $P_{2}$ are $1-\mu$ and $\mu$ respectively; their coordinates are $(x_{1},0,0)=(-\mu,0,0)$, $(x_{2},0,0)=(1-\mu,0,0)$ and their revolution period is $2\pi$. The Hamiltonian $h$ and the Tisserand parameter are related by

[TABLE]

For energy values $E>E_{1}$, where $E_{1}$ is the value of $h$ associated to the Lagrangian point $L_{1}$, the transits from the realm ${\cal S}$ of motions dominated by the Sun to the realm ${\cal P}$ of motions dominated by the planet (i.e. the temporary captures) and back, become possible (see [1], [27], [28], [22] for precise characterizations of the different realms of motion). Respectively for $E>E_{2}$, where $E_{2}$ is the value of $h$ associated to $L_{2}$, the transits from the realms ${\cal E}$ of motions which are external to the binary system to ${\cal P}$ (i.e. the temporary captures from the region external to the planet orbit) and back are possible. In the planar CR3BP (useful to study close encounters with small inclination), [1] has shown that the orbits which perform the transits are contained in two dimensional surfaces of the phase-space, the so-called tube manifolds of $L_{1}$ and $L_{2}$. These tube manifolds are defined as follows. Let us consider for definiteness the tube manifolds at $L_{1}$; for values of $E$ slightly larger than $E_{1}$ there is a periodic orbit of libration around $L_{1}$, the horizontal Lyapunov orbit denoted by $LL_{1}(E)$, whose amplitude increases rapidly as the value of $E$ increases. All the phase-space orbits which are asymptotic in the future (resp. in the past) to $LL_{1}(E)$ form a surface which, close to $LL_{1}(E)$, is topologically a 2-dimensional tube extending on both right and left sides of $LL_{1}(E)$ called the stable tube manifold $W^{s}_{1}(E)$ (resp. the unstable manifold $W^{u}_{1}(E)$) of $LL_{1}(E)$. Analogously, for $E$ slightly larger than $E_{2}$ we have the stable and unstable tube manifolds of $LL_{2}(E)$. The tube manifolds of $LL_{1}(E)$ separate the motions which transit between ${\cal S}$ and ${\cal P}$: for example, an orbit with $E>E_{1}$ which is a temporary satellite of the Sun, when it approaches the Lyapunov orbit $LL_{1}(E)$ transits to the realm ${\cal P}$ only if it is contained in the stable tube $W^{s}_{1}(E)$ of $LL_{1}(E)$, otherwise it bounces back to the realm ${\cal S}$. Therefore, numerical computations of the stable and unstable tube manifolds at various values of $E$ provide the relevant information to understand the transit properties related to the close encounters, e.g. to determine if a comet becomes a temporary satellite of Jupiter. A list of comets which have been identified as potential candidates for temporary captures can be found in [29]. For this reason, the Sun-Jupiter case, with mass ratio, $\mu=\mu_{J}:=9.537{\scriptstyle\times 10}^{[}-4]$, has received particular relevance in the literature (e.g. see [11], [21], [22]). Also in this paper we focus on $\mu=\mu_{J}$, while the methods that we use can be implemented with any other value. In particular, we correlate a property of the tube manifolds to a property of temporary captures which has been little considered in the literature: the number ${\cal N}$ of revolutions performed around Jupiter during the temporary capture and their orientation (clockwise or counter-clockwise) measured in the rotating reference frame. For simplicity, we consider the right-branch of the stable tube of the Lyapunov orbits $LL_{1}(E)$, and we identify three–situations:

• (i)

for $E=E_{1}$, the tube collapses to only one limit orbit (see Fig. 1, top-left panel);

• (ii)

for very small and positive values $E-E_{1}$, the numerical computation of the tube manifolds provides evidence that up to a fixed number ${\cal N}$ of revolutions all the orbits of the tube are not collision orbits, and perform the same number and type of revolutions

around Jupiter as the limit orbit of (i). As a consequence, also the orbits in the interior of the tube share the same properties. In particular, since collisions are excluded within the ${\cal N}$ revolutions, all the orbits in the interior of the stable tube transit to the realm ${\cal S}$ after the ${\cal N}$ revolutions (see Fig. 1, top-right panel, for ${\cal N}=2$).

• (iii)

The question is what happens for increasing values of $E-E_{1}$, which corresponds to increasing amplitudes of the stable/unstable tubes. When these amplitude are large, we find that the tube manifolds become so large and stretched in phase-space that the criterion of using them as separatrices for the transit properties is no more as effective as when their amplitude is small. In this paper we propose a threshold on $E$ based on the appearance of peculiar orbits in the stable tube. This orbit originates from a collision with Jupiter and converges to the Lyapunov orbit $LL_{1}(E)$ before performing a full revolution around the planet. In addition to this peculiar orbit, we find that the stable tube contains also orbits performing either clockwise or counter-clockwise revolutions (see Fig. 1, bottom panels). The same properties are shared by the orbits in the interior of the tube. In particular, due to possible collisions, it is no more granted that all the orbits in the interior of the tube will transit to the realm ${\cal S}$.

The conclusions outlined above have been obtained thanks to a method of computation of the tube manifolds which exploits the well known Levi-Civita regularization of the three-body problem and the method of normalization of an Hamiltonian at a partially hyperbolic equilibrium. The combination of the two techniques allows us to reach values of the energy never investigated before.

There are several methods for computing numerically the stable and unstable manifolds of periodic orbits, such as the flow continuation of the local manifolds, the parametrization method, and the recent method based on chaos indicators. In the first two the manifolds are developed from analytic approximations of the local stable and unstable manifolds (see for example [30]), in the latter the manifolds are obtained as the ridges of a chaos indicator defined from a Hamiltonian normalization ([17], [21], [22]). We remark that Hamiltonian normalizations provide not only high precision computations of the local stable/unstable manifolds, but also of all the orbits in their neighbourhood, suitable for astronomical applications. For this reason, in this paper we push the method of Hamiltonian normalization to its limit, by implementing it in the Levi-Civita regularization of the CR3BP.

Even if the Hamiltonian normalization method has been extensively used to compute the tube manifolds, its limits have to be improved, especially if one aims to extend its application to a broader interval of energies. In the present work, we show that the use of the Levi-Civita regularization is necessary in order to normalize the Hamiltonian at values of the energy for which we detect the transition from situation (ii) to (iii) described above.

In Section 2 we describe all the steps necessary to perform the normalization of the Levi-Civita Hamiltonian. In particular, we find that for $E>E_{1}$ the regularized Hamiltonian still has an equilibrium, not corresponding to an orbit of the CR3BP. Nevertheless the Levi-Civita Hamiltonian can be normalized at these ’fictitious’ equilibria, providing the Lyapunov orbit as well as their stable and unstable manifolds. In Section 3 we discuss the efficiency of the normal form computations of $L_{1}$. In Section 4, we show the computations of the manifolds via the normalized Levi-Civita Hamiltonian, for a wide range of energies such that the amplitude of the corresponding Lyapunov orbit is larger than $\left|1-\mu-x_{L_{1}}\right|$. While for small $E$ the stable tube manifold ${W}^{s}_{1}(E)$ folds around $P_{2}$ exclusively in a clockwise fashion (when integrated backwards in time for ${\cal N}=2$), for large values of $E$ the folding can be either clockwise or counter-clockwise. It is within this range of energy values where we identify the transition between the orbits discussed at ii) and iii).

2 Hamiltonian normalizations

Powerful methods to analyze the dynamics originating at the Lagrangian points $L_{1},L_{2}$ of the CR3BP rely on the Birkhoff normalizations of the Hamiltonian with a large normalization order $N$ ([7], [4], [31]). For the planar problem, this implies the explicit construction, for any value of an integer parameter $N\geq 3$, of a canonical transformation

[TABLE]

conjugating the Hamiltonian of the planar circular restricted three–body problem

[TABLE]

to a normal form Hamiltonian111For convenience, we do not simplify from the Hamiltonian the constant term ${E}_{i}$.

[TABLE]

The Hamiltonian in (4) is analytic in some neighbourhood of the Lagrange equilibrium $L_{i}$, represented by $(Q,P)=(0,0,0,0)$ (the radius of the neighborhood depending on the mass ratio $\mu$ and $N$), $K_{j},{\cal R}_{j}$ are polynomials in $P,Q$ of order $j$, the polynomials $K_{j}$ depend on the $Q_{1},P_{1}$ only through the product $Q_{1}P_{1}$ and

[TABLE]

(alternative reduction methods can be considered, see Section 3 for details). By neglecting the remainder terms ${\cal R}_{j}$ we obtain

• •

The approximated equations of the Lyapunov orbits labeled by $E$

[TABLE]

• •

the approximated equations of the (local) stable and unstable manifolds of ${LL}_{i}(E)$

[TABLE]

• •

the arcs of approximated orbits in a neighbourhood of ${W}^{s}_{i}(E)$, ${W}^{u}_{i}(E)$ which are scattered by $LL_{i}(E)$.

For increasing values of $E$, the periodic orbit ${LL}_{i}(E)$ has increasing libration amplitude. As a consequence, for suitably large values of $E$, we expect a breakdown of the method. A natural limit for this breakdown is given by the singularity in the gravitational potential energy at the position of $P_{2}$ which, for small values of $\mu$, is close to $L_{1},L_{2}$. Since the Lyapunov orbits ${LL}_{i}(E)$ are typically larger in the $y$ coordinate with respect to the $x$ coordinate, a dangerous complex singularity is the one located at $x=x_{L_{j}}$ and $y=\pm i\left|1-\mu-x_{L_{j}}\right|$ (for $\mu=\mu_{J}$, $y=\pm i0.066...$). This singularity severely limits the validity of these methods for libration amplitudes in the $y$ variable of order of $\left|1-\mu-x_{L_{j}}\right|$. We here investigate the possibility to overcome this limit by implementing the normalization methods using the Levi-Civita regularization ([32]) on the secondary body $P_{2}$. Similar approaches has been used in the past for the simpler Hill’s problem, e.g. [33].

Following [32], we first perform the phase-space translation

[TABLE]

on Hamiltonian (3), and then we introduce the Levi-Civita variables $(u_{1},u_{2})$ canonically extended to the momenta $(U_{1},U_{2})$

[TABLE]

and the fictitious time $\tau$

[TABLE]

where $|{u}|^{2}=u_{1}^{2}+u_{2}^{2}=\sqrt{(x-1+\mu)^{2}+y^{2}}$. For any fixed value of $E$ of the Hamiltonian (3), we define the Levi-Civita Hamiltonian

[TABLE]

which is a regularization of the planar three-body problem at $P_{2}$. In fact, ${\cal K}_{E}$ is regular at $u=(0,0)$ which correspons to a collision with $P_{2}$. The solutions $({u}(\tau),{U}(\tau))$ of the Hamilton equations of ${\cal K}_{E}$,

[TABLE]

with initial conditions222The initial conditions of the regularized variables correspond to initial conditions of the original barycentric reference frame such that $h(x(0),y(0),p_{x}(0),p_{y}(0))=E$. satisfying ${u}(0)\neq 0$ and

${\cal K}_{E}({u}(0),{U}(0))=0$, are conjugate, in a neighbourhood of $\tau=0$, via Eq. (8) and

[TABLE]

to solutions $(X(t),Y(t),P_{X}(t),P_{Y}(t))$ of the three-body problem Hamiltonian (see [34]). For definiteness, we describe the normalization of the Levi-Civita Hamiltonian at the Lagrangian point $L_{1}$; a similar method applies to $L_{2}$. Equilibrium points of ${\cal K}_{E}$. The regularized Hamiltonian ${\cal K}_{E}$, for increasing values of $E\geq E_{1}$, has a family of equilibria continued from $L_{1}$. Since the equilibria of the family are characterized by $u_{1},U_{2}=0$ and $U_{1}=-2u_{2}^{3}$, we define the fictitious equilibrium position as

[TABLE]

where $u_{2_{*}}^{\phantom{3}}\!:=u_{2_{*}}^{\phantom{3}}\!(E)$ solves the algebraic equation

$F(u_{2_{*}}^{\phantom{3}}\!(E),E)=0$ with

[TABLE]

In Appendix 1 (see paragraphs A1 and A2) we prove that equation $F(u_{*,2},E)=0$ has a unique solution $u_{2_{*}}^{\phantom{3}}\!$ for all the values of $E$ in an interval $E\in[E_{1},E_{*}]$ where

[TABLE]

The critical value (14) appears in the study of close encounters with $P_{2}$ as a threshold value for considering a close encounters fast or slow (see, for example, Section 3.1 of [15]). With this classification we are considering the regime of slow close encounters. The function $u_{2_{*}}^{\phantom{3}}\!(E)$ is strictly monotone increasing, and the extremal values range from $u_{2_{*}}^{\phantom{3}}\!(E_{1})=-\sqrt{1-\mu-x_{L_{1}}}$ (so that $(u_{1},u_{2})=(0,u_{2_{*}}^{\phantom{3}}\!(E_{1}))$ corresponds to the equilibrium $L_{1}$) to $u_{2_{*}}^{\phantom{3}}\!(E_{*})=0$ (so that $(u_{1},u_{2})=(0,0)$ is the collision point). Except for $E=E_{1}$, we have ${\cal K}_{E}(u_{*},U_{*})\neq 0$, therefore $(u_{*},U_{*})$ does not correspond to a solution of the three–body problem (no equilibria of the CR3BP different from $L_{1}$ exist for $y=0$ and $x\in(-\mu,1-\mu)$). Nevertheless the equilibrium $(u_{*},U_{*})$ can be used to perform a Birkhoff normalization of the Levi-Civita Hamiltonian (10). We remark that for small values of $\mu$ the energy interval $(E_{1},E_{*})$ seems small (for instance, for $\mu=\mu_{J}$, $(E_{1},E_{*})\sim(-1.5193,-1.49809)$) but tiny variations of the energy in this interval produce large variations in the amplitude of the Lyapunov orbits.

In Appendix 1, paragraph A3, we show that for any $\mu\in(0,1/2)$, and any $E\in(E_{1},E_{*})$ the fictitious equilibria $(u_{*},U_{*})$ are of saddle-center type, with two imaginary eigenvalues $\pm i\omega$ and two real eigenvalues $\pm\lambda$, with $\lambda,\omega>0$. Moreover, for $a$ going to zero, both $\lambda,\omega$ tend to zero as $(u_{2_{*}}^{\phantom{3}}\!)^{2}$.

The quadratic part of ${\cal K}_{E}$ expanded at $(u_{*},U_{*})$. For all $E\in(E_{1},E_{*})$, the Levi-Civita Hamiltonian ${\cal K}_{E}$ has a saddle-center equilibrium in $(u_{*},U_{*})$ and, as indicated in Appendix 1 paragraph A4, one can construct a canonical change of variables

[TABLE]

conjugating ${\cal K}_{E}(u,U)$ to the Hamiltonian

[TABLE]

where $K_{0}={\cal K}_{E}(u_{*},U_{*})$,

[TABLE]

and $k_{j}(\xi,\eta,q,p)$ are homogeneous polynomials of degree $j$ expressed as sum of monomials of the type

[TABLE]

Note that the couple of conjugate variables $\xi,\eta$ is related to the hyperbolic behaviour, while the couple of conjugate variables $q,p$ is related to the elliptic behaviour (moreover, $q,p\in{\mathbb{C}}$ but, for real values of $u,U$, they satisfy $q=-ip$). The Hamiltonian in Eq. (16) is the starting point for the algorithm performing the Birkhoff normalization.

Normal forms of the Levi-Civita Hamiltonian at $L_{1}$. We reproduce with the Levi-Civita Hamiltonian the normalization methods which have been introduced in the Cartesian variables ([4], [31]). We perform a normal form scheme that uncouples (up to any arbitrary order) the hyperbolic variables $\xi,\eta$ from the elliptic variables $q,p$. The normalization is a near to the identity canonical transformation ${\cal C}_{N}$ which conjugates the Hamiltonian (16) to the normal form Hamiltonian

[TABLE]

where $N\geq 3$ is an integer called normalization order, ${\cal R}^{(N)}_{j}$ are polynomials of order $j$ expressed as sum of monomials of the type:

[TABLE]

while the $K_{j}$ are polynomials of order $j$ having a specific form defined by a certain strategy (see [35], [36] for an introduction to polynomial normal forms). A traditional strategy is to require that all the $K_{j}$ commute with $K_{2}$. For motivations specific of the three-body problem, such as the optimization of the computational time (which is crucial to perform a large number $N\geq 20$ of normalizations) and due to the bifurcations occurring in the spatial case preventing the integrability on the center manifold, the reduction is usually performed with a weaker normal form. In this paper, for the sake of comparison, we consider the three strategies:

• (a)

the $K_{j}$ depend on the variables $\xi$, $\eta$ only through the product $\xi\eta$ (i.e. the canonical transformation eliminates from the normal form every monomial of order smaller or equal than $N$ for which $k_{2}\neq l_{2}$, see [4]). As a consequence, by neglecting the terms ${\cal R}_{j}$ (which are of order $j\geq N+1$), from the normal form Hamiltonian ${\cal K}^{(N)}_{E}$, we obtain an approximate normal form of type

[TABLE]

where ${\cal J}$ is polynomial. From (21) one computes the center manifold (corresponding to $\xi,\eta=0$) as well as its stable (corresponding to $\xi=0,\eta\neq 0$) and unstable (corresponding to $\eta=0,\xi\neq 0$) manifolds.

• (b)

the $K_{j}$ contain only two type of monomials (18): monomials independent of $\xi,\eta$ and dependent on the $q,p$ only through the product $qp$, as well as monomials at least quadratic in $\xi,\eta$ (i.e. the canonical transformation eliminates every monomial for which $k_{2}+l_{2}=1$, or those which simultaneously satisfy that $k_{2}+l_{2}=0$ and $k_{1}\neq l_{1}$, see [31]). Neglecting the polynomials ${\cal R}_{j}$, we obtain an approximate normal form of type

[TABLE]

where ${\cal J}_{0},\ldots,{\cal J}_{3}$ are polynomial functions. From (22) one computes the center manifold (corresponding to $\xi,\eta=0$), but not its stable and unstable manifolds.

• (c)

the $K_{j}$ contain only two type of monomials (18): monomials independent of $\xi,\eta$ and monomials at least quadratic in $\xi,\eta$ (i.e. the canonical transformation eliminates every monomial for which $k_{2}+l_{2}=1$, to our knowledge this strategy has been never introduced before). Neglecting the polynomials ${\cal R}_{j}$, we obtain an approximate normal form of type

[TABLE]

where ${\cal J}_{0},\ldots,{\cal J}_{3}$ are polynomial functions. From (23) one computes the center manifold (corresponding to $\xi,\eta=0$), but not its stable and unstable manifolds.

While all these strategies allow one to compute the center manifolds, only strategy (a) allows one to compute also its stable and unstable manifolds. On the other hand, for any given normalization order $N$, the normal form (a) has less monomial terms than the normal form (b), which in turn has less monomial terms than the normal form (c). Thus, the computation of the normal form (a) requires the elimination of a larger number of monomials, and is heavier than the computation of the normal form (b) or (c). Therefore, if one is interested only in the computation of the Lyapunov orbits, then the best strategy is to adhere to the procedure that requires the minimum amount of eliminations, i.e. (b) or (c). On the other hand, if one is interested also in the computation of the stable and unstable manifolds $W^{s}_{1}(E),W^{u}_{1}(E)$, and the dynamics in their neighbourhood, then one adheres to the procedure (a).

Algorithm for the construction of the normal forms. For each $N\geq 3$, the canonical transformation ${\cal C}_{N}$ conjugating the Hamiltonian (16) to the normal form Hamiltonian (19) is obtained from the composition of canonical transformations:

[TABLE]

where ${\cal C}_{\chi_{N}}$ is the Hamiltonian flow at time $t=1$ of a suitable generating function defined from the coefficients of ${\cal K}^{(N-1)}_{E}$, and ${\cal C}_{2}$ is the identity. Below we describe the steps required for the computation of each ${\cal C}_{N}$ and ${\cal K}^{(N)}_{E}$ using the Lie series method (for an introduction to the method, see [37], [38]). We assume that the normal form ${\cal K}^{(N-1)}_{E}$ and ${\cal C}_{N-1}$ are known.

Step 1: From ${\cal K}^{(N-1)}_{E}$, we compute the generating function ${\chi_{N}}$, according to the chosen strategy:

– For strategy (a)

[TABLE]

– For strategy (b)

[TABLE]

– For strategy (c)

[TABLE]

Step 2: We compute the canonical transformation

[TABLE]

defined by the Hamiltonian flow at time $t=1$ of $\chi_{N}$ as the Lie series

[TABLE]

where $L_{\chi_{N}}\equiv\{\cdot,\chi_{N}\}$ is the Poisson bracket operator and $\zeta,\tilde{\zeta}$ denote any of the variables $\xi,q,\eta,p$ and $\tilde{\xi},\tilde{q},\tilde{\eta},\tilde{p}$ respectively.

Step 3: We compute the transformed Hamiltonian

[TABLE]

which, by construction, is in normal form up to order $N$, i.e. as in Eq. (19).

In summary, for any given $N$ we obtain the transformation yielding the original variables $(\xi,q,\eta,p)$ in terms of the normalized variables $(\xi^{(N)},q^{(N)},\eta^{(N)},p^{(N)})$, via

[TABLE]

The inverse transformation, i.e. the transformation yielding the normalized variables in terms of the original variables, is given by

[TABLE]

Solutions of Hamilton equations of the normal form Hamiltonian. Let us consider the solutions ($\xi^{(N)}(\tau)$,$q^{(N)}(\tau)$,

$\eta^{(N)}(\tau)$,$p^{(N)}(\tau)$) of the Hamilton’s equations of (19),

[TABLE]

[TABLE]

on the zero-energy level of ${\cal K}_{E}^{(N)}$, transformed back to the variables $u,U$ using the transformations (15) and (30),

[TABLE]

They project to solutions of the three-body problem of energy $E$ via the Levi-Civita transformation. In fact, since the transformations (15) and (30) are canonical, we only have to prove that, for any initial condition ($\xi^{(N)}(0)$,$q^{(N)}(0)$,$\eta^{(N)}(0)$,$p^{(N)}(0)$) satisfying

[TABLE]

transformed back to the variables $(u,U)$, we have

[TABLE]

This is a consequence of the definitions (16) and (19) of ${\cal K}^{(2)}_{E}$ and ${\cal K}^{(N)}_{E}$. In fact we have

[TABLE]

Computation of the Lyapunov orbits and their stable and unstable manifolds. Our aim is to use the normalizing transformations defined above to obtain analytic representations of the Lyapunov orbits as well as their local stable and unstable manifolds. To simplify the notation, from now on, we avoid the superscripts for normalized variables.

For any fixed value of $E>E_{1}$, we normalize the Hamiltonian ${\cal K}_{E}$ with a canonical transformation ${\cal C}_{N}$ (which also depends parametrically on $E$) and obtain the normal form Hamiltonian ${\cal K}^{(N)}_{E}$. The Lyapunov orbit of energy $E$ is represented in the normalized variables by

[TABLE]

and then it is mapped back to the original variables using the inverse of the canonical transformation ${\cal C}_{N}$ and all the transformations needed to pass from the Cartesian variables $(x,y,p_{x},p_{y})$ to the regularized variables $(u,U)$ and to the variables $(\xi,q,\eta,p)$.

The error in computing the set $LL_{1}(E)$ by formula (31), which neglects the remainder terms, depends on the norm of the remainder in a small neighbourhood of the set $LL_{1}(E)$, which is well represented by the quantity

[TABLE]

where the coefficient $\alpha_{N}$ depends on $N$ and $\rho$ denotes the maximum amplitude in the variables $q,p$ on $LL_{1}(E)$. In principle, one wishes to use the largest possible value of $N$ which is compatible with a reasonable CPU time and memory usage of modern computers. However, as it is typical of normalizing transformations, the coefficients $\alpha_{N}$ can increase with $N$ so that it may happen that $\alpha_{N}\rho^{N+1}$ decreases only up to a certain value of $N$, depending on $\rho$. The values of the coefficients $\alpha_{N}$ increase also as $E$ increases, since the values of $\lambda,\omega$, which appear at the denominators of the generating functions $\chi_{i}$, decrease. Moreover, since also $\rho$ increases as $E$ increases, we have an increment of the error terms $\alpha_{N}\rho^{N+1}$ which limits the validity of the method up to a certain value of $E$.

In this paper, we do not quantify a priori these errors, but we set up a numerical method to estimate their effects. In particular, rather than constructing the stable and unstable sets $W^{s}_{1}(E),W^{u}_{1}(E)$, we construct two surfaces, which we call the inner and outer stable/unstable tubes, containing the sets $W^{s}_{1}(E),W^{u}_{1}(E)$.

First, the stable and unstable manifolds $W^{s}_{1}(E),W^{u}_{1}(E)$ are computed from the normalization implemented using strategy (a). In this case, the sets are represented in the normalized variables by

[TABLE]

and then mapped back to the original Cartesian variables $(x,y,p_{x},p_{y})$. As a matter of fact, since the normalizing transformation is valid only in a small neighbourhood of the set $LL_{1}(E)$, we compute the following sections of the stable and unstable manifolds:

[TABLE]

with suitably small values of $\eta_{0}\neq 0,\xi_{0}\neq 0$. Then, the tubes $W^{s}_{1}(E),W^{u}_{1}(E)$ are constructed by computing numerically the orbits with initial conditions in a grid of points of the sets (33).

We quantify the effect of the errors in the computation of the sets $LL_{1},W^{s}_{1}(E),W^{u}_{1}(E)$ as follows:

$\quad(1)\,$ We consider a point in the set $LL_{1}$, obtained as explained above. Then, we compute numerically the Hamilton equations of the Levi-Civita Hamiltonian (10) with such initial condition and we represent it in the Cartesian variables $(x(t),y(t))$. In the case the initial condition is on the Lyapunov orbit, without any error, the numerical integration provides a periodic orbit of period $T$. Otherwise, the numerically integrated orbit does not exactly closes after $T$, and the distance of $(x(T),y(T))$ from $(x(0),y(0))$ provides an estimate of the error in the computation of $LL_{1}$.

$\quad(2)\,$ We consider for definiteness the stable manifold, and its branch which extends on the right of the Lyapunov orbit towards the singularity. The set $W^{s}_{1}(E)$ computed from the first of equations (33) is affected by some error. As done in [22], we can profit of a property of the tube manifolds to construct two surfaces, the inner and outer stable tubes $W^{s,in}_{1}(E)$, $W^{s,out}_{1}(E)$, that contain the true stable manifold $W^{s}_{1}(E)$. As soon as the two surfaces $W^{s,in}_{1}(E),W^{s,out}_{1}(E)$ are very close, the numerical errors in the computation of $W^{s}_{1}(E)$ is small.

The inner and outer stable tubes are defined as by computing the two cycles

[TABLE]

with $0<\xi_{0}<<\eta_{0}$. The parameters $\xi_{0},\eta_{0}$ are chosen so that by computing the initial condition on the two cycles (34) forward in time, we observe that all the orbits with initial conditions on one cycle approach the Lyapunov orbit and then bounce back on the right, while all the orbits with initial conditions on the other cycle approach the Lyapunov orbit and then transit to the left. Since the stable manifold is a separatrix for the motions which approach the Lyapunov orbit and then transit on the left or bounce back to the right, by computing the backward evolution of orbits with initial conditions in the cycles (34) we obtain two tubes containing the tube manifold $W^{s}_{1}(E)$.

In the next Sections we discuss the implementation of this theory and we present numerical demonstrations of the method to the computation of transit orbits.

3 Efficiency of the normal form computations

for the Sun-Jupiter system

Aiming to investigate the limits of the normalizations implemented with the Cartesian variables, we represent the family of Lyapunov orbits $LL_{1}(E)$ at different values of the energy, computed from the normalization of the Hamiltonian expressed with these variables with the three different normalizing strategies (a), (b), (c), see Fig. 2. The starting point of any Hamiltonian normalization in Cartesian variables includes the computation of the Taylor expansion $V_{N}$ at order $N$ of the singular potential energy

[TABLE]

appearing in Hamiltonian (3). As already remarked, for $\mu=\mu_{J}$, the function $V$ has complex singularities located at $x=x_{L_{1}},y=\pm i\sigma\sim 0.06\pm i0.066...$. For the values of the energy $E$ that we consider, we need to perform normalizations for $N\in[20,30]$. In this range of $N$ the complex singularity affects the Taylor expansions for $y>0.066$. For example, by denoting with

[TABLE]

an easy computation shows

[TABLE]

In order to check how this singularity affects the validity of the normal form expansions, in Fig. 2 we show the Lyapunov orbits corresponding to $E\in[-1.5193,-1.502]$ obtained from the normalization of the Cartesian Hamiltonian with the three strategies (a), (b), (c) for $N=30$. The orbits in Fig. 2 are obtained by numerically integrating an initial condition of the approximate Lyapunov orbit obtained from Eq. (6), without applying differential corrections333Since our aim here is to compare the limits of the different methods of normalization, we prefer not to reduce it with a differential correction of the initial conditions.. Due to the abundance of bibliography explaining this type of computations (see [4], [31] and references therein), we limit ourselves to present only the results. In the background we plot in color scale the value of the function $\log_{10}\left|\Delta V_{N}\right|$. With evidence, the largest reliable Lyapunov orbit obtained with the Cartesian variables strategy (a) (Fig. 2, left panel) is within the domain of the plane $x,y$ where the error function $\left|\Delta V_{N}(x,y)\right|$ is small. Instead, we observe that for the largest Lyapunov orbit obtained with the Cartesian variables-strategy (b) (Fig. 2, center panel), the amplitude is well beyond the limit of the complex singularity of $V$. As a matter of fact, the effect of the error due to the difference $\left|\Delta V_{N}\right|$ appears when we consider values of the hyperbolic variables $\xi,\eta\neq 0$. From the specific expression of the normal form $(b)$ (Eq. (22)), for any given $(q_{0},p_{0})$ we have:

[TABLE]

and therefore we expect a quadratic growth of $E-E_{0}(q_{0},p_{0})$ with $\xi$. This behavior is confirmed only for the Lyapunov orbits in black in Fig. 2, center panel. Precisely, by choosing a sample point for each Lyapunov orbit of Fig. 2, given by $q_{0},p_{0},\xi=0,\eta=0$ in the normalized variables, in Fig. 3 we show the values of the energy $E$ computed for $\eta=0$ and $\xi\in[0,1/10]$. Only for the Lyapunov orbits represented in black we have an agreement between the curves of Fig. 3 and the law (36). The same behaviour is found for the normalization method (c) (Fig. 2, right panel). We conclude that the limit of the normalizations implemented with the Cartesian variables is given by $\left|1-\mu-x_{L_{1}}\right|$.

We therefore proceed to check that the normalization of the Hamiltonian expressed with the Levi-Civita variables overcomes this limit. First, in Fig. 4 we show the Lyapunov orbits computed using method (a) in the left panel and method (c) in the right panel. The initial conditions used for each case have been obtained from the normalized Levi-Civita Hamiltonian ${\cal K}^{(N)}_{E}$ in Eq. (19) at the specific energy value, and they were numerically integrated with the non-normalized Hamiltonian ${\cal K}_{E}$ (10). Thanks to the implementation of the normalization method in the regularized Levi-Civita Hamiltonian, with both methods we compute Lyapunov orbits of libration amplitudes larger than the threshold value $0.066...$ of the Cartesian normalizations. Then, the threshold value of the largest Lyapunov orbit depends also on the normalization method. In fact, with the method (c) we reach a larger value for the amplitude with respect to the method (a). The limiting factor for strategy (a) is clearly the large number of terms to eliminate to render the Hamiltonian in normal form. Since for the computation of the stable and unstable manifolds we are forced to use method (a), we have an indication of the improvement introduced by the Levi-Civita variables, which allow us to reach a libration amplitude of about $1.4\ \left|1-\mu-x_{L_{1}}\right|$.

With the sake of explaining the differences between these techniques, we provide a few specific examples of the normal forms computations. We consider here the computation of Lyapunov orbits for two specific values of the energy, the first corresponding to a small amplitude Lyapunov orbit ($E=E_{a}:=-1.516$) and the second corresponding to a moderate amplitude Lyapunov orbit ($E=E_{b}:=-1.509$). For reference, $LL_{1}(E_{a})$ is represented with a dashed red curve in Fig. 4, while $LL_{1}(E_{b})$ corresponds to the dotted blue orbit. We expect that small libration orbits imply small errors in the computation. This is confirmed by the computation of the norm of the remainders $R^{N}$ (evaluated according to Eq. (32)) of such orbits, represented in Table 1, for $N\in\{3,\ldots,24\}$, for the normalizing strategy (c).

Additionally, in Tables 2 and 3, we compare the normal forms obtained with the three strategies for the first two steps of the normalization, giving ${K_{3}}$ and ${K_{4}}$ respectively, for $LL_{1}(E_{b})$. Each entry of the tables refers to a specific monomial appearing in the normal form, for which it is given the integers $k_{1}$, $l_{1}$, $k_{2}$, $l_{2}$, indicated in which normal form it is present (out of the three strategies), and the associated coefficient $\alpha$.

From Table 2, we can see that the less demanding normalization schemes (b) and (c) provide normal forms with many more terms than the strategy (a). At this step of normalization, there are no differences in the coefficients, since the three strategies are applied to the same initial expansion (the terms $K_{0}$ and $K_{2}$ are common to the three). When performing the following step, since now the schemes are applied also on $K_{3}$, the values of the coefficients start differing from strategy to strategy. The coefficients in Table 2 refer to the normalization by strategy (c), and we denote with different symbols the monomials whose coefficient has a different value for the strategy (a) ($\bullet$), the strategy (b) ($\dagger$) or both ($\star$). As we proceed with the normalizations, the resulting $K_{j}$ are different from each other. We also notice the increase in the amount of terms from step $3$ to $4$. This effect is evident also in Table 4, which summarizes the information related to the normal form computation (following strategy (c)) for the set of Lyapunov orbits appearing in Fig. 4 (right panel). The information provided corresponds to: the values of $E$, the order of the initial polynomial expansion $N_{p}$444As rule of thumb, we use an initial polynomial expansion at least three orders larger than the normalization order used, i.e. $N_{p}\geq N+3$., the normalization order $N$, the amount of terms of the normal form ${\cal K}_{E}^{(N)}$ and of the remainder ${\cal R}_{E}^{(N)}$ at order $N$, and the maximum distance that the orbit reaches with respect to the position of $L_{1}$.

As for increasing energies the periodic orbits have larger amplitudes, it is necessary to use larger normalizations orders (see Table 1). Notice, for instance, that a variation of $\sim 1$% in the energy (from $E=-1.518$ to $E=-1.505$) requires to use a normalization order 3 times larger, and the associated normal form is 35 times longer. This exponential growth makes hardly tractable the computation of these normal forms for larger values of the energy.

4 Numerical computation of the stable tubes manifolds

In Fig. 5, 6, 7 and 8, we represent the inner and outer stable tubes computed using the normalized Levi-Civita Hamiltonian as explained in Section 2 (see Eq. (34) and related discussion), for $E=-1.5193$, $-1.5183$, $-1.5167$, $-1.5113$ respectively. For each value of $E$, the inner/outer tubes are obtained first by considering $\xi_{0}=\pm 10^{-4}<<\eta_{0}=10^{-2}$, and by numerically computing a set of points $(q_{j},p_{j})$ ($j=0,\ldots,30$) sampling the level set

[TABLE]

Then, we compute the image of all the points $(\xi_{0},q_{j},\eta_{0},p_{j})$ in the Levi-Civita variables $u,U$. Finally, the tubes are obtained by numerically integrating the Hamilton’s equation of the Levi-CivitaHamiltonian ${\cal K}_{E}(u,U)$ for these initial conditions.

On the top panels of each figure we represent separately the orbits with initial conditions in the inner tube (blue line) and the orbits with initial conditions in the outer tube (green line). When the initial conditions on the outer stable tube are integrated forward in time, as it happens in the top panels, initially they approach the Lyapunov orbit (black curve), but then the orbits of the outer stable tube (orange orbits) bounce back, while the orbits of the inner stable tube (red orbits) transit to the left of the Lyapunov orbit. This mechanism allows us to identify one tube from the other. When integrated backward in time (bottom panels of Fig. 5, 6, 7, 8), the orbits of both inner and outer stable tubes almost overlap, indicating the location of the invariant stable manifold.

Figures 5 and 6 show the traditional behavior for the stable manifold, as described at (ii) in Section

1. When integrated backward in time, both tubes fold around Jupiter clockwise, following the sense of the collapsed invariant manifold emanating from $LL_{1}(E_{1})$ (blue line in the bottom panels), $W^{s}_{1}(E_{1})$. The smaller the value of the energy is, the better $W^{s}_{1}(E_{1})$ represent the behavior of the tube manifold, as it is evident from comparing the two figures. The limit of this representation takes place when the stable tube manifold allows for collisions with the planet. In Fig. 7 we present a stable manifold that, while still folding clockwise, passes very close to $P_{2}$. It is natural to expect that, by increasing the energy, the stable manifold is such that some orbits eventually collide with $P_{2}$. This is clearly visualized in Fig. 8 and, particularly, in 9. For the energy $E=-1.5113$, we show that the manifold embraces the position of the planet, therefore allowing possible collisions with it. In the right panel of Fig. 9 we show some of the orbits belonging to the manifold, indicating with arrows the sense of the orbit when integrated backward in time. We see that, besides passing extremely close to the singularity, the orbits can either revolve around the $P_{2}$ either clockwise or counter-clockwise, indicating that this value of the energy is, already, a good representation of orbits of type (iii) (see description in Introduction).

5 Conclusions

In this work, we studied the transition from non collisional to possibly collisional orbits as a outcome of the change in the structure of the tubes manifolds emanating from the Lyapunov orbits of the CR3BP for the Sun-Jupiter system. To this aim, we constructed a normal form of the Levi-Civita Hamiltonian, at a partially hyperbolic fictitious equilibrium. The normalization of the Levi-Civita Hamiltonian, which is regular at the position of the planet, allows us to reach values of the energy not investigated before via the traditional normalizations of the CR3BP Hamiltonian in Cartesian variables. We found that, for suitably large values of the energy, the tubes manifolds emanating from the horizontal Lyapunov orbits may contain orbits that collide with the secondary body before performing a full circulation around it (when integrated backwards from the Lyapunov orbit), or can revolve either clockwise or counter-clockwise around $P_{2}$. We established that, for $\mu=\mu_{J}$, the transition takes place for values of the energy larger than $E>-1.5167$. We expect this threshold to be relevant for future studies of dynamics of temporary trapped comets and/or for space mission design.

As future work, we plan to extend the present work to the spatial restricted three–body problem. The Levi-Civita regularization is extended to the spatial case by the Kustaanheimo-Stiefel regularization,

[TABLE]

see [39], [40]. Let us denote by ${\cal K}_{E}(u,U)$ the Hamiltonian representing the Kustaanheimo-Stiefel regularization (see for example [34], [41] for a detailed derivation). For values $E>E_{L_{1}}$, the fictitious equilibrium found in the regularized planar problem, $u_{*}=(0,u_{2_{*}}^{\phantom{3}}\!)$, $U_{*}=(-2u_{2_{*}}^{\phantom{3}}\!{}^{3},0)$, extends to an equilibrium of Hamilton equations of the Kustaanheimo-Stiefel Hamiltonian, $u_{*}=(0,u_{2_{*}}^{\phantom{3}}\!,0,0)$, $U_{*}=(-2u_{2_{*}}^{\phantom{3}}\!,0,0,0)$, satisfying the bi-linear equation $\ell(u_{*},U_{*})=0$ where $\ell(u,U)=u_{4}U_{1}-u_{3}U_{2}+u_{2}U_{3}-u_{1}U_{4}$. In principle, this be used to define a normalizing transformation as done in the present work. On the other hand, the internal symmetry of the Kustaanheimo-Stiefel transformation requires specific adaptations, which will be the subject of future work.

Acknowledgments

R.P. was supported by the ERC project 677793 StableChaoticPlanetM (01/07/2018-30/06/2019) and the project MIUR-PRIN 20178CJA2B (01/07/2019-30/06/2020). M.G. also acknowledges the project MIUR-PRIN 20178CJA2B ”New frontiers of Celestial Mechanics: theory and applications”.

Appendix

A1) Computing the equilibria $(u_{*},U_{*})$. The equilibrium points of ${\cal K}_{E}$, satisfying $u_{1}=0$ (which is the condition granting $y=2u_{1}u_{2}=0$, $x=x_{2}+u_{1}^{2}-u_{2}^{2}<x_{2}$) and $u_{2}\in(-1,0)$ ($u_{2}=-1$, $u_{2}=0$ locate the equilibrium on the primary body $P_{1}$ and $P_{2}$ respectively) are the solutions of

[TABLE]

We first notice that the first equation of this system provides $U_{1}=-2u_{2}^{3}$, while the second and the third ones are solved if $U_{2}=0$. Then, by using $U_{1}=-2u_{2}^{3}$ and noticing that since $u_{2}^{2}\in(0,1)$ we have $\sqrt{1-2u_{2}^{2}+u_{2}^{4}}=1-u_{2}^{2}$, the last equation is solved by $u_{2}$ satisfying $F(u_{2},E)=0$, where the function $F$ is defined in (13). A2) Solutions of the equation $F(u_{2},E)=0$. To solve equation $F(u_{2},E)=0$ we first re-write it in the form

[TABLE]

¿From standard calculus, we have

[TABLE]

for all $\mu\in[0,1/2]$ and $u_{2}\in(-1,0)$. As a consequence there is a strictly monotone increasing function $E:=E(u_{2})$ defined for $u_{2}\in(-1,0)$ such that $F(u_{2},E(u_{2}))=0$. Since we know that the Lyapunov orbits exist only for $E\geq E_{1}$, we restrict the domain of $E(u_{2})$ to the interval $[-\sqrt{x-x_{L_{1}}},0)$ and we define the inverse $a(u_{2})$ of $E(u_{2})$:

[TABLE]

[TABLE]

where $E_{*}$ has been defined in (14). A3) Linearization of the Hamilton equations of ${\cal K}_{E}$ at $u_{*},U_{*}$. The Jacobian matrix of the Hamilton vector field of ${\cal K}_{E}$ at $(u_{*},U_{*})$ is

[TABLE]

with

[TABLE]

As a consequence, the matrix $M$ has four eigenvalues:

[TABLE]

whose nature is determined by the values of the functions $A,B$. From standard calculus, for all $\mu\in(0,1/2)$ and $u_{2}\in(u_{2_{*}}^{\phantom{3}}\!,0)$ one proves: $\Gamma<0$, $B>0$, $B>\Gamma+5(u_{2_{*}}^{\phantom{3}}\!)^{4}$. As a consequence, the matrix $M$ has a couple of real eigenvalues $\pm\lambda$ and a couple of purely imaginary numbers $\pm i\omega$. Moreover, for $a$ going to zero, both $\lambda,\omega$ tend to zero as (a constant multiplying) $(u_{2_{*}}^{\phantom{3}}\!)^{2}$. A4) Diagonalization of the quadratic part of ${\cal K}_{E}$. Let ${v}_{i}$, $i=1,\ldots,4$, the four eigenvectors of $M$; matrix $A$ is constructed as :

[TABLE]

where $c_{1},c_{2}\in\mathbb{C}$ are chosen to satisfy the symplectic condition $A^{t}{\mathbb{J}}A={\mathbb{J}}$, where ${\mathbb{J}}$ denotes the standard symplectic matrix.