Abstract

We illustrate the chaotic nature of the circular restricted three-body problem from the perspective of the bifurcation diagram with respect to the mass ratio parameter. Moreover, it is shown that when the frequency ratio in different directions of the classical problem is irrational, it has the quasiperiodic characteristics. In addition, a three-dimensional approximate solution to this problem under two time scales is proposed by using the multiple time scales method.

1. Introduction

As early as the 19th century, plenty of mathematicians, such as Dirichlet and Weierstrass, et al. [1], expected to obtain a series solution of the three-body problem in the following form:

The rate of change of the solution with respect to time in equation (1) appears as a combination of many incommensurable frequencies, , often called quasiperiodic solution.

Almost half a century ago, Farquhar and Kamel [2] proposed an approximation to construct the periodic orbits described above. A few years later, Richardson [3, 4] developed third-order analytic solutions to collinear libration points (CLP) for a class of circular restricted three-body problem (CRTBP) based on a method similar to Lindstedt–Poincaré and successive approximation. The constructed solutions are the basis of determining halo orbits around these points, and they can provide approximate initial values of halo orbits. Furthermore, the solutions can also be developed to investigate spacecraft formation flying. The more accurate the approximate analytical solutions, the less fuel in the system will be consumed. Ruijgrok [5] solved a planar three-body problem for equal mass particles under a particular class of three-body forces. Zagouras and Markellos [6] studied a periodic spatial solution to Hill’s problem, which is the limiting case of a restricted three-body problem (R3BP). Papadakis et al. [7] considered the periodic orbits of the three-body problem generated by CLP and provided the second-order approximate analytic expressions of such orbits. Gómez et al. [8] computed quasihalo orbits in CRTBP semianalytically through the Lindstedt–Poincaré technique. Ten years ago, Lu and Zhao [9] proposed an improved and more precise third-order analytic solution than Richardson’s classical analytic solution. Gidea and Deppe [10] numerically studied the influence of small perturbations on the dynamics of an infinitesimal third body. They investigated the chaotic orbits in an RTBP, which also indicated the complexity and difficulty of obtaining approximate analytic solutions to the problem.

In recent years, Šuvakov and Dmitrašinovic [11] presented some creative results for the periodic orbits of the three-body problem. Thirteen new, zero-angular-momentum, distinct equal mass, planar collisionless periodic orbits are displayed in three new categories. When the two primaries are both oblate spheres, Mittal et al. [12] used the predictor-corrector algorithm to construct periodic orbits. They showed the corresponding periodic orbits under different energy constants, mass ratios, and oblateness factors of the two primaries. Based on photogravitational planar RTBP with oblateness, Pathak et al. [13] studied the seventh-, ninth-, and eleventh-order internal resonance periodic orbits of the Sun-Earth system by using Runge–Kutta–Gill method. For a generalized photogravitational RTBP model, the two primaries are oblate spheroid and are under the gravity of an asteroid belt. Abouelmagd et al. [14] found a secular solution that can be reduced to a periodic one near the triangular libration points when the mass ratio is equal to the critical mass value. Similar results can also be found in [15]. Selim et al. [16] analyzed the stable motion solutions of long-period orbits and short-period orbits when the primary is a triaxial rigid body and the Euler angle of rotation satisfies certain conditions. Besides, for two types of the perturbed RTBP, Gao et al. [17, 18] utilized the Lindstedt–Poincaré perturbation method to give an approximate analytical solution to the periodic orbits near the CLP.

Considering that the dynamical equation of the third body in CR3BP is a time-varying high-dimensional nonlinear system in the inertial frame, it can hardly be solved by using analytical approaches. However, this case will be better in the rotating coordinate system, for the aforementioned governing equations are time-independent, and there is an integrable motion. Using the method of multiple time scales, Nayfeh [19, 20] studied the 3 : 1 and 2 : 1 small-amplitude resonances near the triangular libration points in the plane when the potential function was expanded to the second- and third-order terms of a small parameter, respectively.

Based on the number of terms expanded by the potential function of the problem, the existing time-scale solutions mainly include the three-dimensional single time-scale solution when the potential function expands to the third-order, the three-dimensional (3D) double time-scale solution when the potential function expands to the second-order, and the planar approximate analytic solution. Accordingly, when the third body is also considered to vibrate with a small amplitude in space, the multiple time scales method [21] will be used to construct the approximate 3D CR3BP solution in this paper. The 3D multiple time scales solution has the following form:where , , , and is a small, dimensionless parameter.

2. Dynamic Equations of Classical CRTBP

In this section, the dynamic equations of classical CRTBP in three-dimensional space will be described.

In a rotating framework, it is well known that the governing equations of the classical CRTBP can be characterized by the following set of differential equations:where the potential function

