Impacts of Poynting–Robertson Drag and Dynamical Flattening Parameters on Motion around the Triangular Equilibrium Points of the Photogravitational ER3BP

Using an analytical and numerical study, this paper investigates the equilibrium state of the triangular equilibrium points 
 
 
 
 L
 
 
 4
 ,
  
  
 5
 
 
 
 of the Sun-Earth system in the frame of the elliptic restricted problem of three bodies subject to the radial component of Poynting–Robertson (P–R) drag and radiation pressure factor of the bigger primary as well as dynamical flattening parameters of both primary bodies (i.e., Sun and Earth). The equations of motion are presented in a dimensionless-pulsating coordinate system 
 
 
 
 ξ
 −
 η
 
 
 
 , and the positions of the triangular equilibrium points are found to depend on the mass ratio 
 
 
 
 μ
 
 
 
 and the perturbing forces involved in the equations of motion. A numerical analysis of the positions and stability of the triangular equilibrium points of the Sun-Earth system shows that the perturbing forces have no significant effect on the positions of the triangular equilibrium points and their stability. Hence, this research work concludes that the motion of an infinitesimal mass near the triangular equilibrium points of the Sun-Earth system remains linearly stable in the presence of the perturbing forces.


Introduction
e study of the equilibrium state of an infinitesimal mass (a test particle) with regard to the dynamical system of the restricted three-body problem (R3BP) remains one of the most important and interesting aspects in the study of celestial mechanics and dynamical astronomy. e elliptic restricted three-body problem (ER3BP) deals with the description-study of the motion of an infinitesimal mass (m 3 ), in the vicinity of the gravitational fields of two dominant bodies m 1 and m 2 which are called the primaries (bigger primary and smaller primary, respectively), where m 1 ≥ m 2 ≫ m 3 . ese primary bodies revolve about their common centre of mass in elliptic orbits, under the influence of their mutual gravitational attraction. e ER3BP admits five equilibrium points at which the test particle (infinitesimal mass) would remain fixed if placed there. ree of such points lie on the line joining the two dominant bodies m 1 and m 2 , called the collinear equilibrium points L i (i � 1, 2, 3), and are unstable, while the other two points form equilateral triangles with the two dominant bodies, called the triangular equilibrium points L 4 and L 5 , and are stable for 0 < μ < μ c , where μ is the mass ratio defined by μ � (m 2 /(m 1 + m 2 )) and μ c is the critical mass parameter [1].
In an attempt to have a much more realistic description of the motion of an infinitesimal mass over the decades, the classical R3BP has been modified in the sense that additional dynamical potentials of the system were considered in different approaches [2][3][4][5][6][7][8] and others.
Abouelmagd [9] in his study of the R3BP found out that the locations of the triangular points and their linear stability are affected by the oblateness of the more massive primary, up to the second zonal harmonic J 4 . In his numerical study, he concluded that the existence of J 4 sometimes does not affect the stability of the equilibrium points as in the Earth-Moon, Saturn-Phoebe, and Uranus-Caliban systems. Also, Abouelmagd et al. [10] examined the effects of oblateness of the three participating bodies together with small perturbations in Coriolis and centrifugal forces. Afterwards, Abouelmagd et al. [11] studied the effect of the first two even zonal harmonic coefficients of both the primaries as well as the periodic orbits around the equilibrium points. ey found out that the triangular equilibrium points are linearly stable for 0 < μ < μ C and unstable for μ C ≤ μ ≤ 0.5, while the collinear equilibrium points remain linearly unstable.
Singh and Tyokyaa [12] examined the stability of the triangular points in the ER3BP with oblateness up to the second zonal harmonic J 4 of both primaries. ey concluded that the location and linear stability of the triangular points are affected by the oblateness of the primary bodies, eccentricity of the orbits of the primaries, and the semi-major axis of the system and both destabilized the system. In the same year, Singh et al. [13] investigated the influence of the zonal harmonics (J 2 and J 4 ) of the primary and the radiation pressure of the secondary on the locations and linear stability of the triangular points. ey however claimed that the parameters involved in the system affect the positions of the equilibrium points and destabilize the system as well.
e Doppler shift and absorptions and subsequent reemission of incident radiation, that is, the so-called Poynting-Robertson's (P-R) drag, are often neglected by many researchers in the estimation of light radiation force. Poynting [14] while studying radiation in the Solar System stated that radiation affects temperature and small bodies. He asserted that particles such as cosmic dust grain or small meteors are affected considerably by gravitational and light radiation force as they approach luminous celestial bodies. Furthermore, infinitesimal bodies in solar orbits suffer a gradual loss of angular momentum and ultimately spiral into the Sun. Later, Robertson [15] in a modified theory of Poynting considered only terms of the first order in the ratio of velocity of the particles to that of light. He investigated the dynamical effects of drag in the Solar System and derived the expression for the times of fall from circular orbits. us at a cosmically rapid rate, the P-R effect sweeps small particles of the Solar System into the sun.
Researchers like Burns et al. [16], Murray [17], Singh and Simeon [18], Alhussain [19], Chakraborty and Narayan [20], Amuda et al. [21], and others studied the R3BP by taking into account the P-R drag in different views. Mishra et al. [22] examined the stability of triangular points under the assumption that the bigger primary is a source of radiation with the incident P-R drag while the smaller is an oblate spheroid in the frame of the ER3BP. ey concluded that the triangular points are unstable. In a recent study, Singh and Amuda [23] investigated the linear stability around L 4, 5 of a test particle in the field of post-AGB binary system with the effective P-R drag force. ey asserted that P-R drag and the mass parameter μ contribute in shifting the locations of the triangular points and the triangular points are unstable in the linear sense due to the presence of complex conjugate roots. e aim of the present paper is to further investigate the effects of radiation pressure and P-R drag of the bigger primary on L 4, 5 in the ER3BP, taking into account the effects of dynamical flattening parameters of the primary bodies.

