A Solution of Kepler ’ s Equation

The present study deals with a traditional physical problem: the solution of the Kepler’s equation for all conics (ellipse, hyperbola or parabola). Solution of the universal Kepler’s equation in closed form is obtained with the help of the two-dimensional Laplace technique, expressing the universal functions as a function of the universal anomaly and the time. Combining these new expressions of the universal functions and their identities, we establish one biquadratic equation for universal anomaly ( ) χ for all conics; solving this new equation, we have a new exact solution of the present problem for the universal anomaly as a function of the time. The verifying of the universal Kepler’s equation and the traditional forms of Kepler’s equation from this new solution are discussed. The plots of the elliptic, hyperbolic or parabolic Keplerian orbits are also given, using this new solution.


Introduction
In the Keplerian problem, a body of mass 2 m follows a conic orbit, for which the focus is identified to the center of attracting body (with mass 1 m ).The Keplerian motion is described by the fundamental differential equation of the physical two-body problem: where r is the vector position of the moving body related to the attraction center and µ the gravitational pa- J. N. Tokis rameter defined by ( ) > (where G the Newtonian gravitational constant) (see [1], Equation (6.2.3)).
The traditional form of Kepler's equation, which can be obtained directly from Equation (1) (see [1], Section 6.3 and [2]), is normally written as: For elliptic orbits ( ) for hyperbolic orbits ( ) where ε is the eccentricity, E the eccentric anomaly and M the mean anomaly, which is defined as ( ) (see [2], Equation (11) and [3], Equations (4.35) and (4.51)).Remark that the time t is measured from peri- center and a is the semimajor axis, which is positive for ellipses, negative for hyperbolas, and infinite for pa- rabolas; also, q is the pericenter distance of the orbit.The case of 0 ε = leads to a circular orbit and the sim- ple solution e E M = (cf., Equation (2a)), so that we will regard 0 ε > hereafter in the present work.
Johannes Kepler announced the relevant laws of above equation early in 1609 and 1619 [4].He has used physics as a guide in this discovery [5].For four centuries, the Kepler's problem is to solve the nonlinear Kepler's Equation (2) for the eccentric anomaly.
Early analytical solution of Kepler's equation was considered in a comprehensive study of Tisserand [6].Recently, analytical works of the solution and use of Kepler's equation have been proposed by various authors (see, e.g., [7]- [12]).
In virtually every decade from 1650 to the present, there have appeared papers devoted to the solution of this Kepler's equation.Its exact analytical solution is unknown, and therefore, efficient procedures to solve it numerically have been well discussed in many standard text books of Celestial Mechanics and Astrodynamics as well as in a large number of papers.Colwell [13] contains extensive references to the Kepler problem in his book.During last two decades, studies were carried out by several investigators of the present problem [2] [14]- [18].In these studies, they used numerical or approximations methods for solution of the Kepler's equation.Hence, it appears that an analytical solution of the Kepler's equation will be of great interest.
In the current study, an analytical investigation of the Kepler's equation real roots in closed form is presented.In Section 2, we will establish the general form of Kepler's equation and will clear up the useful identities of the universal functions.In Section 3, using the two-dimensional Laplace transform technique, we will present an analytical solution for the universal Kepler's equation, obtaining the universal functions ( ) , 0,1, 2,3 n U n = as function of the universal anomaly ( ) χ and the time ( ) In Section 4, in the first step we will establish one new biquadratic equation for universal anomaly χ for all conics with the help of Newman's equation (cf.Eq- uation (23)) and some identities of the new expressions of universal functions.Then, the solution = of the present problem has been obtained, solving this biquadratic equation for all conics.
Finally, discussion of the results, thus obtained, is presented in Section 5; the new solution of the problem will prove that verifies the traditional form of Kepler's equations for elliptic, hyperbolic or parabolic orbits.The elliptic, hyperbolic or parabolic Keplerian motion is easily plotted, using this new solution.
Working for the Kepler's Equation (2), we consider an object following a path of same eccentricity ε about the center of attracting body; the object is at time t in (vector) position ( ) r t with (vector) velocity 0 v .The time t is measured from the pericenter passage; so, when 0 0 t = this object was at position 0 r (with 0 0 r = r ) of the pericenter with velocity 0 v (with 0 0 υ = v ) and eccentric anomaly 0 0 E = .We emphasize that the vec- tors r and 0 r originate at the center of attraction.Then, we introduce the universal anomaly χ , which is de- fined by Sundman transformation: (see [3], Equation (4.71)) and related to the classical eccentric anomaly by where E is the eccentric anomaly angle of elliptic or hyperbolic orbit and D the parabolic eccentric anomaly of parabolic orbit with dimension 1 2   L (see [9], Equation ( 18)).The α denotes the reciprocal of the semimajor axis a , namely Depending on the sign of α or the value of the eccentricity ε , the type of the orbit is determined such that: 0 α > (or 1 ε < ) for elliptic orbits; 0 α < (or 1 ε > ) for hyperbolic orbits and 0 α = (or 1 ε = ) for para- bolic orbits.Note that the universal anomaly χ is a new independent variable with dimension 1 2 L (see [3], Equation (4.70)).
From the initial condition 0 0 E = and the known relations: for hyperbolic orbits (see [3], Sections 4.3 and 4.4), we have also for the present problem ( ) where q stands for the pericenter distance of the orbit related to the parameter p with the relation Remark that the p is a non-negative parameter and the pericenter distance q may be positive or zero; both of them have dimensions of length (see [3], Section 4.1).Now, using the universal functions defined by with their following useful properties: ( ) ( ) ; ; U aU (see [3], Equation (9.73)), the two forms of Kepler's Equation ( 2) are incorporated in one universal equation which is a standard form of the traditional Kepler's Equations (2) with the epoch at pericenter passage (see [9], Equation ( 23) and [2], Equation ( 29)).The general formula (10) is valid for all values of α and ε ; in particu- lar, it is good for parabolic orbits where 1 0 a = = a .
To find out the expression of many orbital quantities, e.g. the magnitude of the position vector ( ) r t , we must solve the standard universal form (10) of the Kepler's equation for the universal anomaly as function of the time, namely