For more details on the establishment of equation (3), please refer to [22].

It is easy to find that in equation (4) is the only parameter in the classical CRTBP. In this paper, we define it as the ratio of the mass of one primary to the mass sum of the two primaries, so it is called the mass ratio parameter, and its value range is between 0 and 1. Under the initial conditions [−0.001, 0.2, 0.1, −0.3, 0.87, 0.01], we divide the range of into 1000 equal parts, and then, respectively, return its corresponding coordinates of three different planes (x-y, x-z, y-z) for each value of , and the bifurcation diagram with respect to the parameter based on ode45 algorithm shows that any small change of will lead to the third body’s apparent chaotic dynamic behavior (see Figure 1). We find that any value of corresponds to an infinite number of points rather than a single point, which means that the probability of obtaining periodic solutions of the three-body problem by numerical simulation is almost zero. Moreover, even if there are some results that are seem to be periodic, part of the existing conclusions based on numerical simulation will no longer hold when the step size continues to shrink. When the iteration accuracy is gradually improved, it is difficult for the orbit of the three-body system with chaotic characteristics to come back after iteration from a particular initial point, and researchers do not even know what will happen when the accuracy is less than 10⁻¹⁶. Therefore, we expect that the periodic solution which is beneficial to the actual mission can be realized by “construction.”

(a)

(b)

(c)

(a)
(b)
(c)

Figure 1 

(a–c) Bifurcation diagrams of mass ratio parameters in µ-x-y, µ-x-z, and µ-y-z frames, respectively.

Next, we first deal with the complex potential function by introducing Legendre polynomials , which satisfies Rodrigues formula:

Then, we have [22]where and .

Substituting equations (4) and (6) into equation (3), it yieldswhere , for libration points and and , for libration point , , where represents the angular velocity of relative motion between the two primaries.

Then, equation (7), up to the third approximation, can be written in the following form:

3. Construction of Approximate 3D Multiple Time Scales Solution

Consider the small amplitude of nonhomogeneous terms of equation (8). Then,

We write the solution of equation (9) as follows:where .

Then, time derivatives becomewhere .

Substituting equations (10)–(12) into equation (9), we find

Now, from equations (13a)–(13c), equating the coefficients of , , and to zero, it leads to the following three sets of equations.

Order :

Order :

Order :

Consider that the last equation of (14) is decoupled from the first two equations of (14). We employ the following transformations:

Then, the first two equations of (14) can be represented by the following equivalent form:

Thus, the coefficient matrix of equation (18) can be denoted as

Therefore, the characteristic equation can be obtained aswhere

It is noticed that and are two equal and opposite real roots, and and are a pair of pure imaginary roots. Hence, the general solution of equation (14) assumes the following form:where coefficients and are all constants and is the mould of and .

Since the motion of the third body may be unbounded due to the influence of exponential terms, appropriate initial conditions are selected such that , which implies thati.e.,where

If the value of in equation (24) is an irrational number, according to the time history diagrams in Figure 2, this linearization motion of the system appears to make the periodic motion in the and components, respectively. However, it is not periodic but Lissajous-type quasiperiodic in the x-y plane, and this fact can be verified by the integral curves with respect to and (see Figure 3) and its projections in the x-y plane, as well as the portraits in the frame of x-y-z (see Figure 4). The significance of Figures 2–4 is that if the frequency of the third body in all directions is not commensurable, it will exhibit a quasiperiodic motion. Therefore, the value of is a rational number that performs a critical role in the motion law of the system. In general, if the value of is an irrational number, the time history diagram can only serve as a reference but a standard to judge if the system makes the periodic motion or not. On the contrary, the motion law of the system can be judged by phase portraits and time history diagrams.

(a)

(b)

(a)
(b)

Figure 2 

Time history diagrams with respect to (a) x and (b) y, respectively, when .

Figure 3 

Quasiperiodic orbit of Lissajous type in x-y-t space for .

(a)

(b)

(a)
(b)

Figure 4 

Quasiperiodic orbit of Lissajous type in the (a) x-y plane and (b) x-y-z space, respectively, when .

In addition, note that the halo-type periodic orbits play particular importance in the actual mission ACE, ISEE-3/ICE, MAP, SOHO, Genesis, etc., that is, the in-plane and out-of-plane motion makes those characteristic frequencies equal within the sufficiently large motion region. Therefore, without loss of generality, we consider here the case of . Then, equation (24) can be represented as

Substituting equation (26) into equation (15), it yields

Let the secular term be equal to zero, and we find that

Hence,where and are constants.

Substituting equation (29) into equation (26), and set constant , we obtain

Then, equation (27) can be expressed as

Similarly, solutions to equation (31) can be denoted as

Substituting equations (26) and (32) into equation (16), it yields

We define

Thus,