Equations of Motion
Let m 1 , m 2 , and m 3 be the three masses, m 1 and m 2 being the dominant bodies having an elliptic orbit about their common centre of mass, while m 3 being the infinitesimal mass which moves in the same plane with the dominant bodies under the influence of their force-field without influencing their motion. Let (x 1 , y 1 , z 1 ), (x 2 , y 2 , z 2 ), and (x 3 , y 3 , z 3 ) denote the coordinates of m 1 , m 2 , and m 3 , respectively, in the sidereal coordinates. Using Newton's law, the equations of motion of an infinitesimal mass m 3 in the sidereal coordinates system can be represented as where V � m 3 k 2 ((m 1 /R 1 ) + (m 2 /R 2 )); R 1 and R 2 are the distances of an infinitesimal mass from m 1 and m 2 , respectively, and are defined as ; k 2 is the Gaussian constant of gravitation and the dot indicates differentiation with respect to time t.
Equations (1) can be rewritten as [3] € where W � (V/m 3 ). Now, we choose and rotate the synodic coordinates uniformly with a positive unit rate, which has the same origin at the centre of mass of the two dominant bodies with the sidereal coordinates. e direction of the x− axis is chosen such that the two dominant bodies always lie on it. erefore, the equations of motion in sidereal coordinates are related to the equations of motion in synodic coordinates with respect to the true anomaly. By transforming the equations of motion from true anomaly to eccentric anomaly in a dimensionlesspulsating (rotating) coordinate system (ξ, η, ζ), we have where the prime ( ′ ) represents differentiation with respect to eccentric anomaly (E) and Ω(ξ, η, ζ) � (1/(n 2 ρ))[(n 2 /2) , e and n are the eccentricity of the orbits and mean motion of the dominant bodies, while are the distances of an infinitesimal mass from the bigger primary and smaller primary, respectively, and are defined as , and a is the semi-major axis of m 2 around m 1 . See [3,24,25]. Now, integrating equations (3) with respect to the eccentric anomaly (E) and averaging, we get where Ω is the potential-like function defined by us, the defined dynamical system in equations (4) is the required equations of motion for the ER3BP.
Let q be the radiation pressure factor of the bigger primary which is given by [26], where F g1 is the gravitational force of m 1 and F p1 is the radiation pressure of m 1 .
Considering the potential theory, the external gravitational potential due to a body that has axial symmetry can be written in terms of Legendre polynomials as V 0 � − ((k 2 m 0 )/r 0 )[1 − ∞ n�2 J n P n (cos θ)(R 0 /r 0 ) n ] (see de Pater and Lissauer [27] for more details), where m 0 denotes the mass of the body; r 0 denotes the radial distance from the centre of the particle to the centre of any other body; θ denotes the angle between the body's symmetry axis and vector to the particle; R 0 denotes the mean radius of the body; J n denotes the dimensionless coefficient that characterizes the degree of nonspherical components of the potential, J n is zero for odd n, and when n is even, J n is called a zonal harmonic coefficient; and the term P n (cos θ) denotes the Legendre polynomials of degree n and is defined by P n (x) � (1/(2 n n!))(d n /dx n )(x 2 − 1) n . We assume that the infinitesimal mass moves in the same plane of motion as the dominant bodies and this plane coincides with the equatorial plane (i.e. θ � 90 0 ). erefore, the potential energy of the infinitesimal mass under the effects of the dynamical flattening parameters of both dominant bodies can be written as where 2 (i � 1, 2) represent the dynamical flattening parameters of the bigger and smaller primaries, respectively, J A 2i and J B 2i (i � 1, 2) are the zonal harmonics coefficients; r i (i � 1, 2) is defined as in equation (3); and R 1 and R 2 are the mean radii of the dominant masses m 1 and m 2 , respectively. e distance between the dominant bodies is r � a(1 − e cos E) in the elliptic orbit, where a, e , and E are semi-major axis between the dominant bodies, common eccentricity of the dominant bodies, and eccentric anomaly, respectively.
And the mean distance between them is as follows: Assuming that the dominant bodies are in elliptical motion with constant angular velocity n(mean motion), then the orbits of m 1 and m 2 with respect to the centre of mass, with semi-major axes, would be a 1 � (m 2 /(m 1 + m 2 ))a and a 2 � (m 1 /(m 1 + m 2 ))a, respectively, having the same eccentricity; thus, the motion of the bigger and smaller primary can be written as [1,28] where n, r, and k are the mean motion, distance between the dominant bodies, and Gaussian constant of gravitation, respectively. e distance between the dominant bodies r is defined to be the semi-major axis a of the orbit (i.e., r � a), since the dominant bodies are in elliptic orbits. Hence, adding equations (8) together, we obtain Assume that the sum of the masses of the dominant bodies and the semi-major axis between them are the units of mass and length, respectively; i.e., m 1 + m 2 � 1, and a 1 + a 2 � a � 1. Also, the unit of time is chosen so as to make the Gaussian constant, k 2 � 1. Hence, equation (9) becomes Using equations (4), (5), (6), and (10), the equations of motion of an infinitesimal mass in the frame of the ER3BP can be modified, taking into account the dynamical flattening parameters of both dominant bodies together with the radiation pressure as well as P-R drag due to the bigger primary in a dimensionless-pulsating (rotating) coordinate system (ξ, η) as where Advances in Astronomy and c d is the dimensionless speed of light [18,22,29].

Locations of the Triangular Equilibrium Points L 4, 5
To obtain the equilibrium positions of an infinitesimal mass, the equations ξ ′ � η ′ � ξ ″ � η ″ � 0 must be satisfied in the equations of motion (11); i.e., they are the solutions of the equations U ξ � U η � 0, and thus resulting in which can be rewritten as

Stability of the Triangular Equilibrium Points
To examine the stability of the triangular equilibrium points, we place the infinitesimal mass at one of the equilibrium points and give it a small velocity. e point is stable for the oscillatory solutions with small amplitude and unstable for exponentially diverging solutions.
Assume that (ξ 0 , η 0 ) are the coordinates of the equilibrium points under consideration and let (x, y) be the small displacements from these coordinates of the equilibrium points. ese can be written as Advances in Astronomy Now, using equations (25) and holding only the linear terms of Taylor's theorem as the second and higher power of x and y being very small, equation (11) become where the superscript 0 of equations (26) indicates that the partial derivatives are evaluated at the equilibrium points Suppose x � Ae λt and y � Be λt are the trial solution of equations (27) (variational equations). en, by using these values of the trial solutions in equations (27), we get

Equations (28) have a nontrivial solution if
Hence, the required characteristic equation of this dynamical system is given by [3,22] where Now, the second partial derivatives of the modified potential-like function U at triangular equilibrium points are

Advances in Astronomy
Hence, the general expressions for the roots of the characteristic equation (30) are , where

Numerical Application
In this section, we study numerically the locations and stability of the triangular equilibrium points L 4, 5 of the Sun-Earth system, by taking into account the dynamical flattening parameters of both the Sun (m 1 ) and Earth (m 2 ), the radiation pressure factor, and P-R drag of the Sun. For the purpose of computation in this paper, the astrophysical data of the Sun-Earth system are borrowed from NASA ADS, Ragos et al. [29], Mecheri et al. [30], and Singh and Umar [31]. e first two even zonal harmonics of the Sun are J A 2 ∼ 2.2 × 10 − 7 & J A 4 ∼ − 4.5 × 10 − 9 and those of the Earth are J B 2 ∼ 1.0 × 10 − 3 & J B 4 ∼ − 1.6 × 10 − 6 . Also, the orbital eccentricity of the Earth is e ∼ 0.0167. Now, the dynamical flattening parameters of the primary bodies are given by e radiation pressure factor of the Sun q is defined as q � 1 − α such that α can be expressed as α � (L ⊙ /(2πGm 1 cκ)) [28], where L ⊙ is the luminosity of the Sun, G is the gravitational constant, c is the speed of light, and κ is the mass per unit area. By using Stefan-Boltzmann's law, the luminosity of the primary can be expressed as where σ is the Stefan-Boltzmann's constant (see [28]). Also, the dimensionless velocity of light and the mass ratio of the Sun-Earth system are given by c d ∼ 10064.84 and μ ∼ 3.0035 × 10 − 6 , respectively [29].

Discussion
e modified equations of motion of an infinitesimal mass in the framework of the elliptic restricted three-body problem under the effects of dynamical flattening parameters of both primaries, radiation pressure factor, and P-R drag of the bigger primary (i.e., the Sun) are given in equations (11), while the locations and characteristic equation of the triangular points L 4, 5 are given in equations (24) and (34), respectively. Figure 1 shows the effects of the perturbing forces involved in the problem under consideration for the three different cases and classical case as well on the locations of triangular equilibrium points L 4, 5 . Graph (a) is the classical case, graph (b) shows the effects of radiation pressure factor together with P-R drag, and graph (c) shows the effects of dynamical flattening parameters while graph (d) shows the combine effects of radiation pressure factor, P-R drag, and dynamical flattening parameters. is clearly shows that the perturbing forces under consideration have no significant effect on the locations of the triangular equilibrium points L 4, 5 in the vicinity of the Sun-Earth system. ese effects can only be seen in the table (see Table 1). Table 2 shows numerical roots of the characteristic equation (30) for the classical case, as well as three other cases. In all cases, the characteristic roots reveal that all the roots are purely imaginary. is shows that the perturbing forces under consideration have no significant effect on the stability of the triangular equilibrium points L 4, 5 . Hence, the motion of an infinitesimal mass near the triangular equilibrium points L 4, 5 of the Sun-Earth system is stable in the linear sense under the influence of these perturbations.
For A 2 1 � B 2 1 � A 1 B 1 � 0 and W � 0, the present results of the triangular points in the circular case are in conformity with Singh and Taura [5] for p 2 � M b � 0 in their results and for α � W � 0, which also agrees with those of Abouelmagd et al. [11]; the difference in configuration of the primary bodies is responsible for the difference in sign. eir results are [5,11] Advances in Astronomy  Table 1: Effects of radiation pressure factor, P-R drag, and dynamical flattening parameters on the locations of triangular equilibrium points L 4, 5 of the Sun-Earth system. Table 2: Effects of radiation pressure factor, P-R drag, and dynamical flattening parameters on the stability of triangular equilibrium points L 4, 5 of the Sun-Earth system.
is is because we have used the modified mean motion (n)(equation (10)) which does not contain the eccentricity (e) and semi-major axis (a). However, by substituting e � 0 and a � 1, their results in the triangular case coincide fully with ours upon relaxing some parameters in our problem. e same applies to the second partial derivatives. eir results are [12,14] e P-R drag finds importance in the investigation of the stability of zodiacal cloud, orbital evolution of cometary meteor steams, asteroidal particles, and dust rings around planets. is model is applicable not only to the Sun-Earth system but also to other systems in both the solar and stellar systems as well.

Data Availability
e data used to support the findings of this study are available from the corresponding author upon request.

Disclosure
It is part of work/research for the partial fulfilment of Ph.D. in Mathematics by the authors, respectively.

Conflicts of Interest
e authors declare that they have no conflicts of interest.