Periodic Orbits Around the Triangular Points with Prolate Primaries

The restricted problem of three bodies is considered one of the most important and famous problems of dynamics. This is due to its wide extremely important applications in the field of space dynamics, it describes accurately many real-world problems. In the restricted three-body problem, a body of negligible mass moves under the influence of the gravitational fields of two massive bodies. These two primary bodies rotate in circular or elliptic orbits about their common center of mass. Having negligible mass, the force exerted on the two primaries by the third body may be neglected.

The dynamic system of the restricted three-body problem is characterized by the presence of five equilibrium points. In this system, the gravitational and the centrifugal forces on a spacecraft mass cancel each other out. These fixed points are called equilibrium points. Three of these points are collinear, and two of them are triangular. These points rotate at the same frequency as the massive bodies, and thus the spacecraft mass’s position relative to the primaries is constant. This makes them very important for research and space operations (Marsola et al., 2021); (Reiff et al., 2022).

Furthermore, the periodic orbits around these equilibrium points acquired great attention and interest due to the crucial need for space orbits in the proximity of one of the collinear or triangular equilibrium points, (Abd El-Salam, 2019). Also, periodic orbits can be utilized to explore small solar system bodies, including comets and asteroids.

Different methodologies have been used to address the restricted three-body problem. In general, quantitative methods, either analytical or numerical, give precise and accurate information on the evolution of differential systems. However, this information is usually limited to the solution of interest and to a small vicinity. Also, in most cases, the accuracy decreases as time increases. In the current work, to obtain the required accuracy of the actual space mission orbit, we combined an analytical perturbed solution with a qualitative method. This technique gives partial but also rigorously demonstrated properties that are valid at least for long periods of time. Moreover, it deals with questions of existence, integrals of motion, uniqueness, periodic orbits, stability, etc.

Over the years, many researchers have investigated the issue of the restricted problem from various aspects, such as locations, stability of stationary points, and the periodic orbits, to mention some (Abouelmagd et al., 2016) (Burgos et al., 2019) (Pathak et al., 2019).

Recently, Poddar and Sharma, (2021) studied the equations of motion for the problem, which are regularized in the neighborhood of one of the finite masses. Further, the authors studied the existence of periodic orbits in a three-dimensional coordinate system when the reduced mass equals zero. Radwan and Abd El Motelp,(2021) investigated the linear stability of the restricted three-body problem when both of the massive primaries are triaxial. Also, they studied the periodic orbits in the vicinity of the triangular points. The authors showed that the shape of periodic orbits changed because of the triaxiality of the primary bodies. (Alrebdi et al., 2022) investigated how the mass ratio μ and the transition parameter influence the stationary points of the pseudo-Newtonian planar circular restricted problem. The authors also, showed how these parameters influence the networks of simple symmetric periodic orbits.

In the current work, we study the periodic orbits around the triangular points in the elliptic restricted three-body problem frame of work. To obtain a more realistic representation, the problem is generalized in the sense that bigger and smaller primaries are modeled as prolate spheroids. Also, we study in detail the variations in the angular frequencies for the long and short periodic orbits due to the shape of the primaries. Moreover, we compute explicit expressions for the eccentricities of the ellipses and determine the orientations of the principal axes for the ellipses that represent periodic orbits.

It is well known in the field of space science that most celestial bodies are often irregular in shape. In the original version of the restricted problem, the massive primaries are supposed to be spherical and symmetrical bodies (Szebehely, 1967). However, when studying various problems, the irregular shapes of these bodies must be taken into account in order to obtain highly efficient solutions. In some cases, considering the two primaries as point mass is not sufficient to describe the dynamic problem.

Over the past decades, several modifications have been proposed to include different additional parameters in the effective potential, such as the oblateness, the triaxiality, or the radiation of the two massive primaries (AbdulRaheem and Singh, 2008) (Beatty and Chaikin, 1999) (Radwan and Abd El Motelp, 2021) (Sharma and Subba, 1975) (Zahra et. al, 2017), and (Zotos, 2020). The mentioned reasons motivated us to study the dynamics of the problem under the influence of the real shape of the primaries. Furthermore, periodic orbits give more insights into a better understanding of the complex dynamical system of the restricted problem. Therefore, the crucial need for periodic orbits motivated us to study these orbits when both primaries are prolate spheroids.