Solution of the Universal Kepler's Equation
In order to obtain the analytical solution of the present problem, we shall solve first the universal Kepler's equation Equation ( 10), obtaining the universal functions n U , 0,1, 2,3 n = as a function of the universal anomaly and the time.For this purpose, we will use the double Laplace transformation technique, which was analytically studied by Aghili and Salkhordeh-Moghaddam [19] and by Valkó and Abate [20].
The universal functions 2 U : For this case, we introduce a new variable ( ) (see Equations ( 9)) and the universal Kepler's Equation ( 10) becomes From the initial condition 0 E = or 0 χ = and Equation (9a), we have the corresponding initial conditions to the Equation (12a) The application of double Laplace transform (with respect to χ and t ) to the Equations ( 12) gives the so- lution ( ) , 1 s qs where s and ξ are the transform variables of χ and time t , respectively.For the solution of partial diffe- rential Equation (12a), we use the Appendix.Now, the universal function 2 U can be obtained by taking the inverse transform of Equation ( 13) (cf., Ap- pendix).So, we get ( ) where we have abbreviated ( ) A L = .The universal functions 3 U : For this case, we will use a new variable ( ) (see Equations (9c)) and the universal Kepler's Equation (10) becomes For 0 E = or 0 χ = and Equation (9a), we have the initial conditions to the Equation (17a) Similarly as in the case of 2 U , using the double Laplace transform (with respect to χ and t ) for the Equa- tions (17), we obtained the solution ( ) (cf., Appendix).Inverting ( ) where we have defined the non-dimensional function ( ) Then, substituting the results of Equation ( 14) and ( 19) into the relations where A and B are given from Equations ( 15) and ( 20), respectively.

