Johannes Kepler studied positional data for the planet Mars and came up with three Laws between 1609 and 1619 which describe how a body orbits around another. All of Kepler’s calculation were done using geometry as it was the only mathematical tool available then.

In the late 17th century calculus was discovered independenly by Sir Isaac Newton and Gottfried Leibnitz. There was much controversy over which man discovered it first with Newton accusing Leibnitz of plagiarising his work. Calculus became a very powerful tool which enabled Newton to derive Kepler’s laws from first principles.

*The orbit of every planet is an ellipse with the Sun at
one of the two foci.*

The variables are:

- \(a\) is the semi major axis distance.
- \(b\) is the semi minor axis distance.
- \(e\) is the eccentricity of the ellipse.
- \(p\) is the semi latus rectum distance which passes through a focus.
- \(\theta\) is the real anomaly which is the angle the planet has orbited from perihelion.
- \(r\) is the distance from the planet to the sun.

The values \(a\), \(b\) and \(e\) are related by the equation:

\[b^2 = a^2(1-e^2)\]

The area of the ellipse \(A = \pi ab\).

The cartesian equation of an ellipse centred on its focus is:

\[\begin{equation} \frac{(x + ae)^2}{a^2} + \frac{y^2}{b^2} = 1 \end{equation}\]

In polar coordinates this is:

\[r = \frac{a(1 - e²)}{1 + e\cos\theta}\]

At the semi latus rectum \(x = 0\) and \(y = p\). Substituting these values and multiplying by \(b^2\) gives:

\[\begin{equation} e^2b^2 + p^2 = b^2 \implies p^2 = b^2(1 - e^2) \implies p^2 = a^2(1 - e^2)^2 \implies p = a(1 - e^2) \end{equation}\]

This gives the equation of the ellipse to be:

\[r = \frac{p}{1 + e\cos\theta}\]

Newton stated that the energy of a body orbiting another is the difference of the kinetic energy and the gravitational potential energy. This is a constant due to the conservation of energy:

\[E = \frac{m}{2}(\dot{r}^2 + r^2\dot{\theta}^2) - \frac{GMm}{r}\]

Where \(G\) is the gravitational constant, \(m\) is the mass of the orbiting body and \(M\) is the mass of the body being orbited. The angular momentum \(L\), which is also constant due to the conservation of angular momentum, is given by:

\[L = m r^2\dot{\theta}\]

Combining the two equations and integrating gives the equation of an ellipse.

\[r = \frac{p}{1 + e\cos\theta}\]

Where the semi latus rectum \(p\) is:

\[p = \frac{L^2}{GMm^2}\]

The eccentricity e is:

\[e^2 = 1 + \frac{2E}{GMm}p\]

*A line joining a planet and the Sun sweeps out equal
areas during equal intervals of time.*

Kepler determined how to describe the true anomaly in terms of time. This is not a straightforward process.

The mean anomaly \(M\) sweeps out the same area of the containing circle for equal intervals of time. The mean anomaly can be expressed as a function of time \(t\) given the orbital period \(P\) by:

\[M = \frac{2\pi t}{P}\]

The area swept out is:

\[A(t)=\frac{a^2M}{2}\]

Kepler defined an intermediate angle called the eccentric anomaly \(E\) to the point A on the diagram. It is chosen such that the area of the section \(\widehat{PSA}\) on the diagram is the same area as the area swept out by the mean anomaly. The area of \(\widehat{PSA}\) is the area of the sector \(\widehat{PZA}\) minus the area of the triangle \(\widehat{ASZ}\). Multiplying by two and dividing by \(a^2\) gives Kepler’s equation:

\[M = E - e\sin E\]

Kepler’s equation cannot be solved analytically. As the eccentricity is small for most planets, it can be approximated using a power series expansion. It can also be solved numerically using iteration.

The final step is to be able to derive the value of the true anomaly \(\theta\) from the eccentric anomaly \(E\). From the diagram we can see that:

\[a\cos E = ae + r\cos\theta\]

The radius \(r\) can be calculated from the eccentric anomaly by combining the ellipse equation with the true anomaly equation.

\[r = a(1 - e\cos E)\] \[\cos\theta = \frac{\cos E - e}{1 - e\cos E}\]

The angular momentum \(L\) of an orbiting particle is:

\[L = I\frac{d\theta}{dt}\]

Where \(I = r^2m\) is the moment of inertia of a particle of mass \(m\). The angular momentum is constant as there is no torque.

\[L = mr^2\frac{d\theta}{dt}\]

The area swept out \(A\) is given by:

\[\frac{dA}{dt} = \frac{1}{2}r^2\frac{d\theta}{dt}\]

Combining the equations gives:

\[\frac{dA}{dt} = \frac{L}{2m}\]

So, the rate of change of area with respect to time is constant.

*The square of the orbital period \(T\) is directly proportional to the
cube of the semi-major axis \(a\) of its orbit.*

Kepler’s equation for the third law is:

\[T^2 \propto a^3\]

Newton’s form of Kepler’s Third Law is:

\[T^2 = \frac{4\pi^2}{G(M + m)}a^3\]

Where \(M\) and \(m\) are the masses of the two bodies. In the case of our solar system the mass of the Sun is considerably greater than the masses of the planets. In general if \(M >> m\), then the \(m\) term can be ignored. For our solar system if \(T\) is measured in years and \(a\) is measured in AU then for any body orbiting the Sun:

\[T^2 = a^3\]