The current dynamical system contains an infinitesimal mass that rotates in the orbital plane of the two massive bodies, the primary m₁ and the secondary mass m₂. The third infinitesimal one is considered to act as a test particle while the two primaries are prolate triaxial and circulate about their common centre of mass. The motion of the infinitesimal body doesn’t have any dynamic impact on the motion of the main bodies, due to its insignificant mass. In order to remove the time dependence from the equations of motion, it is better to use a synodic-rotating frame that rotates with constant angular velocity about the z-axis. The origin of the reference frame is centered at the barycentre of the system, and the x-axis lies on the line joining the two primary bodies. For convenience, we use a units system where the constant of gravity G and the distance between the centers of the two primaries are both equal to unity. Utilizing the reduced mass μ=m1m1+m2$$\mu = {{{m_1}} \over {{m_1} + {m_2}}}$$, we can express the dimensionless masses of the two primaries as m₁ = 1 – μ and m₂ = μ. Following the notations of Szebehely, (1967), the equations of motion of the tiny object in the dimensionless rotating-synodic frame are given by 1 x¨−2ny˙=∂U∂x,y¨+2nx˙=∂U∂y, $$\matrix{ {\mathop x\limits^{\unicode{x00A8}} - 2\,n\,\dot y = {{\partial U} \over {\partial x}},} & {\,\,\,\,\,\,\,\mathop y\limits^{\unicode{x00A8}} + 2\,n\,\dot x = } \cr } {{\partial U} \over {\partial y}},$$ where the amended potential function U can be written as 2 U=n22(x2+y2)+(1−μ)r1(1+Aσ2r12)+μr2(1+Aγ2r22) $$U = {{{n^2}} \over 2}({x^2} + {y^2}) + {{(1 - \mu )} \over {{r_1}}}(1 + {{{A_\sigma }} \over {2\,r_1^2}}) + {\mu \over {{r_2}}}(1 + {{{A_\gamma }} \over {2\,r_2^2}})$$ and 3 r1=(x+μ)2+y2,r2=(x+μ−1)2+y2, \begin{array}{*{35}{l}} {{r}_{1}}=\sqrt{{{\left( x+\mu \right)}^{2}}+{{y}^{2}},} \\ {{r}_{2}}=\sqrt{{{\left( x+\mu -1 \right)}^{2}}+{{y}^{2}},} \\ \end{array}\

The perturbed mean motion of the primaries is given by 4 n=1a(1+32(Aγ+Aσ)(1+e2)), $$n = \sqrt {{1 \over a}(1 + {3 \over 2}({A_\gamma } + {A_\sigma })(1 + {e^2}))} ,$$ where r₁ and r₂ are the distances of the two massive bodies from the infinitesimal third body. A_γ and A_σ represent the prolateness coefficients. a and e are the semi-major axis and eccentricity of either primary, respectively.

The locations of the equilibrium triangular points L₄ and L₅ can be obtained by setting all relative velocity and relative acceleration components equal to zero and solving the resulting system of equations U_x = U_y = 0. The first derivatives of the potential function can be written as 5 Ux=n2x−(3Aγ+2r22)μ(−1+x+μ)2r25−(−1+μ)(x+μ)[ 3Aσ2r15−1r13 ] \begin{array}{*{35}{l}} {{U}_{x}}={{n}^{2}}x-\frac{\left( 3\,{{A}_{\text{ }\!\!\gamma\!\!\text{ }}}+2\,r_{2}^{2} \right)\mu \left( -1+x+\mu \right)}{2\,r_{2}^{5}}-\left( -1+\mu \right) \\ \,\,\,\,\,\left( x+\mu \right)\left[ \frac{3\,{{A}_{\sigma }}}{2\,r_{1}^{5}}-\frac{1}{r_{1}^{3}} \right] \\ \end{array}\ 6 Uy=y[ n2+(−1+μ)(3Aσ2r15+1r13)−μ(3Aγ2r25+1r23) ] $${U_y} = y\left[ {{n^2} + \,\,( - 1 + \mu )\left( {{{3\,{A_\sigma }} \over {2\,r_1^5}} + {1 \over {r_1^3}}} \right) - \mu \left( {{{3\,{A_\gamma }} \over {2\,r_2^5}} + {1 \over {r_2^3}}} \right)} \right]$$

Since the perturbations considered in the present work are small, i.e., the prolateness coefficients are much smaller than unity, therefore, we can ignore its values (i.e., r₁ = r₂ = 1). Then it may be reasonable here to suppose that the locations of the triangular points L_4,5 are the same as given by classical restricted problem but perturbed by terms δ₁, δ2=OAγ,Aσ$${\delta _2} = {\cal O}\left( {{A_\gamma },\,{A_\sigma }} \right)$$). In this case, the solution of the classical restricted problem can be written as 7 ri=1+δi,δi<<1,(i=1,2). $$\matrix{ {{r_i} = 1 + {\delta _i},} & {{\delta _i} < < 1,} & {(i = 1,\,\,2)} \cr } .$$

Using equations (5) and (6) and solving for x and y up to order one in the involved small quantities δ₁, δ₂, we obtain 8 x=12(2δ1−2δ2−2μ+1),y=±321+43(δ1+δ2), $$\matrix{ {x = {1 \over 2}(2\,{\delta _1} - 2\,{\delta _2} - 2\,\mu + 1),} & {y = \pm } \cr } {{\sqrt 3 } \over 2}\sqrt {1 + {4 \over 3}\,\,({\delta _1} + {\delta _2})} ,$$