Analytical Solution of the Problem
In order to obtain a solution for the universal anomaly ( ) , we use the explicit expression for χ : where . This relation was discovered by C. M. Newman and its χ does not involve any of the universal functions n U (see [3], Equation (4.86)).The σ can be also given by the known equation (see [3], Substituting Equation (24) (with 1 U given by ( 21)) into Equation (23), the following relation is obtained To find out two more relations between χ , A and B , similar to Equation (25), we will use the basic rela- tion ( ) ( ) sin cos 1 q q χ χ + = and the definitions of A and B Equations ( 15) and ( 20), respectively.Thus, we can easily obtain the relation Further, we can find one more relation using the basic identity ( ) [3], Equation (4.93)); indeed, with the help of Equations ( 14), ( 21) and ( 22), we obtain the relation ( ) The three Equations ( 25), ( 26) and (27) are a system of the three unknowns: χ , A and B .Solving this system, we get from two Equations ( 25) and (26) the following relations ( ) ( ) Finally, substituting Equations (28) into Equation ( 27), we obtain the following biquadratic equation for universal anomaly χ where we have abbreviated The solution of the biquadratic Equation (29) gives the relation between the universal anomaly χ and the time t for all conics: ellipse, hyperbola or parabola.The solution of this equation can be obtained using the standard formula of the solution for biquadratic equation.
Solving the new biquadratic Equation (29), we get the solution of the present problem for the universal anomaly χ at time 0 t ≥ as shown below ( ) where we have abbreviated particularly, for elliptic orbits ( ) for hyperbolic orbits ( ) and for parabolic orbits ( ) Remark that the discriminant of the biquadratic Equation ( 29) is ( ) Consequently, the solution (31) is real in the cases of hyperbolic and parabolic cases since the corresponding ( ) t ϕ , given from Equations (32c,d), are always real-valued.In the case of the elliptic Keplerian orbits, the Equ- ation (31) is real only in the special case for which ( ) t ϕ , given from Equation (32b), becomes real, namely for the case ( ) ( ) where e M is the mean anomaly defined by Equation (2c).The upper limit of this mean anomaly e M is de- fined from the relation (32f) for every completed trip of the orbiting body in its elliptic orbit about the center of the attracting body; this limit is useful for determination of each Keplerian ellipse (cf., the applications in the next section).
In the case of parabolic orbits where the limiting case → ∞ a (or 0 α = ) and 1 ε → corresponds to 0 2 r q p = = and , the equation ( 29) and its solution (31) are reduced to The Equation (31) is the solution of the present problem for all conics (ellipse, hyperbola or parabola) and expresses the relation between the universal anomaly χ and the time t.
Knowing the solution of the universal anomaly ( ) t χ , we establish the exact expressions of the universal functions n U , n = 0, 1, 2, 3 as functions of the time t .Indeed, using Equation (31), the Equation ( 28) are ob- tained as functions of the time t : Then, the universal functions ( 19), ( 14), ( 21) and ( 22) are expressed as functions of the time t as show below ( ) ( ) The magnitude of the position vector r of the orbiting body is (see [3], Equation (4.82) and [14], Equation ( 8)).Substituting into Equation (36b) the new expressions of 2 U , given by Equation (35b), we obtain the time-dependent distance r of the orbiting body from the center of attraction as Furthermore, if we work in the orbital reference system with the origin at the attracting center (or focus), we chose the ( ) , X Y plane to be the plane of motion with X -axis pointing toward pericenter and Y -axis in the direction for which the true anomaly ( ) f is 90˚; in this way the Z -axis is parallel to the angular momentum.
Then, if X and Y are the coordinates of a point ( ) , P X Y on the conic orbit, we have cos X r f = and r p X where p is the non-negative quantity given by Equation ( 7) (see [3], Equation (4.3)).Thus, the X and Y coordinates can be expressed as function of the time t , using Equations (38) and the identities: where 1 U and 2 U are given from Equations (35) (see [14], Equations ( 6)-( 7)).

