Sfdwhf

Kepler's equation solutions for five different eccentricities between 0 and 1

In orbital mechanics, Kepler's equation relates various geometric properties of the orbit of a body subject to a central force.

It was first derived by Johannes Kepler in 1609 in Chapter 60 of his Astronomia nova,^[1][2] and in book V of his Epitome of Copernican Astronomy (1621) Kepler proposed an iterative solution to the equation.^[3][4] The equation has played an important role in the history of both physics and mathematics, particularly classical celestial mechanics.

Equation

Kepler's equation is

M=E−esin⁡E{displaystyle M=E-esin E}

where M is the mean anomaly, E is the eccentric anomaly, and e is the eccentricity.

The 'eccentric anomaly' E is useful to compute the position of a point moving in a Keplerian orbit. As for instance, if the body passes the periastron at coordinates x = a(1 − e), y = 0, at time t = t₀, then to find out the position of the body at any time, you first calculate the mean anomaly M from the time and the mean motion n by the formula M = n(t − t₀), then solve the Kepler equation above to get E, then get the coordinates from:

x=a(cos⁡E−e)y=bsin⁡E{displaystyle {begin{array}{lcl}x&=&a(cos E-e)\y&=&bsin Eend{array}}}

where a is the semi-major axis, b the semi-minor axis.

Kepler's equation is a transcendental equation because sine is a transcendental function, meaning it cannot be solved for E algebraically. Numerical analysis and series expansions are generally required to evaluate E.

Alternate forms

There are several forms of Kepler's equation. Each form is associated with a specific type of orbit. The standard Kepler equation is used for elliptic orbits (0 ≤ e < 1). The hyperbolic Kepler equation is used for hyperbolic trajectories (e ≫ 1). The radial Kepler equation is used for linear (radial) trajectories (e = 1). Barker's equation is used for parabolic trajectories (e = 1).

When e = 0, the orbit is circular. Increasing e causes the circle to become elliptical. When e = 1, there are three possibilities:

• a parabolic trajectory,


• a trajectory going in or out along an infinite ray emanating from the centre of attraction,


• or a trajectory that goes back and forth along a line segment from the centre of attraction to a point at some distance away.

A slight increase in e above 1 results in a hyperbolic orbit with a turning angle of just under 180 degrees. Further increases reduce the turning angle, and as e goes to infinity, the orbit becomes a straight line of infinite length.

Hyperbolic Kepler equation

The Hyperbolic Kepler equation is:

M=esinh⁡H−H{displaystyle M=esinh H-H}

where H is the hyperbolic eccentric anomaly.
This equation is derived by redefining M to be the square root of −1 times the right-hand side of the elliptical equation:

M=i(E−esin⁡E){displaystyle M=ileft(E-esin Eright)}

(in which E is now imaginary) and then replacing E by iH.

Radial Kepler equation

The Radial Kepler equation is:

t(x)=sin−1⁡(x)−x(1−x){displaystyle t(x)=sin ^{-1}({sqrt {x}})-{sqrt {x(1-x)}}}

where t is proportional to time and x is proportional to the distance from the centre of attraction along the ray. This equation is derived by multiplying Kepler's equation by 1/2 and setting e to 1:

t(x)=12[E−sin⁡(E)].{displaystyle t(x)={frac {1}{2}}left[E-sin(E)right].}

and then making the substitution

E=2sin−1⁡(x).{displaystyle E=2sin ^{-1}({sqrt {x}}).}

Inverse problem

Calculating M for a given value of E is straightforward. However, solving for E when M is given can be considerably more challenging. There is no closed-form solution.

One can write an infinite series expression for the solution to Kepler's equation using Lagrange inversion, but the series does not converge for all combinations of e and M (see below).

Confusion over the solvability of Kepler's equation has persisted in the literature for four centuries.^[5] Kepler himself expressed doubt at the possibility of finding a general solution.

Inverse Kepler equation

The inverse Kepler equation is the solution of Kepler's equation for all real values of e{displaystyle e}:

E={∑n=1∞Mn3n!limθ→0+(dn−1dθn−1((θθ−sin⁡(θ)3)n)),e=1∑n=1∞Mnn!limθ→0+(dn−1dθn−1((θθ−esin⁡(θ))n)),e≠1{displaystyle E={begin{cases}displaystyle sum _{n=1}^{infty }{frac {M^{frac {n}{3}}}{n!}}lim _{theta to 0^{+}}!{Bigg (}{frac {mathrm {d} ^{,n-1}}{mathrm {d} theta ^{,n-1}}}{bigg (}{bigg (}{frac {theta }{sqrt[{3}]{theta -sin(theta )}}}{bigg )}^{!!!n}{bigg )}{Bigg )},&e=1\displaystyle sum _{n=1}^{infty }{frac {M^{n}}{n!}}lim _{theta to 0^{+}}!{Bigg (}{frac {mathrm {d} ^{,n-1}}{mathrm {d} theta ^{,n-1}}}{bigg (}{Big (}{frac {theta }{theta -esin(theta )}}{Big )}^{!n}{bigg )}{Bigg )},&eneq 1end{cases}}}

Evaluating this yields:

E={s+160s3+11400s5+125200s7+4317248000s9+12137207200000s11+15143912713500800000s13+⋯ with s=(6M)1/3,e=111−eM−e(1−e)4M33!+(9e2+e)(1−e)7M55!−(225e3+54e2+e)(1−e)10M77!+(11025e4+4131e3+243e2+e)(1−e)13M99!+⋯,e≠1{displaystyle E={begin{cases}displaystyle s+{frac {1}{60}}s^{3}+{frac {1}{1400}}s^{5}+{frac {1}{25200}}s^{7}+{frac {43}{17248000}}s^{9}+{frac {1213}{7207200000}}s^{11}+{frac {151439}{12713500800000}}s^{13}+cdots {text{ with }}s=(6M)^{1/3},&e=1\\displaystyle {frac {1}{1-e}}M-{frac {e}{(1-e)^{4}}}{frac {M^{3}}{3!}}+{frac {(9e^{2}+e)}{(1-e)^{7}}}{frac {M^{5}}{5!}}-{frac {(225e^{3}+54e^{2}+e)}{(1-e)^{10}}}{frac {M^{7}}{7!}}+{frac {(11025e^{4}+4131e^{3}+243e^{2}+e)}{(1-e)^{13}}}{frac {M^{9}}{9!}}+cdots ,&eneq 1end{cases}}}

These series can be reproduced in Mathematica with the InverseSeries operation.

InverseSeries[Series[M - Sin[M], {M, 0, 10}]]

InverseSeries[Series[M - e Sin[M], {M, 0, 10}]]

These functions are simple Maclaurin series. Such Taylor series representations of transcendental functions are considered to be definitions of those functions. Therefore, this solution is a formal definition of the inverse Kepler equation. However, E is not an entire function of M at a given non-zero e. The derivative

dM/dE=1−ecos⁡E{displaystyle dM/dE=1-ecos E}

goes to zero at an infinite set of complex numbers when e<1. There are solutions at E=±icosh−1⁡(1/e),{displaystyle E=pm icosh ^{-1}(1/e),} and at those values

M=E−esin⁡E=±i(cosh−1⁡(1/e)−1−e2){displaystyle M=E-esin E=pm ileft(cosh ^{-1}(1/e)-{sqrt {1-e^{2}}}right)}

(where inverse cosh is taken to be positive), and dE/dM goes to infinity at these points. This means that the radius of convergence of the Maclaurin series is cosh−1⁡(1/e)−1−e2{displaystyle cosh ^{-1}(1/e)-{sqrt {1-e^{2}}}} and the series will not converge for values of M larger than this. The series can also be used for the hyperbolic case, in which case the radius of convergence is cos−1⁡(1/e)−e2−1.{displaystyle cos ^{-1}(1/e)-{sqrt {e^{2}-1}}.} The series for when e = 1 converges when m < 2π.

While this solution is the simplest in a certain mathematical sense,^[which?], other solutions are preferable for most applications. Alternatively, Kepler's equation can be solved numerically.

The solution for e ≠ 1 was found by Karl Stumpff in 1968,^[7] but its significance wasn't recognized.^{[8][clarification needed]}

One can also write a Maclaurin series in e. This series does not converge when e is larger than the Laplace limit (about 0.66), regardless of the value of M (unless M is a multiple of 2π), but it converges for all M if e is less than the Laplace limit. The coefficients in the series, other than the first (which is simply M), depend on M in a periodic way with period 2π.

Inverse radial Kepler equation

The inverse radial Kepler equation (e = 1) can also be written as:

x(t)=∑n=1∞[limr→0+(t23nn!dn−1drn−1(rn(32(sin−1⁡(r)−r−r2))−23n))]{displaystyle x(t)=sum _{n=1}^{infty }left[lim _{rto 0^{+}}left({frac {t^{{frac {2}{3}}n}}{n!}}{frac {mathrm {d} ^{,n-1}}{mathrm {d} r^{,n-1}}}!left(r^{n}left({frac {3}{2}}{Big (}sin ^{-1}({sqrt {r}})-{sqrt {r-r^{2}}}{Big )}right)^{!-{frac {2}{3}}n}right)right)right]}

Evaluating this yields:

x(t)=p−15p2−3175p3−237875p4−18943931875p5−329321896875p6−241809262077640625p7− ⋯ |p=(32t)2/3{displaystyle x(t)=p-{frac {1}{5}}p^{2}-{frac {3}{175}}p^{3}-{frac {23}{7875}}p^{4}-{frac {1894}{3931875}}p^{5}-{frac {3293}{21896875}}p^{6}-{frac {2418092}{62077640625}}p^{7}- cdots {bigg |}{p=left({tfrac {3}{2}}tright)^{2/3}}}

To obtain this result using Mathematica:

InverseSeries[Series[ArcSin[Sqrt[t]] - Sqrt[(1 - t) t], {t, 0, 15}]]

Numerical approximation of inverse problem

For most applications, the inverse problem can be computed numerically by finding the root of the function:

f(E)=E−esin⁡(E)−M(t){displaystyle f(E)=E-esin(E)-M(t)}

This can be done iteratively via Newton's method:

En+1=En−f(En)f′(En)=En−En−esin⁡(En)−M(t)1−ecos⁡(En){displaystyle E_{n+1}=E_{n}-{frac {f(E_{n})}{f'(E_{n})}}=E_{n}-{frac {E_{n}-esin(E_{n})-M(t)}{1-ecos(E_{n})}}}

Note that E and M are in units of radians in this computation. This iteration is repeated until desired accuracy is obtained (e.g. when f(E) < desired accuracy). For most elliptical orbits an initial value of E₀ = M(t) is sufficient. For orbits with e > 0.8, an initial value of E₀ = π should be used. A similar approach can be used for the hyperbolic form of Kepler's equation.^[9]:66–67 In the case of a parabolic trajectory, Barker's equation is used.

Fixed-point iteration

A related method starts by noting that E=M+esin⁡E{displaystyle E=M+esin {E}}. Repeatedly substituting the expression on the right for the E{displaystyle E} on the right yields a simple fixed-point iteration algorithm for evaluating E(e,M){displaystyle E(e,M)}. This method is identical to Kepler's 1621 solution.^[4]

function E(e,M,n)

    E = M

    for k = 1 to n

        E = M + e*sin E

    next k

    return E

The number of iterations, n{displaystyle n}, depends on the value of e{displaystyle e}. The hyperbolic form similarly has H=esinh⁡H−M{displaystyle H=esinh H-M}.

This method is related to the Newton's method solution above in that

En+1=En−En−esin⁡(En)−M(t)1−ecos⁡(En)=En+(M+esin⁡En−En)(1+ecos⁡En)1−e2(cos⁡En)2{displaystyle E_{n+1}=E_{n}-{frac {E_{n}-esin(E_{n})-M(t)}{1-ecos(E_{n})}}=E_{n}+{frac {(M+esin {E_{n}}-E_{n})(1+ecos {E_{n}})}{1-e^{2}(cos {E_{n}})^{2}}}}

To first order in the small quantities M−En{displaystyle M-E_{n}} and e{displaystyle e},

En+1≈M+esin⁡En{displaystyle E_{n+1}approx M+esin {E_{n}}}.

See also

• Kepler's laws of planetary motion


• Kepler problem


• Kepler problem in general relativity


• Radial trajectory

References

External links

Wikimedia Commons has media related to Kepler's equation.

• 
  Danby, J. M.; Burkardt, T. M. (1983). "The solution of Kepler's equation. I". Cel. Mech. 31: 95–107. Bibcode:1983CeMec..31...95D. doi:10.1007/BF01686811.
  


• 
  Conway, B. A. (1986). An improved algorithm due to Laguere for the solution of Kepler's equation. doi:10.2514/6.1986-84.
  


• 
  Mikkola, Seppo (1987). "A cubic approximation for Kepler's equation". Cel. Mech. 40 (3). Bibcode:1987CeMec..40..329M. doi:10.1007/BF01235850.
  


• 
  Nijenhuis, Albert (1991). "Solving Kepler's equation with high efficiency and accuracy". Cel. Mech. Dyn. Astr. 51 (4): 319–330. Bibcode:1991CeMDA..51..319N. doi:10.1007/BF00052925.
  


• 
  Fukushima, Toshio (1996). "A method solving kepler's equation without transcendental function evaluations". Cel. Mech. Dyn. Astron. 66 (3): 309–319. Bibcode:1996CeMDA..66..309F. doi:10.1007/BF00049384.
  


• 
  Charles, E. D.; Tatum, J. B. (1997). "The convergence of Newton-Raphson iteration with Kepler's equation". Cel. Mech. Dyn. Astr. 69 (4): 357–372. Bibcode:1997CeMDA..69..357C. doi:10.1023/A:1008200607490.
  


• 
  Stumpf, Laura (1999). "Chaotic behaviour in the newton iterative function associated with kepler's equation". Cel. Mech. Dyn. Astr. 74 (2): 95–109. doi:10.1023/A:1008339416143.
  


• 
  Palacios, M. (2002). "Kepler equation and accelerated Newton method". J. Comput. Appl. Math. 138: 335–346. Bibcode:2002JCoAM.138..335P. doi:10.1016/S0377-0427(01)00369-7.
  


• 
  Boyd, John P. (2007). "Rootfinding for a transcendental equation without a first guess: Polynomialization of Kepler's equation through Chebyshev polynomial equation of the sine". Appl. Num. Math. 57 (1): 12–18. doi:10.1016/j.apnum.2005.11.010.
  

• Kepler's Equation at Wolfram Mathworld