Substituting the values of x, y, r₁, and r₂ into equations (5) and (6), and expanding the resulting equations, we can retained only first order terms in δ₁, δ₂. Therefore, we get 9 δ1=13−13a(1+e2+32Aγ(1+e2)+Aσ(1+e2)),δ2=13−13a(1+e2+32Aσ(1+e2)+Aγ(1+e2)). \begin{array}{*{35}{l}} {{\delta }_{1}}=\frac{1}{3}-\frac{1}{3\,a}\left( 1+{{e}^{2}}+\frac{3}{2}{{A}_{\text{ }\!\!\gamma\!\!\text{ }}}\left( 1+{{e}^{2}} \right)+{{A}_{\sigma }}\left( 1+{{e}^{2}} \right) \right), \\ {{\delta }_{2}}=\frac{1}{3}-\frac{1}{3\,a}\left( 1+{{e}^{2}}+\frac{3}{2}{{A}_{\sigma }}\left( 1+{{e}^{2}} \right)+{{A}_{\text{ }\!\!\gamma\!\!\text{ }}}\left( 1+{{e}^{2}} \right) \right). \\ \end{array}\

Substituting the values of δ₁, δ₂ into equations (8) yields the coordinates of the equilibrium triangular points 10 x=12−μ−Aγ6a(1+e2)+Aσ6a(1+e2),y=±318[ 13−(1+e2)[ 4a+5a(Aγ+Aσ) ] ] \begin{array}{*{35}{l}} x=\frac{1}{2}-\mu -\frac{{{A}_{\text{ }\!\!\gamma\!\!\text{ }}}}{6\,a}\left( 1+{{e}^{2}} \right)+\frac{{{A}_{\sigma }}}{6\,a}\left( 1+{{e}^{2}} \right), \\ y=\pm \frac{\sqrt{3}}{18}\left[ 13-\left( 1+{{e}^{2}} \right)\left[ \frac{4}{a}+\frac{5}{a}\left( {{A}_{\text{ }\!\!\gamma\!\!\text{ }}}+{{A}_{\sigma }} \right) \right] \right] \\ \end{array}\

Note that if we ignore the involved perturbations, equations (10) will lead to the corresponding classical one.

It is well known that periodic orbits are of great importance, and they represent the backbone of studying the behavior of dynamic systems in the field of celestial mechanics. Let the locations of the equilibrium points be given as (x_L4,5, y_L4,5). Let us give the equilibrium points a small displacement (ξ₀, η₀), i.e., ξ₀, η₀ ≪ 1. We have 11 x=xL4,5+ξ0,y=yL4,5+η0 \begin{matrix} x={{x}_{{{L}_{4,5}}}}+{{\xi }_{0}}, & y={{y}_{{{L}_{4,5}}}}+{{\eta }_{0}} \\ \end{matrix}\