Discussion
Using the standard form of the universal Kepler's equation (10) with the epoch at pericenter passage, we have derived a new biquadratic equation ( 29) for universal anomaly χ as a function of the time t .This equation governs the motion of an orbiting body following a path with position vector ( ) r t from the center of attracting body and with velocity vector ( ) The new solution (31) of the present problem was obtained solving Equation ( 29) with initial-value conditions for the orbiting body at time 0 t = passages from the pericenter, with minimum position vector 0 r , velocity vector 0 v and eccentric anomaly 0 0 E = .The solution (31) is a solution of the present problem.Indeed, the new expressions of the universal functions 3 U , and 1 U given by (35a,c) verify universal Kepler's Equation (10).Moreover, the solution (31) is a solution of the Equation ( 29), since it can be easily verified by substitution of (31) into (29).
The solution (31) verifies also the traditional forms of Kepler's Equation ( 2).Particularly, in the case of elliptic orbits ( ) the solution (31) of the problem is reduced for the eccentric anomaly (cf.Equation (4a)) to the form ( ) where e M is the mean anomaly given by Equation ( ( ) ( ) Similarly, in the case of hyperbolic orbits ( ) < the solution (31) of the problem is reduced, for the eccentric anomaly (cf.Equation (4b)), to ( ) ( ) where h M is the mean anomaly given by Equation ( Finally, in the case of parabolic orbits ( ) 1 ε = the solution (33b) gives for the parabolic eccentric anomaly (cf.Equations (4c) and (2e)): ( ) ( ) ( ) From the other hand, the standard and hyperbolic trigonometric functions of Equations ( 2) are expressed as ( ) ( ) where we have use Equation ( 35) and the relationship between the function U 1 and the standard and hyperbolic trigonometric functions of Equations ( 2): ( ) ( ) ( ) and hyperbola, respectively (see [3], Problem 4 -21 and [9], Equation ( 24)).Now, we will prove that the Equations ( 40) and (42) represent the solutions of the traditional forms of Kepler's Equations (2).Indeed, it can be shown that the left-hand sides of Equations ( 2) are reduced to the righthand sides, namely ( ) ( ) in accordance with the Equations ( 40), ( 42) and (45).
It should be pointed out that our solutions for eccentric anomaly (cf., Equations ( 40) and ( 42)) are ready for physical applications in the corresponding Keplerian orbits.
In addition to above Keplerian orbits, the new solution (33b) of the present problem for the case of parabola ( 1 ε = and 0 α = ) verifies the traditional Barker's equation for parabolic orbits Remark that the parabolic Keplerian equation is called Barker's equation (see [3], Equation (4.24)).Indeed, from the definition (8) for the universal functions we have ( ) Using the solution (33b), these universal functions are expressed from Equations (35a,c), so that ( ) ( ) In order to study the Keplerian orbits with the help of the new solutions, we use also the cartesian coordinates x and y with the origin at the center of the ellipse or hyperbola.In this case the X and Y coordinates of the orbital ( ) , X Y plane system given by Equations (39) are related to the system ( ) , x y with the relations X x ε = − a and Y y = (see [3], Equation (4.4)); so, for the case 1 ε ≠ , these new coordinates can be obtained in the explicit forms where we have defined the non-dimensional relation ( ) in accordance to the Equations (39).Then, we introduce the non-dimensional coordinates The new expressions (52) verify the following equations of ellipse ( ) and hyperbola ( ) where ζ is defined by Equation (51).Remark that the equations (53) are the non-dimensional forms of the el- lipse and hyperbola, respectively, in the cartesian system ( ) , x y with the origin at the center of the ellipse or hyperbola (see [1], Equations (A.2.1) and (A.4.2)).
In the other hand, for the case of parabola ( ) 1 ε = , the coordinates X and Y can be also obtained, from Equations (39), as following ( ) Then, we have, from Equations (54), ( ) ( ) The last Equation ( 55) is the equation of parabola, which passes through its pericenter with coordinates ( ) with the non-dimensional coordinates ( ) In addition to above results for the non-dimensional coordinates we have (for physical application) the corresponding expressions: For elliptic orbits ( ) for hyperbolic orbits ( ) for parabolic orbits ( )  2c,d,e)).These mean anomalies can be varied from 0 to 2 π k , where k is an integer.Note that the non- dimensional coordinates of attractive center (or focus) F in the non-dimensional cartesian system ( ) given by ( ) 2 ,0 F ε for all Keplerian orbits (ellipse, hyperbola or parabola).In order to get a physical insight into the new solution of the Kepler's problem, we apply the above results for the system Earth-Moon.For this system the eccentricity of the Moon is 0.0549 ε = and the upper limit for the mean anomaly can be obtained from relation (32f) as 1.94432 ) and using Equations (58a) and (52b), the elliptic Keplerian motion of the Moon about the Earth can be easily plotted in the non-dimensional cartesian system ( ) , ζ η (Figure 1).Remark that the non-dimen- sional coordinates of the Earth are obtained as ( ) 0.003, 0 .Note that the use of the upper limit of the mean anomaly, given from the relation ( )

Conclusions
This work presents a solution to the well known Keplerian two body physical problem.From the investigation relation (32f)), is import for the plotting of all elliptic Keplerian orbits.To confirm that we give two more examples: (a) We consider an object following an elliptic orbit with eccentricity 0.5 ε = about the center of attracting body; the upper limit of the mean anomaly, from relation (32f), is1.4142 e M < ; then, varying the mean anomaly e M (with 1.414 e M ≤ ) and using Equations (58a) and (52b), we plot the elliptic Keplerian orbit in the non-dimensional cartesian system ( ) , ζ η (Figure2). (b) We consider another object in an elliptic orbit with

Figure 4 .
Figure 4.The parabola of an orbiting body (with 1 ε = ) about the attractive