Then the corresponding characteristic equation of the current problem is given by Szebehely, (1967) 12 λ4+(4n2−UxxL4,5−UyyL4,5)λ2+UxxL4,5UyyL4,5−(UxyL4,5)2=0 $${\lambda ^4} + \,\,(4\,{n^2} - \,\,U_{xx}^{{L_{4,5}}} - \,\,U_{yy}^{{L_{4,5}}}){\lambda ^2} + \,\,U_{xx}^{{L_{4,5}}}\,U_{yy}^{{L_{4,5}}} - \,\,{(U_{xy}^{{L_{4,5}}})^2} = 0$$ where 13 UxxL4,5=−12+54a(1+e2)+Aσ[ −32+338a(1+e2)+3μ2−11μ4a(1+e2) ]+Aγ[ 118a(1+e2)−3μ2+−11μ4a(1+e2) ], [\begin{array}{*{35}{l}} U_{xx}^{_{{{L}_{4,5}}}}=\frac{-1}{2}+\frac{5}{4\,a}\left( 1+{{e}^{2}} \right)+{{A}_{\sigma }}\left[ \frac{-3}{2}+\frac{33}{8\,a}\left( 1+{{e}^{2}} \right)+\frac{3\,\mu }{2}-\frac{11\,\mu }{4\,a}\left( 1+{{e}^{2}} \right) \right] \\ \,\,\,\,\,\,\,\,\,\,+{{A}_{\text{ }\!\!\gamma\!\!\text{ }}}\left[ \frac{11}{8\,a}\left( 1+{{e}^{2}} \right)-\frac{3\,\mu }{2}+-\frac{11\,\mu }{4\,a}\left( 1+{{e}^{2}} \right) \right], \\ \end{array}\ 14 UyyL4,5=12+74a(1+e2)+Aσ[ −32+598a(1+e2)+3μ2−17μ4a(1+e2) ]+Aγ[ 258a(1+e2)−3μ2+17μ4a(1+e2) ], \begin{array}{*{35}{l}} U_{yy}^{{{L}_{4,5}}}=\frac{1}{2}+\frac{7}{4\,a}\left( 1+{{e}^{2}} \right)+{{A}_{\sigma }}\left[ \frac{-3}{2}+\frac{59}{8\,a}\left( 1+{{e}^{2}} \right)+\frac{3\,\mu }{2}-\frac{17\,\mu }{4\,a}\left( 1+{{e}^{2}} \right) \right] \\ \,\,\,\,\,\,\,\,\,\,+{{A}_{\text{ }\!\!\gamma\!\!\text{ }}}\left[ \frac{25}{8\,a}\left( 1+{{e}^{2}} \right)-\frac{3\,\mu }{2}+\frac{17\,\mu }{4\,a}\left( 1+{{e}^{2}} \right) \right], \\ \end{array}\ and 15 UxyL4,5=33 [−12+114a(1+e2)+(1−112a(1+e2))μ+Aγ( 298a(1+e2)+(5 −352a(1+e2)μ ) )+Aσ(−5+1118a(1+e2)+(5−352a(1+e2)μ)) ]. \begin{array}{*{35}{l}} U_{xy}^{{{L}_{4,5}}}=\frac{\sqrt{3}}{3}\left[ \frac{-1}{2}+\frac{11}{4\,a}\left( 1+{{e}^{2}} \right)+\left( 1-\frac{11}{2\,a}\left( 1+{{e}^{2}} \right) \right)\,\,\mu +{{A}_{\text{ }\!\!\gamma\!\!\text{ }}}\left( \frac{29}{8\,a}\left( 1+{{e}^{2}} \right)+\left( 5 \right. \right. \right. \\ \left. \left. \left. \,\,\,\,\,\,\,\,\,\,-\frac{35}{2\,a}\left( 1+{{e}^{2}} \right)\mu \right) \right)+{{A}_{\sigma }}\left( -5+\frac{111}{8\,a}\left( 1+{{e}^{2}} \right)+\left( 5-\frac{35}{2\,a}\left( 1+{{e}^{2}} \right)\,\,\mu \right) \right) \right]. \\ \end{array}\ UxxL4,5$$U_{xx}^{{L_{4,5}}}$$, UyyL4,5$$U_{yy}^{{L_{4,5}}}$$ and UxyL4,5$$U_{xy}^{{L_{4,5}}}$$ are the second partial derivatives of the amended potential function evaluated at the triangular points. The roots of the characteristic polynomial λ_i, of the present system, in the range 0 ≤ μ ≤ μ_critical, are purely imaginary. Therefore, the motion about the triangular equilibrium points L_4,5 is stable and composed of two harmonic motions governed by the variations 16 ξ=C1coss1t+D1sins1t+C2coss2t+D2sins2t,η=C¯1coss1t+D¯1sins1t+C¯2coss2t+D¯2sins2t \begin{array}{*{35}{l}} \xi ={{C}_{1}}\cos \,{{s}_{1}}t+{{D}_{1}}\sin \,{{s}_{1}}t+{{C}_{2}}\cos \,{{s}_{2}}t+{{D}_{2}}\sin \,{{s}_{2}}t, \\ \eta ={{{\bar{C}}}_{1}}\cos {{s}_{1}}t+{{{\bar{D}}}_{1}}\sin {{s}_{1}}t+{{{\bar{C}}}_{2}}\cos {{s}_{2}}t+{{{\bar{D}}}_{2}}\sin {{s}_{2}}t \\ \end{array}\ where s₁, s₂ are the frequencies for long and short periodic orbits, respectively. The coefficients C₁, D₁, C¯1$${\bar C_1}$$, and D¯1$${\bar D_1}$$ are the long periodic terms, while the coefficients C₂, D₂, C¯2$${\bar C_2}$$, and D¯2$${\bar D_2}$$ are the short periodic terms. The frequencies s₁, and s₂ are given up to order μ² as 17 s1=13226944 [ 9821a(1+e2)(−5488+47432μ+157542μ2)−7Aσ( 16464( −42 +2521a )−392(−1764+257321a)μ+732408621aμ2+21ae2( 411600 −1008616μ+7324086μ2 ) )+7Aγ( −7134421a+392(1764+833921a) μ+2197872621aμ2+21ae2(−71344+3268888μ+21978726μ2) )−6a21(−1+e2)( Aσ(44688+165592μ−13661634μ2)+14( −784+1624μ+128330μ2+Aγ(−2352−429352μ+2817990μ2)+14 [ 8( 686(421+21a)+14406ae2 +862421μ+16839921μ2−Aσ [ −8545621+633864a+56( 933521− 14406μa )+948114621μ2−57624ae2(−11+14μ) ]+8Aγ( 147( −821 +147a )−14(6421+7203a)μ+107747721μ2−7203ae2(−3+14μ) ) ]. \begin{array}{*{35}{l}} {{s}_{1}}=\frac{1}{3226944}\left[ 98\frac{\sqrt{21}}{a}\left( 1+{{e}^{2}} \right)\left( -5488+47432\,\mu +157542{{\mu }^{2}} \right)-7{{A}_{\sigma }}\left( 16464\left( -42 \right. \right. \right. \\ \left. \,\,\,\,\,+\frac{25\,\sqrt{21}}{a} \right)-392\left( -1764+2573\frac{\sqrt{21}}{a} \right)\,\,\mu +7324086\frac{\sqrt{21}}{a}{{\mu }^{2}}+\frac{\sqrt{21}}{a}{{e}^{2}}\left( 411600 \right. \\ \left. \left. \,\,\,\,\,\,-1008616\mu +7324086\,{{\mu }^{2}} \right) \right)+7{{A}_{\text{ }\!\!\gamma\!\!\text{ }}}\left( -71344\frac{\sqrt{21}}{a}+392\left( 1764+8339\frac{\sqrt{21}}{a} \right) \right. \\ \left. \,\,\,\,\,\,\mu +21978726\frac{\sqrt{21}}{a}{{\mu }^{2}}+\frac{\sqrt{21}}{a}{{e}^{2}}\left( -71344+3268888\,\mu +21978726\,{{\mu }^{2}} \right) \right)-6\,a\sqrt{21} \\ \,\,\,\,\,\left( -1+{{e}^{2}} \right)\left( {{A}_{\sigma }}\left( 44688+165592\,\mu -13661634\,{{\mu }^{2}} \right)+14\left( -784+1624\,\mu +128330\,{{\mu }^{2}} \right. \right. \\ \,\,\,\,\,+{{A}_{\text{ }\!\!\gamma\!\!\text{ }}}\left( -2352-429352\,\mu +2817990\,{{\mu }^{2}} \right)+14\left[ 8\left( 686\left( 4\sqrt{21}+\frac{21}{a} \right)+\frac{14406}{a}{{e}^{2}} \right. \right. \\ \,\,\,\,\,+8624\sqrt{21}\,\mu +168399\,\sqrt{21}\,{{\mu }^{2}}-{{A}_{\sigma }}\left[ -85456\,\sqrt{21}+\frac{633864}{a}+56\left( 9335\sqrt{21}- \right. \right. \\ \left. \left. \,\,\,\,\,14406\frac{\mu }{a} \right)+9481146\,\sqrt{21}\,{{\mu }^{2}}-\frac{57624}{a}{{e}^{2}}\left( -11+14\,\mu \right) \right]+8\,{{A}_{\text{ }\!\!\gamma\!\!\text{ }}}\left( 147\left( -8\,\sqrt{21} \right. \right. \\ \left. \left. \left. \,\,\,\,\,+\frac{147}{a} \right)-14\left( 64\,\sqrt{21}+\frac{7203}{a} \right)\,\,\mu +1077477\,\sqrt{21}\,{{\mu }^{2}}-\frac{7203}{a}{{e}^{2}}\left( -3+14\,\mu \right) \right) \right]. \\ \end{array}\ and 18 s2=13226944 [ −9821a(1+e2)(−5488+47432μ+157542μ2)+7Aσ( 16464 (42+2521a)−392(1764+257321a)μ+732408621aμ2+21ae2 ( 411600−1008616μ+7324086μ2 )−7Aγ( −7134421a+392( −1764+8339 21a )μ+2197872621aμ2+21ae2(−71344+3268888μ+21978726μ2) )+6a21(−1+e2)( Aσ(44688+165592μ−13661634μ2)+14( −784+1624μ +128330μ2 )+Aγ (−2352−429352μ+2817990μ2) )−14 [ 8( 686(421−21a) −14406ae2+862421μ+16839921μ2−Aσ [ −392(21821+1617a) +56( 933521+144061a )μ+948114621μ2+57624ae2(−11+14μ) ] +8Aγ( −147( 821+147a )−14(6421−72031a)μ+107747721 μ2+7203ae2(−3+14μ) ) ]. \begin{array}{*{35}{l}} {{s}_{2}}=\frac{1}{3226944}\left[ -98\frac{\sqrt{21}}{a}\left( 1+{{e}^{2}} \right)\left( -5488+47432\,\mu +157542{{\mu }^{2}} \right)+7{{A}_{\sigma }}\left( 16464 \right. \right. \\ \,\,\,\,\,\,\left( 42+\frac{25\,\sqrt{21}}{a} \right)-392\left( 1764+2573\frac{\sqrt{21}}{a} \right)\,\,\mu +7324086\frac{\sqrt{21}}{a}{{\mu }^{2}}+\frac{\sqrt{21}}{a}{{e}^{2}} \\ \left. \,\,\,\,\,\,\left( 411600 \right.-1008616\mu +7324086\,{{\mu }^{2}} \right)-7{{A}_{\text{ }\!\!\gamma\!\!\text{ }}}\left( -71344\frac{\sqrt{21}}{a}+392\left( -1764+8339 \right. \right. \\ \left. \,\,\,\,\,\,\left. \frac{\sqrt{21}}{a} \right)\mu +21978726\frac{\sqrt{21}}{a}{{\mu }^{2}}+\frac{\sqrt{21}}{a}{{e}^{2}}\left( -71344+3268888\,\mu +21978726\,{{\mu }^{2}} \right) \right) \\ \,\,\,\,\,+6\,a\,\sqrt{21}\,\,\left( -1+{{e}^{2}} \right)\left( {{A}_{\sigma }}\left( 44688+165592\,\mu -13661634\,{{\mu }^{2}} \right)+14\left( -784+1624\,\mu \right. \right. \\ \left. \,\,\,\,\,+128330\,{{\mu }^{2}} \right)+{{A}_{\text{ }\!\!\gamma\!\!\text{ }}}\left. \left( -2352-429352\,\mu +2817990\,{{\mu }^{2}} \right) \right)-14\left[ 8\left( 686\left( 4\,\sqrt{21}-\frac{21}{a} \right) \right. \right. \\ \,\,\,\,\,-\frac{14406}{a}{{e}^{2}}\,+8624\,\sqrt{21}\,\mu +168399\,\sqrt{21}\,{{\mu }^{2}}-{{A}_{\sigma }}\left[ -392\,\left( 218\,\sqrt{21}+\frac{1617}{a} \right) \right. \\ \left. \left. \,\,\,\,\,+56\,\,\left( 9335\,\sqrt{21}+ \right.14406\frac{1}{a} \right)\mu +9481146\,\sqrt{21}\,{{\mu }^{2}}+\frac{57624}{a}{{e}^{2}}\left( -11+14\,\mu \right) \right] \\ \left. \,\,\,\,\,+8\,{{A}_{\text{ }\!\!\gamma\!\!\text{ }}}\left( -147\left( 8\,\sqrt{21} \right. \right.+\frac{147}{a} \right)-14\left( 64\,\sqrt{21}-7203\frac{1}{a} \right)\,\,\mu +1077477\sqrt{21} \\ \left. \left. \,\,\,\,\,\,{{\mu }^{2}}+\frac{7203}{a}{{e}^{2}}\left( -3+14\,\mu \right) \right) \right]. \\ \end{array}\

It can be seen from equations (17) and (18) that the frequencies of the orbit of both short and long periodic motions are affected by the prolateness coefficients of the primaries, the mass ratio, the semi-major axis a, and the eccentricity e. as can be seen in the following illustrative graphs:

Figs. 1a and 1b illustrate the variations of the two frequencies s₂ and s₁ for different values of the prolateness coefficients A_σ, A_γ and e = 0.07, a = 0.94. Figs. 2a and 2b depict the variations of the short- and long-periodic frequencies s₂ and s₁ with the mass ratio μ for different values of the eccentricity of either primary (e = 0.05, 0.09, 0.4, A_σ = −0.06, A_γ = −0.04, and a = 0.94). It is observed that, in the above mentioned curves, the short-period frequency s₂ is a decreasing function, while the long-period frequency s₁ is an increasing one. Figs. 3a and 3b depict the variations of the long- ,and short periodic frequencies s₁ and s₂ with the mass parameter μ, for different values of the semi-major axis (a = 0.90, 0.95, 0.99). The figures show that the long period frequency, s₁, is an increasing function, while the short-period frequency, s₂ is a decreasing function. Figs. 4a, and 4b depict the variations in angular frequencies s₁, and s₂ under the effect of the perturbation considered in comparison with the classical case. It can be seen from both figures the effect of the perturbing forces on the behavior of the curves representing the angular frequencies. The perturbing forces cause these curves to depart from the classical case.

Figure 1.a.

Figure 1.b.

Figure 2.a.

Figure 2.b.

Figure 3.a.

Figure 3.b.

Figure 4.a.

Figure 4.b.

The expansion of the amended potential function U about the triangular equilibrium points L_4,5 is 19 U=UL4,5+UxxL4,5ξ2+UyyL4,5η2+UxyL4,5ξη+O(3) $$U = {U^{{L_{4,5}}}} + \,\,U_{xx}^{{L_{4,5}}}{\xi ^2} + \,\,U_{yy}^{{L_{4,5}}}{\eta ^2} + \,\,U_{xy}^{{L_{4,5}}}\xi \eta + {\cal O}(3)$$

As we can see equation (19) is quadratic, thus, the periodic orbits around the libration points L_4,5 are elliptical, since the Hessian UxxUyy−Uxy2>0$${U_{xx}}{U_{yy}} - U_{xy}^2 > 0$$.

6.2.

Eccentricities of the ellipses

In order to obtain the eccentricities of the ellipses, we use the equations Szebehely, (1967) 27 e1=(1−α12)12 \[{{e}_{1}}={{\left( 1-\alpha _{1}^{2} \right)}^{\frac{1}{2}}}\] 28 e2=(1−α22)12 \[{{e}_{2}}={{\left( 1-\alpha _{2}^{2} \right)}^{\frac{1}{2}}}\] and 29 αi=2sisi2+λ¯1=si2+λ¯22si \[{{\alpha }_{i}}=\frac{2{{s}_{i}}}{s_{i}^{2}+{{{\bar{\lambda }}}_{1}}}=\frac{s_{i}^{2}+{{{\bar{\lambda }}}_{\text{2}}}}{2\,{{s}_{i}}}\]] where λ¯1$${{\rm{\bar \lambda }}_1}$$ and λ¯2$${{\rm{\bar \lambda }}_2}$$ are the roots of the characteristic equation. For i = 1, 2, with similar expression for α₂, we have 30 α1=−13687936a350981+22977173 [ 392−3+21 (1+e2)( 28( 5838−127421+(−903+19721)Aγ+(151431−330521) Aσ )+( −4091010+89273821+(−17969553+3292131321)Aγ+( −8382927 +182931121 ) Aσ ) μ )+48−3+21( −2415+52721a2(−1+e2)(196−798Aσ−406μ−2957Aσμ+Aγ(42+7667μ))+212a52(−2+3e2)( 336(746413−488647)Aγ+( 107279400553− 70230856177 )Aγμ+(−1524368755913+997933630297)Aσμ−28( −14668643+9602887+676427923Aσμ−442826047Aσ− 775851233μ+507913897μ ) )+92a72(−2+5e2)( Aγ( 29997212 3−196377727+526998946533μ−345001798987μ )+28( 4895492 3−32048527+4227629213μ−2767633087μ )+Aσ( −9040057628 3+59181070207−105263356577733​μ+6891118085707μ ) )+112−3+21a( 56(−2415+52721)(7+22μ)+Aσ( −6938988+ 151421221+25797723μ−562955521μ )+Aγ( 405720−8853621− 2944578μ+64255421μ ) )+2942a32(−2+e2)( Aγ( 14490283− 9486127+2355337253μ−1541930287μ )+Aσ( −502969323 +329270767−23105682253μ+15126220567μ )+28( 707003− 462847+121(130973−85747)μ ) ) ]. \[\begin{array}{*{35}{l}} {{\alpha }_{1}}=\frac{-1}{\text{3687936}\,a}\sqrt{\text{350981}+\text{229771}\sqrt{\frac{7}{3}}}\left[ 392\sqrt{-3+\sqrt{21}} \right.\left( 1+{{e}^{2}} \right)\,\,\left( 28 \right. \\ \,\,\,\,\,\left( \text{5838}-\text{1274}\,\sqrt{21}+\left( -903+197\,\sqrt{21} \right){{A}_{\gamma }}+\left( \text{151431}-\text{3305}\,\,\sqrt{21} \right)\left. {{A}_{\sigma }} \right)+ \right. \\ \,\,\,\,\,\left( -\text{4091010}+\text{892738}\,\sqrt{21}+\left( -\text{17969553}+\text{32921313}\,\sqrt{21} \right){{A}_{\gamma }}+\left( -\text{8382927} \right. \right. \\ \left. \,\,\,\,\,+\text{1829311}\,\sqrt{21} \right)\left. {{A}_{\sigma }} \right)\left. \mu \right)+48\,\sqrt{-3+\sqrt{21}}\left( -\text{2415}+\text{527}\,\sqrt{21}\,{{a}^{2}}\left( -1+{{e}^{2}} \right) \right. \\ \,\,\,\,\,\left( 169-798{{A}_{\sigma }}-406\,\mu -2957{{A}_{\sigma }}\mu +{{A}_{\gamma }}\left( 42+\text{7667}\mu \right) \right)+21\,\sqrt{2}\,{{a}^{\frac{5}{2}}} \\ \,\,\,\,\,\left( -2+3\,{{e}^{2}} \right)\left( 336\left( \text{74641}\,\sqrt{3}-\text{48864}\,\sqrt{7} \right) \right.{{A}_{\gamma }}+\left( \text{10727940055}\,\sqrt{3}- \right. \\ \left. \,\,\,\,\,\text{7023085617}\sqrt{7} \right){{A}_{\gamma }}\,\mu +\left( -\text{152436875591}\,\sqrt{3}+\text{99793363029}\,\sqrt{7} \right){{A}_{\sigma }}\,\mu - \\ \,\,\,\,\,28\left( -\text{1466864}\,\sqrt{3}+\text{960288}\,\sqrt{7}+\text{67642792}\,\sqrt{3}\,{{A}_{\sigma }}\,\mu -\text{44282604}\,\sqrt{7}{{A}_{\sigma }}- \right. \\ \left. \left. \,\,\,\,\,\text{77585123}\sqrt{3}\,\mu +\text{50791389}\,\sqrt{7}\mu \right) \right)+9\sqrt{2}\,{{a}^{\frac{7}{2}}}\left( -2+5\,{{e}^{2}} \right)\left( {{A}_{\gamma }}\left( \text{29997212} \right. \right. \\ \left. \,\,\,\,\,\sqrt{3}-\text{19637772}\,\sqrt{7}+\text{52699894653}\,\sqrt{3}\,\mu -\text{34500179898}\,\sqrt{7}\mu \right)+28\left( \text{4895492} \right. \\ \left. \,\,\,\,\,\sqrt{3}-\text{3204852}\,\sqrt{7}+\text{422762921}\,\sqrt{3}\,\mu -\text{276763308}\,\sqrt{7}\,\mu \right)+{{A}_{\sigma }}\left( -\text{9040057628} \right. \\ \left. \left. \,\,\,\,\,\sqrt{3}+\text{5918107020}\,\sqrt{7}-\text{10526335657773}\,\sqrt{3}\,\mu +\text{689111808570}\,\sqrt{7}\,\mu \right) \right)+ \\ \,\,\,\,\,112\,\sqrt{-3+\sqrt{21}}\,a\left( 56\left( -\text{2415}+\text{527}\,\sqrt{21} \right)\left( 7+22\,\mu \right)+{{A}_{\sigma }}\left( -\text{6938988}+ \right. \right. \\ \left. \,\,\,\,\,\text{1514212}\,\sqrt{21}+\text{25797723}\,\mu -\text{5629555}\,\sqrt{21}\,\mu \right)+{{A}_{\gamma }}\left( \text{405720}-\text{88536}\,\sqrt{21}- \right. \\ \left. \left. \,\,\,\,\,\text{2944578}\mu +\text{642554}\,\sqrt{21}\,\mu \right) \right)+294\,\sqrt{2}\,{{a}^{\frac{3}{2}}}\left( -2+{{e}^{2}} \right)\left( {{A}_{\gamma }}\left( \text{1449028}\,\sqrt{3}- \right. \right. \\ \left. \,\,\,\,\,\text{948612}\,\sqrt{7}+\text{235533725}\,\sqrt{3}\,\mu -\text{154193028}\,\sqrt{7}\,\mu \right)+{{A}_{\sigma }}\left( -\text{50296932}\sqrt{3} \right. \\ \left. \,\,\,\,\,+\text{32927076}\,\sqrt{7}-\text{2310568225}\,\sqrt{3}\,\mu +\text{1512622056}\,\sqrt{7}\,\mu \right)+28\left( \text{70700}\,\sqrt{3}- \right. \\ \left. \left. \left. \,\,\,\,\,\text{46284}\,\sqrt{7}+121\left( \text{13097}\,\sqrt{3}-\text{8574}\,\sqrt{7} \right)\,\,\mu \right) \right) \right]. \\ \end{array}\]

Equation (26) determines the orientation of the orbits with respect to the rotational coordinate system. It is observed that the orientation of the orbits is affected by the involved perturbations. Equation (27) depicts, for i = 1, 2, the eccentricities of the short- and long-periodic orbits around the triangular points L_4,5.

We can observe that, from Fig. 5a, and 5b, the eccentricity of the long-period orbit decreases under the effect of the perturbations, while the eccentricity of the short-period one increases. Also, we see from the figures that the perturbed case are shifted from the classical case because of the influence of the disturbing forces. Ignoring all the perturbations considered in the present work, our results will be the same as those obtained by Szebehely, (1967).

Figure. 5.a.

Figure 5.b.

In this work, we have investigated the periodic orbits around the triangular libration points L_4,5, in the range 0 < μ < μ_c. We formulated the problem in a more general way and used a more complex mathematical model than previously published papers that considered the classical case (Abouelmagd and Mostafa, 2015). The prolateness coefficients of both primaries are taken into account as a perturbing force. We investigated the variations of the angular frequencies for the long and the short periodic orbits. The variation of both frequencies is represented graphically versus the mass parameter μ for distinct values of the included perturbations. It is found that for small mass ratio μ, an increment in the perturbing forces results in a decrease in the frequency of the short-period orbit, while an increment in the same parameters will increase the frequency of the long-period one. Both frequencies coincide at the critical value of the mass parameter μ_c. In addition, we derived explicit expressions for the eccentricities e₁ and e₂ of the long and short-period orbits. We represented graphically both eccentricities versus the mass parameter. It is found that the eccentricities e₁ and e₂ of the long and short-period orbits are decreasing and increasing functions, respectively. Furthermore, we studied the orientation of the principle axes of the ellipses. It is observed that the included perturbing forces influence the orientation of the principal axes. The perturbing forces result in a change in the inclination angle of the orbits. Finally, in our opinion, we believe that the current research has special importance to space science applications to send spacecraft into stable regions in planetary systems.