Senior Project

Andy Erickson
May 15, 2009

1 Abstract
In this paper, Kepler’s three Laws of Planetary Motion are proven using New-
ton’s Law of Universal Gravitation. In addition, several results pertaining to
the orbital period of a satellite are derived. An equation for the velocity of a
satellite, as well as the minimum and maximum velocities necessary for a satel-
lite to stay in orbit are also derived. Finally, the anomalous orbit of Mercury is
examined using Newton’s Law of Universal Gravitation and Einstein’s Theory
of Relativity. This section assumes a basic familiarity with General Relativity
though a knowledge of tensor calculus is not required to follow the analysis of
Mercury’s orbit.

2 Introduction
For the past two millennia, people have endeavored to mathematically model
the motions of the planets. The goal of predicting these motions began with
the Greeks around 4 BC. Over the course of the last 2000 years, the model of
the universe has changed from geocentric to heliocentric to a universe in which
neither space nor time are constant and our Solar System is of little consequence
(Linton, 1).
The first significant challenge to the 1500 year old geocentric model of the
universe occured in 1543 when Nicholas Copernicus published On the Revolutions
of the Heavenly Spheres (Linton, 119). Thus began a period of rapid develop-
ment in mathematical astronomy. At the beginning of the 17th century, Jo-
hannes Kepler, using Tycho Brahe’s observational data and the Copernican
model of a heliocentric solar system, derived his laws of planetary motion (Lin-
ton, 177). These three laws are:
(1) A planet revolves in an elliptic path with the sun as one of the foci of the
ellipse
(2) The radius vector from the sun to a planet sweeps out equal areas in equal
intervals of time.
(3) The squares of the periods of revolution of the planets around the sun are
proportional to the cubes of their mean distances from the sun.

1
 The next significant development in mathematical astronomy was Issac New-
ton’s formulation of a law of attraction governing massive objects. Newton’s
masterpiece, the Mathematical Principles of Natural Philosophy, known as the
Principia was published in three volumes over a nearly 40 year period. In 1687,
Newton published the first volume which contained his three laws of motion:

(1) Every body perseveres in its state of being at rest or of moving uniformly
straight forward, except insofar as it is compelled to change its state by
forces impressed upon it. (Inertia)
(2) A change is proportional to the motive force impressed and takes place
along the straight line in which that force is impressed. (F = ma)

(3) To any action there is always an opposite and equal reaction; in other
words, the actions of two bodies upon each other are always equal and
always opposite in direction. (Conservation of Motion) (Linton, 263).
In 1726, Newton published the third book of his Principia, which contained
his statement of the law of universal gravitation and the resulting motion of
planetary bodies in the solar system. Newton’s Law of Universal Gravitation
states
GM m
F =
r2
where G is the gravitational constant, M and m are the mass of each body, and
r is the distance between the centers of each mass. From this force law and his
laws of motion, Newton was able to prove Kepler’s laws of planetary motion
(Linton, 272).
Over the next 150 years, Newton’s Law of Universal Gravitation was shown,
with few exceptions, to adequately model the motions of the planetary bodies in
our solar system. Cheif among the problems with Newton’s Theory of Gravity,
was the discrepancy between the predicted and observed motion of Mecury. By
1850, tables of the motion of Mercury were still shockingly inaccurate relative to
the tables of the other planets, and the anomolous motion of Mercury became
a major problem in astronomy.
During this time, Urbain Jean Joseph Le Verrier began to work on the
problem of Mercury’s orbit. It was not until 1859 that Le Verrier was able to
correctly model the orbit of Mercury, although he could not explain the anomaly.
Le Verrier discovered that the error in predictions of Mercury’s orbit was a result
of the advance of its perihelion, the closest point on the orbit to the Sun. By
analyzing the effect of the other planets on the perihelion advance of Mercury,
Le Verrier discovered that the advance of the perihelion predicted by Newton’s
Theory of Gravity was less than the observed perihelion advance, although he
was unable to explain the cause of this anomaly (Roseveare, 20-37).
Though Le Verrier had successfully corrected the tables of Mercury’s orbit,
explaining the cause of this inconsistency would take another 50 years and a new
Theory of Gravity. In developing his Theories of Special and General Relativity,
Einstein was not immediately concerned with solving the anomalies of planetary

2
motion. However, after developing the Theory of General Relativity, Einstein
realized that it explained the anomalous advance of Mercury’s perihelion. This
observational verification of General Relativity, as well as gravitational lensing
of light and gravitational redshift of light, proved valuable in cementing General
Relativity as a viable theory. In correcting the difference between the predicted
and observed advance of Mercury’s perihelion, Einstein brought to a close nearly
100 years of mathematical quarreling over the cause of Mercury’s anomalous
perihelion and changed our conception of the universe (Roseveare, 147-186).
In this paper, we replicate the work of Newton in proving Kepler’s three
Laws of Planetary Motion using Newton’s Law of Universal Gravitation. In
addition, we will also verify other minor results involving satellite orbit. We
also examine the differences in planetary orbits as predicted by Newton’s Law
of Universal Gravitation and Einstein’s Theory of Relativity.

3 Newton’s Law of Universal Gravitation

As previously mentioned, Newton’s Law of Universal Gravitation states that
the force governing two massive bodies is inversely proportional to the square
of the distance between their centers of mass and directly proportional to the
product of their masses. In mathematical terms, we have
GM m
F = (1)
r2
where G is the gravitational constant, M and m are the masses of each body,
and r is the distance between the centers of each mass.

3.1 Kepler’s Laws

For the purpose of verifying Kepler’s Laws, we will find it useful to choose our
reference frame such that one of the objects is stationary. We choose the larger
mass M to be stationary and the smaller mass m to be in orbit around the
larger mass. In Cartesian coordinates, the position of mass m is given by (x, y).
The gravitational force is directed towards the origin, and has magnitude,
GM m GM m
|F~ | = = . (2)
x2 + y 2 r2

3.1.1 Dynamic Equations of Motion

We have that the gravitational force is given by
GM m GM m
F~ = − r̂ = − ~r,
r2 r3

3
where r̂ is the unit vector in the direction of ~r. If we resolve the force into the
x and y directions, we have
GM m
Fx = − x and (3a)
r3
GM m
Fy = − y. (3b)
r3
We apply Newton’s second law (F = ma) to each component, obtaining the
dynamic equations of motion for the satellite,
d2 x GM m
m 2
=− x and (4a)
dt r3
d2 y GM m
m 2 =− y. (4b)
dt r3
We cancel the mass of the satellite in both equations to obtain
d2 x GM
2
= − 3 x and (5a)
dt r
d2 y GM
= − 3 y. (5b)
dt2 r

3.1.2 Kepler’s First Law

We now prove that the shape of a closed orbit of a satellite around a planet is
elliptic. From equations 5a and 5b, we derive a formula for angular momentum.
We will show that this formula is equal to a constant, which implies that angular
momentum is conserved through a planet’s orbit. We multiply equation 5a by
y and 5b by x and subtract 5b from 5a yielding
d2 y d2 x
 
d dy dx GM GM
x 2 −y 2 = x −y = 3 xy − 3 xy = 0. (6a)
dt dt dt dt dt r r
Integrating with respect to t yields
dy dx
x −y =H (6b)
dt dt
for some constant H. This result is the conservation of angular momentum
where H is the conserved value of angular momentum per unit mass of the
dy
orbiting motion. We now multiply 5a by dx
dt and 5b by dt and add them together
which yields
dx d2 x dy d2 y
       
dx GM dy GM
+ = − 3 x + − 3 y .
dt dt2 dt dt2 dt r dt r
We can now rearrange terms on the right side of the equation
dx d2 x dy d2 y
     
GM dx dy
+ =− 3 x +y .
dt dt2 dt dt2 r dt dt

4
For simplicity, in the right side of the equation, we can convert to polar coor-
dinates, such that x = r cos θ and y = r sin θ. Likewise, dx = cos θdr − r sin θdθ
and dy = sin θdr + r cos θdθ. Therefore, we have

dx d2 x dy d2 y
      
GM dr dθ
+ = − r cos θ cos θ − r sin θ
dt dt2 dt dt2 r3 dt dt
 
dr dθ
+ r sin θ sin θ + r cos θ
dt dt
 
GM dr dr
= − 3 r cos2 θ + r sin2 θ .
r dt dt

This simplifies to

d2 x d2 y
   
dx dy GM dr
+ =− (7a)
dt dt2 dt dt2 r2 dt

Next, we rearrange the left side of Equation 7a, yielding

   
dx d dx dy d dy GM dr
+ =− 2 .
dt dt dt dt dt dt r dt
We recognize that the left side of the equation is the derivative of the sum of
the squares of the first derivatives of x and y. Using the Chain Rule, we find
we need a factor of 12 , and we now have
"   2 #
2
1 d dx dy GM dr
+ =− 2 .
2 dt dt dt r dt

We can now integrate with respect to time, which produces

"   2 #
2
1 dx dy GM
+ − = E, (7b)
2 dt dt r

for some constant E. We can see that the first terms of equation 7b is a statement
of kinetic energy per unit mass. In addition, the second term is a statement of
potential energy per unit mass. Thus the constant E represents the total energy
per unit mass. Thus total energy is conserved in the system. Because there are
no external forces on the two body system, conservation of energy is expected.
For convenience, we now convert equations 6b and 7b to polar coordinates.
The transformation of coordinates is given by

x = r cos θ and y = r sin θ (8a)

with the derivatives given by

dx = cos θdr − r sin θdθ and dy = sin θdr + r cos θ. (8b)

5
 We begin with equation 6b, we substitute for x and y, and rearrange terms.
We then have
   
dr dθ dr dθ
r cos θ sin θ + r cos θ − r sin θ cos θ − r sin θ =H
dt dt dt dt
dθ dθ
r2 cos2 θ + r2 sin2 θ =H
dt dt
dθ
r2 = H. (9a)
dt
Next we substitute for x and y in equation 7b which yields
" 2  2 #
1 dr dθ dr dθ GM
cos θ − r sin θ + sin θ + r cos θ − =E
2 dt dt dt dt r
"  2  2 #
1 2 dr dr dθ 2 2 dθ
cos θ − 2r cos θ sin θ + r sin θ
2 dt dt dt dt
"  2  2 #
1 dr dr dθ dθ GM
+ sin2 θ + 2r cos θ sin θ + r2 cos2 θ − =E
2 dt dt dt dt r
"   2 #
2
1 dr 2 dθ GM
+r − = E. (9b)
2 dt dt r

We then take the initial position of the satellite to be r0 and the initial
velocity to be v0 . We also define the angle between the initial position vector
and the initial velocity vector to be α.

Figure 1: Initial Conditions of Launching (Kwok)

Thus, for conservation of angular momentum, we have

 
dθ
H = r02 = r0 v0 sin α (10a)
dt t=0
and for conservation of energy, we have
v02 GM
E= − . (10b)
2 r0

6
 Substituting these values into equations 9a and 9b, we find
dθ
r2 = r0 v0 sin α, and (11a)
dt
 2  2  
dr dθ 1 1
+ r2 = v02 + 2GM − . (11b)
dt dt r r0
By eliminating dependence on t by using dr dr dθ
dt = dθ dt , equation 11b can be
transformed to
 2  2  2
dr dθ 2 dθ 2GM 2GM
+r = v02 + − .
dθ dt dt r r0
dθ
We substitute Equation 11a for dt which yields
 2  2  2
dr r0 v0 sin α 2 r0 v0 sin α 2GM 2GM
+r = v02 + − .
dθ r2 r2 r r0

dr 2


We then rearrange terms and solve for dθ to get.
2
r4 r2 v 2 sin2 α
  
dr 2GM 2GM
= v02 − + − 0 0 2 .
dθ r0 v0 sin2 α
2 2 r0 r r
dr
We take the square root of both sides and simplify to solve for dθ giving
s
dr r2 2GM 2GM r2 v 2 sin2 α
= v02 − + − 0 0 2 . (12a)
dθ r0 v0 sin α r0 r r

We now have a first order separable differential equation. We can begin the
arduous process of solving for θ in terms of r. It is advantageous to substitute
r = z1 and thus dr = − dz
z 2 . We then have
s 
1 dz 1 2GM
− 2 = v02 − + 2GM z − v02 r02 sin2 αz 2 .
z dθ v0 r0 z 2 sin α r0

Solving for dθ yields

v0 r0 sin αdz
−dθ = r  . (12b)
2 2GM 2 2 2 2
v 0 − r0 + 2GM z − v0 r0 sin αz

We now complete the square in the denominator in order to substitute a trigono-

metric function before integrating. The complete derivation of this is included
in Appendix A, but can be ignored for continuity. From completing the square
in the denominator of equation 12b, we have

7
 v0 r0 sin αdz
−dθ = r  2 . (12c)
(G2 M 2 −2GM v02 r0 sin2 α+v04 r02 sin2 α) GM
v02 r02 sin2 α
− v0 r0 sin αz − v0 r0 sin α

We rearrange terms and find a suitable trigonometric substitution in order to

integrate. We rewrite equation 12c as
v0 r0 sin αdz
−dθ = q r   2 .
(G2 M 2 −2GM v02 r0 sin2 α+v04 r02 sin2 α) v02 r02 sin2 α GM
v02 r02 sin2 α
1− G2 M 2 −2GM v02 r0 sin2 α+v04 r02 sin2 α
v0 r0 sin αz − v0 r0 sin α

We use the trigonometric substitution

 
v0 r0 sin α GM
q v0 r0 sin αz − = cos U
G2 M 2 − 2GM v02 r0 sin2 α + v04 r02 sin2 α v0 r0 sin α

and consequently,
q
G2 M 2 − 2GM v02 r0 sin2 α + v04 r02 sin2 α
dz = − sin(U )dU.
v02 r02 sin2 α

We then substitute into equation 12c, which yields

sin(U )dU
dθ = p = dU.
1 − cos2 (U )

Integrating from θ0 to θ yields

U = θ − θ0 .

We now reverse substitute for U. We know

 
 
v 0 r0 sin α GM
U = arccos  q v0 r0 sin αz − .
2 2 2 2 4 2 2
G M − 2GM v0 r0 sin α + v0 r0 sin α v 0 r0 sin α

Taking the cosine of both sides yields

v02 r02 sin2 α

 
GM
cos(θ − θ0 ) = q z− .
G2 M 2 − 2GM v02 r0 sin2 α + v04 r02 sin2 α v02 r02 sin2 α

We then solve for z. Rearranging terms yields

s
GM 2GM G2 M 2
v0 r0 sin αz = + v02 − + 2 2 2 cos(θ − θ0 ).
v0 r0 sin α r0 v0 r0 sin α

8
We substitute z = 1/r and solve for r yielding

v0 r0 sin α
r= q .
GM 2GM G2 M 2
v0 r0 sin α + v02 − r0 + v02 r02 sin2 α
cos(θ − θ0 )

Finally, we can simplify and find

v02 r02 sin2 α/GM

r= q . (13)
v0 r0 sin α 2GM G2 M 2
1+ GM v02 − r0 + v02 r02 sin2 α
cos(θ − θ0 )

We now have the orbital path in polar form, r as a function of θ. To understand

this orbital path geometrically, we would like to associate the solution with the
equation of a conic section. In polar form, the equation of a conic section is
pe
r= , (14)
1 + e cos θ
where e is the eccentricity of the conic and p is the distance from the focus to
the directrix. The conic is an ellipse when e < 1, a parabola when e = 1, and a
hyperbola when e > 1.

Figure 2: Elliptical Orbit of a Planet about the Sun

Comparing equations 13 and 14 and taking θ0 = 0 (the geometric interpre-

tation of this is that the initial position vector lies along the x-axis), we find
that
s
v 2 r2 sin2 α 2GM
 
e = 1 − 0 02 2 − v02 , (15a)
G M r0
v02 r02 sin2 α
pe = (15b)
G2 M 2
We now have that the solution of the orbital path of a satellite is a conic section
with the type of orbit (i.e. elliptic, parabolic, or hyperbolic) governed by the
eccentricity e. Of these orbits, the only closed orbit is elliptic, and thus we have
Kepler’s First Law: a planet revolves in elliptic path with the sun as one of the
foci of the ellipse.

9
3.1.3 Kepler’s Second Law
We now prove Kepler’s Second Law, which states that the radius vector from
the sun to a planet sweeps out equal areas in equal intervals of time. As we
seen in Figure 1, the physical manifestation of this property is that the satellite
has the greatest velocity when it is closest to the central mass and least velocity
when it is furthest from the central mass.

Figure 3: Kepler’s Second Law (Wikipedia)

We have already defined ~r as the vector from the central mass to the satellite
and additionally we say r = |~r|. We now let d~r be the vector tangent to the
orbit over which the satellite moves in time dt. Thus,
1
dA = |~r × d~r|.
2
or equivalently
1 dθ
dA = |~r × ~r dt|.
2 dt
We simplify to get
dA 1 dθ 1 dθ
= |~r × ~r| = |~r||~r| .
dt 2 dt 2 dt
Thus,
dA 1 dθ
= r2 .
dt 2 dt
But from Equation 9a, we know
dθ
r2 = H.
dt
Thus,
dA H
=
dt 2
and as a satellite moves around an object, the position vector ~r sweeps out equal
areas in equal times. In addition, we see that Kepler’s Second Law is simply
another way of describing conservation of angular momentum.

10
3.1.4 Kepler’s Third Law
Finally, we prove Kepler’s Third Law, which states the squares of the periods of
revolution of the planets around the sun are proportional to the cubes of their
mean distances from the sun. Using Kepler’s Second Law and equation 10a, we
have
dA 1 dθ H 1
= r2 = = r0 v0 sin α. (16)
dt 2 dt 2 2
In Section 2.2, we showed that the path of a revolving satellite is elliptic.
We also know that the area of an ellipse is πab, where a and b are the lengths
of the semi-major axis and semi-minor axis respectively. To find the period of
revolution T , we integrate equation 16 from t = 0 to t = T yielding
Z T
dA
πab = A(T ) − A(0) = dt
0 dt
Z T
1
= r0 v0 sin αdt
0 2
1
= r0 v0 sin αT. (17)
2
From our initial conditions for deriving the orbital path and Equation 15b,
we have that rmin = pe/(1 + e) occurs when θ = 0 and the maximum value
rmax = pe/(1 − e) occurs when θ = π. Therefore the length of the semi-major
axis, a is
max + min pe GM
a= = = 2GM (18a)
2 1 − e2 r − v02 0

and similarly, the length of the semi-major axis is

max − min p v0 r0 sin α
b= = a 1 − e2 = q . (18b)
2 2GM
− v2 r0 0

We can now substitute equations 18a and 18b into equation 17, and obtain the
period of revolution of the elliptic path of a orbiting satellite,

2πGM 2πa3/2
T = 3/2 = √ (19)
2GM 2 GM
r0 −v

or, to match the wording of the Third Law,

4π 2 a3
T2 = .
GM
Thus we have Kepler’s Third Law: the square of the period of revolution of a
planet around the sun is proportional to the cube of the mean distance from
the sun.

11
3.1.5 Orbital Period
We now examine the orbital period of a satellite. We wish to express the time t in
terms of the angular displacement θ. We cover the major steps of this derivation
here, while the complete derivation is included in Appendix B. We begin with
the statement of conservation of angular momentum in polar coordinates using
the initial conditions. We have
dθ
r2 = r0 v0 sin α.
dt
We then use the equation for position r in terms of θ,

v02 r02 sin2 α/GM

r= q .
v0 r0 sin α 2GM G2 M 2
1+ GM v02 − r0 + v02 r02 sin2 α
cos(θ − θ0 )

Substituting for r and rearranging terms yields

v03 r03 sin3 αdθ

dt =  r 2 .
2 2 2
 
v0 r0 sin α 2GM
2 2
G M 1+ 1− G2 M 2 r0 − v02 cos(θ − θ0 )

We integrate from 0 to θ and take θ0 = 0, and α = 0, meaning the initial

position is at the perihelion. We have
θ
v03 r03 dθ
Z
t=  r 2 .
0
 
v02 r0 sin α
2 2
2GM
G2 M 2 1+ 1− G2 M 2 r0 − v02 cos(θ)

From our definition of eccentricity, e, we have

θ
v03 r03 dθ
Z
t= 2.
0 G2 M 2 (1 + e cos(θ))
(1−e2 )
In order to evaluate this integral, We begin by multiplying by (1−e2 ) and
our integral becomes
θ
v03 r03 (1 − e2 )dθ
Z
t= 2.
G M 2 (1 − e2 )
2
0 (1 + e cos(θ))
We now utilize the Weierstrass Substitution
 
θ
U = tan .
2
Using properties of trigonometric functions, we find

1 − U2 2U
cos θ = 2
and sin θ = .
1+U 1 + U2

12
We also find
2
dθ = dU.
1 + U2
We substitute these into our integral to get

v03 r03 2(1 − e2 )dU

Z
i2 .
G2 M 2 (1 − e2 )
h 
2
1 + e 1−U
1+U 2 (1 + U 2)

This integral is equivalent to

v03 r03 −2e(1 + e) + 2e(1 − e)U 2 v03 r03

Z Z
2
2 2 2 2 2
dU + 2 2 2
dU.
G M (1 − e ) [(1 + e) + (1 − e)U ] G M (1 − e ) (1 + e) + (1 − e)U 2

We integrate the first term by reversing the Weierstass substitution and

subsequently using integration by parts. The second term is integrated by using
another trigonometric substitution. We are then left with
θ
v03 r03
Z
dθ
t= 2
G2 M 2 0 (1 + e cos(θ))
" "r  # #
v03 r03 −e sin θ 2 1−e θ
= 2 2 +√ tan−1 tan + nT ,
G M (1 − e2 ) 1 + e cos θ 1 − e2 1+e 2
where T is the period of the orbit, given by

2πa3/2
T = √
GM
and n is the number of complete revolutions the planet has made, or for φ > 2π,
φ = θ + 2πn. Thus, the nT term accounts for the planet moving through more
than one complete orbit. As previously mentioned, the complete derivation of
this equation is presented in Appendix B.

3.2 The Cosmic Velocities

Now that we have proven Kepler’s Laws, we would like to explore some other
results governing orbits. We begin by finding an equation for the velocity of a
satellite at any position r. In our orbital system, we know that the total energy
E, given in the equation,
"   2 #
2
1 dr 2 dθ GM
+r − = E,
2 dt dt r

is constant by the principle of conservation of energy. The first term of this

equation is the kinetic energy of the system, and the second term is the potential
energy. We can now show that the type of orbit is dependent on the total energy
of the system. We begin by recognizing that the square of the velocity is equal

13
to the sum of the squares of the radial and angular components of the velocity.
We have  2  2
dr 2 dθ
+r = v2 ,
dt dt
or equivalently,
2GM
v2 = + 2E.
r
We then recall equation 15a
s
v 2 r2 sin2 α
 
2GM
e= 1 − 0 02 2 2
− v0 .
G M r0

Substituting our equation for velocity into the previous equation, we have
s
2v 2 r2 sin2 α
e = 1 + 0 02 2 E.
G M

We can now compare the three conditions, E = 0, E < 0, and E > 0. When
E = 0, we have e = 1, and thus the orbit is parabolic. The condition in
which E = 0 means that the kinetic energy is equal to the potential energy,
which in physical terms means that the object has the minimum amount of
energy necessary to escape orbit. Because every term in the coefficient of E is
squared, we know the coefficient is positive. Thus, for E < 0, we have e < 1
and thus the orbit is elliptic. The object does not have sufficient kinetic energy
to escape orbit. Similarly, for E > 0, we have e > 1 and therefore the orbit is
hyperbolic. Because the kinetic energy is greater than the potential energy, the
object escapes a closed orbit.

We now move on to orbital speed. We have previously established that

 2  2  
dr 2 dθ 1 1
+r = v02 + 2GM − .
dt dt r r0

This is equivalent to  
2 1 1
v = v02 + 2GM −
r r0
which gives velocity in terms of the initial conditions and the current position.
We have also previously shown that
GM
a= 2GM
.
r0 − v02

Rearranging terms yields,

2GM GM
v02 = − .
r0 a

14
Substituting this statement of our initial conditions into our equation for velocity
in terms of r yields
 
2 2GM GM 2GM 2GM 2 1
v = − + − = GM − ,
r0 a r r0 r a

or equivalently, s  
2 1
v= GM − .
r a
Thus we have an expression for the velocity of a satellite at any position r. We
can also see that if a = r, a circular orbit, that
GM
v2 = .
r
We now demonstrate that the initial propulsion speed v0 of a satellite launched
from the Earth’s surface has to fall within a certain range in order for the satel-
lite to remain in orbit, which is to say, neither hit the surface of the earth,
nor move out of the earth’s gravitational field. We know from the definition of
conic sections, that the size of e in the equation 14 determines the shape of the
orbit. We can see from equation 15b that the sign of quantity (v02 − 2GM/r0 )
determines whether e is less than, equal to, or greater than 1. We can then
define v0∗ as the velocity when e = 1, or
r s r
∗ 2GM 2gR2 Rp
v0 = = = 2gR, (20)
r0 r0 r0

where g = GM R2 . From equation 15a we can see that when the velocity is less
than the minimum initial velocity, v < v0∗ , then e < 1 and the object will either
fall back to Earth or orbit elliptically. If v = v0∗ , then e = 1 and the object
will follow a parabolic path, having just enough energy to escape the Earth’s
gravitational field. Finally, if v > v0∗ , then e > 1 and the orbit will be hyperbolic,
with the object also escaping the Earth’s gravitational
q √ field.
We now know that v0 must be less than rR0 2gR in order for a satellite to
stay within the gravitational field of the earth. We derive the minimum velocity
necessary for the satellite to stay in orbit. Thus, we are limited by rmin = r0 ,
the point at which the satellite collides with the surface of the earth. We know
rmin = pe/(1 + e) and, using equations 15a and 15b, we have

1 v02 r02 sin2 α

r0 = , (21a)
1+e GM
which is equivalent to
s
v02 r02 sin2 α v02 r02 sin2 α
 
2GM
−1= 1− − v02 . (21b)
GM G2 M 2 r0

15
The left hand side of equation 21b is always nonnegative, so we have
r r r
GM GM Rp
v0 ≥ csc α ≥ = gR. (22)
r0 r0 r0
Thus we have that the initial speed v0 of a projectile launched from the
surface of the earth must fall within the range,
r r
Rp Rp
gR ≤ v0 ≤ 2gR, r0 > R. (23)
r0 r0
Finally, if we square both sides of equation 21b, we find that sin α = 1 and
thus α = π/2. This means that at rmin , the angle between the position vector
and the velocity vector is π/2. Thus, the point where the satellite is closest to
the earth occurs along the major axis. If the satellite is launched with minimum
velocity from this point, the satellite will be further from the earth than rmin
for all other points along its orbit.
We now examine an alternative method for deriving the second cosmic ve-
locity, that is the maximum initial velocity an object can be launched from the
surface of a planet and remain in orbit. We previously established that for an
object launched from the planet’s surface to remain in orbit,
p p
gR < v0 < 2gR.

The alternative approach is to find the minimum velocity such that the object
will leave the planet and never return, which is to say that as r increases without
bound, the velocity remains nonnegative.
We know from F = ma, that

d2 y(t) GM
=− 2 .
dt2 r
If we define r = R + y(t) such that R is the radius of the planet and y(t) is
the altitude of the object above the surface of the planet, then our previous
equation becomes
d2 y(t) R2 g
2
=− . (24a)
dt [R + y(t)]2
We can now rearrange terms and solve for the velocity of the object in terms
of its height above the surface of the planet. We can rewrite equation 24a as

d[v(t)] R2 g
=− . (24b)
dt [R + y(t)]2

Rearranging terms and multiplying each side of equation 24b by v yields

vdt
vdv = −R2 g · . (25a)
[R + y(t)]2

16
We can now use the substitution U = R + y(t) and thus dU = vdt to get

dU
vdv = −R2 g · , (25b)
U2
or Z v Z y(t)−R
dU
vdv = −R2 g · . (25c)
v0 R U2
Evaluating this integral, we have

1 2 1 1 y(t)
(v − v02 ) = R2 g = R2 g  . (26a)

2 U R + y(t) 0

This simplifies to
 
2 1 1
v(t) = v02 + 2R g 2
− . (26b)
R + y(t) R

As previously stated, the minimum velocity such that the object never returns
to the planet occurs when the velocity reaches zero as y(t) → ∞. In equation
26b, if we let y(t) → ∞, we have
 
2 1
v∞ = v02 + 2R2 g 0 − ,
R

where v∞ is the velocity when y(t) is arbitrarily large. We simplify the previous
equation to
2
v∞ = v02 − 2Rg.
If v∞ > 0 at this point, the object will leave earth. Thus we have

0 < v02 − 2Rg,

or p
v0 > 2Rg.
√
Thus we have verified that the second cosmic velocity is v0 = 2Rg.

4 Einstein’s Theory of General Relativity

4.1 The Anomalous Advance of Mercury’s Perihelion
For nearly two hundred years after Newton’s Principia was published, the in-
verse square law was the accepted method by which to model the motions of the
planets and their satellites. Newton’s Law of Universal Gravitation was even
used by Urbain Jean Joseph Le Verrier to predict the existence of Neptune.
By the mid 19th century, however, problems began to arise, most notably with
predictions of the orbit of Mercury. As Le Verrier began his work on the orbit
of Mercury, the errors in the prediction of Mercury’s orbit were so large that

17
there was nearly an hour between the predicted and actual start of transit, the
time at which Mercury passes in front of the sun and so can be clearly seen and
an accurate position recorded. By 1845, Le Verrier had reduced this error to 16
seconds, which was still unacceptable in comparison to the predicted orbits of
other planets (Linton, 437).
Le Verrier discovered that the error in predictions of Mercury’s orbit was
a result of the advance of its perihelion, the closest point on the orbit to the
Sun. By analyzing the effect of the other planets on the perihelion advance,
Le Verrier discovered that Mercury’s perihelion advanced by 565” per century
rather than the 527” predicted by Newtonian mechanics. Despite this discovery,
Le Verrier was unable to explain the cause of the anomalous 38”. (Roseveare,
20-37).

Figure 4: Perihelion Advance (Cornell)

Though the problem of Mercury’s orbit was problematic for Newton’s Law of
Universal Gravitation, the scientific community was not yet willing to question
the dogma of Newtonian mechanics. Instead other causes of Mercury’s perihe-
lion advance were investigated. Based on the success of Le Verrier’s prediction
of the existence and position of Neptune, similar techniques were employed in
hopes of discovering a planet between the sun and mercury that could be effect-
ing the orbit. This hypothetical planet, deemed Vulcan, was “observed” many
times during the second half of the 19th century, though each discovery was
subsequently discredited (Linton, 441). Alternative hypotheses included that
the sun was oblate, as well as effects from the ether, the medium through which
light was thought to travel.
Though models of Mercury’s orbit were now mathematically accurate, deter-
mining the cause of the anomaly would require a shift in the way we conceived of
space and time. Shortly after the start of the 20th century, Einstein had already
published his theory of special relativity, which held that “the laws of physics
are of the same form in all inertial frames, and that in any given inertial frame,
the speed of light is the same whether the source is at rest or in uniform motion”
(Linton, 454). In response, the mathematician Hermann Minkowski translated
special relativity into a geometric framework, and as such created the concept
of four-dimensional spacetime. Within special relativity, which is to say in the

18
absence of a massive object, spacetime, as a four-dimensional Euclidean space,
can be thought of as ‘flat.’ Any object moving through this flat spacetime takes
the shortest path and thus travels in a straight line. In developing general rel-
ativity, Einstein postulated that mass had the effect of deforming Minkowskian
spacetime. An object moving through this deformed space time will still take
the shortest path possible, leading to the elliptical orbit of a planet (Linton,
463). The two-dimensional analog of this can be seen in Figure 3, mass deforms
flat space time and causes the trajectory of the moving mass to curve.

Figure 5: Mass deforms space time (Carroll)

Unfortunately the mathematical underpinnings of these concepts require a

knowledge of tensor calculus and are thus beyond the scope of this paper. Rather
we will begin with the general relativistic analog of the classical angular mo-
mentum and energy equation (previously defined as 11b) derived by simplifying
the spherically symmetric Swcharzchild solution to Einstein’s field equations.
After developing his gravitational theory, Einstein soon realized that it ex-
plained the anomalous advance of Mercury’s perihelion. In correcting the dif-
ference between the predicted and observed advance of Mercury’s perihelion,
Einstein brought to a close nearly 100 years of mathematical quarreling over
the cause of Mercury’s anomalous perihelion and changed our conception of the
universe (Roseveare, 147-186).

4.2 Relativistic Angular Momentum and Energy

In order to investigate the effect of General Relativity on the orbit of Mercury, we
must investigate the general-relativity analog of the previously stated classical
angular momentum and energy equation given by
 2  2  
dr dφ 1 1
+ r2 = v02 + 2GM − . (11b)
dt dt r r0

19
We let E be a constant related to the energy of the orbit and given by

v02 GM
E= − .
2 r0
Substituting this into Equation 11b yields
 2  2  2
dr dφ dφ 2GM
+ r2 = 2Eh2 + . (27a)
dφ dt dt r

We also have h, the angular momentum per unit mass given by

 
2 dφ
h=r .
dt

Dividing by h2 , we have
 2
1 dr 1 2E 2GM/r
+ = 2 +  2 . (27b)
r4 dφ r 2 h
r4 dφ
dt

2E
Because energy and angular momentum are conserved, h2 is a constant which,
for simplicity, we will call E0 and thus we have
 2
1 dr 1 2GM/r
+ = E0 +  2 . (27c)
r4 dφ r 2
r4 dφ
dt

We then introduce the substitution u = 1r , and thus du = − dr

r 2 . This yields
 2
du 2GM
+ u2 = E0 + u. (27d)
dφ h2

This is the classical angular momentum and energy equation. In solving Ein-
stein’s field equations, which is beyond the scope of this project, a general
relativistic analog of the previous equation is derived. We find that the gen-
eral relativistic equation is the classical equation with another term, equal to
2GM u3 /c2 . We expect that this term will change the orbit of a satellite from
pure Keplerian motion. We then have the equation
2
2GM u3

du 2GM
+ u2 = E0 + u + . (28a)
dφ h2 c2

The quantity 2GM/c2 is small compared to the radius of planetary orbits. We

can denote this quantity as ε, and we have
 2
du 2GM
+ u2 = E0 + u + εu3 . (28b)
dφ h2

20
Aphelion and Perihelion occur where du/dφ = 0, and thus we have
 
2GM
εu3 − u2 + u + E0 = 0. (29)
h2

This cubic equation then has three roots, let them be denoted u1 , u2 , and u3 .
Because ε is small, unless u is large, the roots will be close to the roots of the
analogous classical equation
 
2 2GM
u − u − E0 = 0.
h2

Thus, we say u1 and u2 are close to the roots of the Newtonian model and thus
we let u1 be the aphelion and u2 be the perihelion. Thus we know u1 ≤ u ≤ u2 .
We also know that u1 + u2 + u3 = 1/ε. For a proof of this see appendix C.
Because 1/ε is large, u3 is also large and thus has no physical meaning in our
model. Substituting our solution for equation 29 into equation 28a, we have
du 1
= [ε(u − u1 )(u2 − u)(u3 − u)] 2 . (30a)
dφ
Rearranging terms, we have
du 1
= [[(u − u1 )(u2 − u)(ε(u3 − u))] 2 .
dφ
We then rewrite this as
du 1
= [[(u − u1 )(u2 − u)(ε(u1 − u1 + u2 − u2 + u3 − u))] 2 .
dφ
We then group terms such that
du 1
= [[(u − u1 )(u2 − u)(ε(u1 + u2 + u3 ) − ε(u1 + u2 + u)))] 2 .
dφ
We can cancel ε and u1 + u2 + u3 . This yields
du 1
= [[(u − u1 )(u2 − u)(1 − ε(u1 + u2 + u)))] 2 .
dφ

In terms of dφ/du, we have

dφ 1
= [1 − ε(u1 + u2 + u)]−1/2 . (30b)
du [(u − u1 )(u2 − u)]1/2

Using the first order approximation from the Taylor expansion, (1 + x)k =
1 + kx, we have
dφ 1 + 12 ε(u1 + u2 + u)
≈ . (30c)
du [(u − u1 )(u2 − u)]1/2

21
 We now let α = 12 (u1 + u2 ) and β = 12 (u2 − u1 ). We then have

dφ 1 + 12 ε(2α + u)
= .
du [−u2 + (u1 + u2 )u − u1 u2 ]1/2
This is equivalent to
dφ 1 + 12 ε(2α + u)
= 1 2 .
du [ 4 (u2 − 2u1 u2 + u21 ) − u2 + (u1 + u2 )u − 14 (u21 + 2u1 u2 + u22 )]1/2
This reduces to
dφ 1 + 1 ε(2α + u)
= 2 2 . (31a)
du [β − (u − α)2 ]1/2
Integrating equation 31 with respect to u from u1 to u2 allows us to find the
angle between an aphelion and the subsequent perihelion. We have
1 + 12 ε(2α + u)
Z u2
∆φ = 2 2 1/2
du. (31b)
u1 [β − (u − α) ]

Rearranging terms, we have

u2 1 3
2 ε(u − α) + 1 + 2 εα
Z
∆φ = du.
u1 [β 2 − (u − α)2 ]1/2
We can split the integral into two terms such that,
1 + 32 εα
Z u2 Z u2
1 (u − α)
∆φ = ε 2 2 1/2
du + 2 2 1/2
du.
2 u1 [β − (u − α) ] u1 [β − (u − α) ]

The first term is relatively easy to integrate and we simply have

Z u2
1 (u − α) 1
ε 2 2 1/2
du = − ε(β 2 − (u − α)2 )1/2 |uu21 .
2 u1 [β − (u − α) ] 2
This yields
Z u2 2
(u2 − u1 )2

1 (u − α) 1 1
ε 2 2 1/2
du = − ε − (u2 − (u2 + u1 ))
2 u1 [β − (u − α) ] 2 4 2
2
(u2 − u1 )2

1 1
+ ε − (u1 − (u2 + u1 )) .
2 4 2

This simplifies to
u2
(u − α)
Z
1
ε du = 0.
2 u1 [β 2 − (u − α)2 ]1/2
For the second term, we have
u2
1 + 23 εα
Z u2  
3 u−α 3
2 2 1/2
du = (1 + εα) arcsin = (1 + εα)π.
u1 [β − (u − α) ] 2 β u1 2

22
Doubling ∆φ and subtracting 2π gives the angle between successive perihelions.
We then have that each orbit advances the perihelion by
 
3GM π 3GM π 1 1
φ = 3εαπ = (u 1 + u2 ) = + , (32)
c2 c2 r1 r2

where r1 and r2 are the values of r at aphelion and perihelion. In the case of
Mercury, r1 and r2 are small compared to orbital radii of the other planets.
As a result, Mercury has a larger perihelion advance then the other planets..
In addition, because the orbit of Mercury is more elliptical than the orbits of
most of the other planets, aphelion and perihelion are easier to observe. For
Mercury, the advancement of perihelion works out to about 43” per century.
While small, it is enough to be observed, and previous to Einstein’s Theory of
General Relativity, was unexplainable.

5 Conclusion
In this paper, we began by verifying Kepler’s Laws of Planetary Motion using
Newton’s Law of Universal Gravitation. We then derived several results per-
taining to the orbital period of a satellite, including solving for orbital period in
terms of angular displacement. We then derived an equation for the velocity of a
satellite, as well as the cosmic velocities, the minimum and maximum velocities
necessary for a satellite to stay in orbit. Finally, we examined the anomalous
orbit of Mercury using Newton’s Law of Universal Gravitation and Einstein’s
Theory of Relativity.
In our use of classical mechanics in exploring planetary motion, we assumed
several things. First, we assumed that we have a two-body system, using the
example of a satellite orbiting Earth. In reality, a massive body’s orbit is per-
turbed by the gravitational attraction to any other mass in its vicinity. Thus,
the orbits of the planets are not simply two-body systems with a given planet
orbiting the Sun, rather the orbit of each planet is affected by the other plan-
ets. Unfortunately, this multiple-body problem is very complicated and must
be modeled numerically. Fortunately, the Sun is vastly more massive than the
planets and thus the two-body model is a reasonable approximation. Our second
assumption was that the Earth is is a homogeneous sphere rather than a hetero-
geneous ellipsoid. The ellipsoidal shape of the Earth causes small perturbations
of the elliptic path of a satellite.
Additional considerations include that the perihelion advance of Mercury’s
orbit is not simply due to the effects of General Relativity. The gravitational
effects of other planets, as well as the fact that we are observing Mercury from
a moving platform, contribute to Mercury’s perihelion advance. In our analysis
we simply derived the additional orbital advance that was discovered by Le
Verrier and explained by General Relativity.
Our analysis also assumes no knowledge of tensor calculus in the deriva-
tion of Mercury’s perihelion advance from the relativistic angular momentum

23
and energy equation. Given an understanding of tensor calculus, the angular
momentum and energy equation can be derived from Einstein’s field equations.
Reasonable areas for further work would be deriving additional results per-
taining to the orbits of individual satellites as well as analysis of the three-body
problem. A useful resource for additional work in this area is Solar System
Dynamics by Carl D. Murray and Stanley F. Dermott. Given an understanding
of tensor calculus, the angular momentum and energy equation could also be
derived, as well as other results in General Relativity. A useful resource for addi-
tional work in this area is A Short Course in General Relativity as referenced in
the bibliography. Further information on the historical context of the mathemat-
ics covered in the paper can be found in From Eudoxus to Einstein: A History
of Mathematical Astronomy as referenced in the bibliography.

24
6 Appendices
6.1 Appendix A
On page 7, we used the standard process of completing the square to transform
equation 12b into equation 12c. The full derivation is given here. We begin
with equation 12b
v0 r0 sin αdz
−dθ = r  . (12b)
2 2GM 2 2 2 2
v 0 − r0 + 2GM z − v0 r0 sin αz

The denominator is
 
2 2GM
(v02 r02 2
sin α)z − (2GM )z + 2
− v0 .
r0

This is equivalent to
 2
GM
v0 r0 sin αz − + C,
v0 r0 sin α

where C is some constant that we will determine. We then expand the above
statement, yielding

G2 M 2
(v02 r02 sin2 α)z 2 − (2GM )z + + C.
v02 r02 sin2 α

Equating this statement with the square of the denominator of equation 12b
and simplifying we have

G2 M 2 2GM
2 2 2 +C = − v02 .
v0 r0 sin α r0

Solving for C, we have

2GM G2 M 2
C= − v02 − 2 2 2 .
r0 v0 r0 sin α

We then find a common denominator, which gives us

2GM v02 r0 sin2 α − v04 r02 sin2 α − G2 M 2

C= .
v02 r02 sin2 α

Thus equation 12b becomes

v0 r0 sin αdz
−dθ = r  2 . (12c)
(G2 M 2 −2GM v02 r0 sin2 α+v04 r02 sin2 α) GM
v 2 r 2 sin2 α
− v0 r0 sin αz − v0 r0 sin α
0 0

25
6.2 Appendix B
In section 2.1.5, we examined the orbital period of a satellite and found an
equation for the time t in terms of the angular displacement θ. The complete
derivation is presented here. We begin with the statement of conservation of
angular momentum in polar coordinates using the initial conditions. We have
dθ
r2 = r0 v0 sin α.
dt
We then use the equation for position r in terms of θ,

v02 r02 sin2 α/GM

r= q .
v0 r0 sin α 2GM G2 M 2
1+ GM v02 − r0 + v02 r02 sin2 α
cos(θ − θ0 )

Substituting for r yields

 
4
v04 r04 sin α/G M 2 2
 dθ
= r0 v0 sin α.

2 
dt
 q
v0 r0 sin α 2GM G2 M 2
1+ GM v02 − r0 + v02 r02 sin2 α
cos(θ − θ0 )

Rearranging terms, we have

v03 r03 sin3 αdθ

dt =  q 2 .
v0 r0 sin α 2GM G2 M 2
G2 M 2 1 + GM v02 − r0 + v02 r02 sin2 α
cos(θ − θ0 )

Again we rearrange terms, which yields

v03 r03 sin3 αdθ

dt =  q 2 .
2 2 v04 r02 sin2 α 2v02 r0 sin2 α
G M 1+ G2 M 2 − GM + 1 cos(θ − θ0 )

This is equivalent to

v03 r03 sin3 αdθ

dt =  r 2 .
2 2 2
 
v0 r0 sin α 2GM
G2 M 2 1 + 1 − G2 M 2 r0 − v02 cos(θ − θ0 )

For θ0 = 0, and α = 0, meaning the initial position is the perihelion, we have

Z θ
v03 r03 dθ
t=  r 2 .
0
 
2 2 v02 r02 sin2 α 2GM 2
G M 1 + 1 − G2 M 2 r0 − v0 cos(θ)

From our definition of eccentricity, e, we finally have

Z θ
v03 r03 dθ
t= 2.
2 2
0 G M (1 + e cos(θ))

26
 (1−e2 )
In order to evaluate this integral, We begin by multiplying by (1−e2 ) such
that our integral becomes
θ
v03 r03 (1 − e2 )dθ
Z
t= 2.
G2 M 2 (1 − e2 ) 0 (1 + e cos(θ))
We now utilize the Weierstrass Substitution such that
 
θ
U = tan .
2
We can use properties of trigonometric functions to find expressions for
cos(θ), sin(θ), and dθ. We begin with cos θ and we know that
1 − cos θ
U2 = .
1 + cos θ
We now rewrite this equation, and we have

U 2 + U 2 cos θ = 1 − cos θ
2
cos θ + U cos θ = 1 − U2
(1 + U 2 ) cos θ = 1 − U2
1 − U2
cos θ =
1 + U2
To find an expression for sin θ, we begin with

sin2 θ = 1 − cos2 θ
2
1 − U2

sin2 θ = 1−
1 + U2
1 + 2U 2 + U 2 1 − 2U 2 + U 2
sin2 θ = −
(1 + U 2 )2 (1 + U 2 )2
2
4U
sin2 θ =
(1 + U 2 )2

and we have
2U
sin θ = .
1 + U2
To find dθ, we can use the derivative of sin θ. We then have

(1 + U 2 ) · 2 − (2U )(2U )
cos θdθ = dU.
(1 + U 2 )2

We then have
(1 + U 2 ) · 2 − (2U )(2U ) 1 + U 2
dθ = dU.
(1 + U 2 )2 1 − U2

27
This simplifies to
2
dθ = dU.
1 + U2
We can now substitute into our integral. This yields

v03 r03 2(1 − e2 )dU

Z
i2 .
G2 M 2 (1 − e2 )
h 
2
1 + e 1−U
1+U 2 (1 + U 2)

Rearranging terms, we have

v03 r03 2(1 − e2 )dU

Z
i2 .
G2 M 2 (1 − e2 )
h
(1+U 2 )+e(1−U 2 ) 2)
1+U 2 (1 + U

This simplifies to

v03 r03 2(1 − e2 )(1 + U 2 )dU

Z
.
G2 M 2 (1 − e2 ) [(1 + e) + (1 − e)U 2 ]2

Expanding the numerator, we have

v03 r03 2(1 − e2 + U 2 − e2 U 2 )dU

Z
.
G2 M 2 (1 − e2 ) [(1 + e) + (1 − e)U 2 ]2

We can rewrite the numerator, yielding

v03 r03 2(1 − e2 ) + 2(1 − e2 )U 2

Z
2 2 2
dU.
G M (1 − e ) [(1 + e) + (1 − e)U 2 ]2

Factoring yields

v03 r03 2(1 + e)(1 − e) + 2(1 + e)(1 − e)U 2

Z
dU.
G2 M 2 (1 − e2 ) [(1 + e) + (1 − e)U 2 ]2

We now rewrite the numerator as

v03 r03 −2e(1 + e) + 2(1 + e) + 2e(1 − e)U 2 + 2(1 − e)U 2
Z
2 2 2
dU.
G M (1 − e ) [(1 + e) + (1 − e)U 2 ]2

We can now separate this into two terms

v03 r03 −2e(1 + e) + 2e(1 − e)U 2 v03 r03 2(1 + e) + 2(1 − e)U 2
Z Z
dU + dU.
G2 M 2 (1 − e2 ) [(1 + e) + (1 − e)U 2 ]2 G2 M 2 (1 − e2 ) [(1 + e) + (1 − e)U 2 ]2

The second term reduces, and we have

v03 r03 −2e(1 + e) + 2e(1 − e)U 2 v03 r03

Z Z
2
dU + 2 2 dU.
G2 M 2 (1 − e2 ) [(1 + e) + (1 − e)U 2 ]2 G M (1 − e2 ) (1 + e) + (1 − e)U 2

28
 We begin by integrating the first term,

v03 r03 −2e(1 + e) + 2e(1 − e)U 2

Z
2 2 2
dU.
G M (1 − e ) [(1 + e) + (1 − e)U 2 ]2

We rewrite this as
v03 r03 −2e[(1 + e) − (1 − e)U 2 ] 1 + U2
Z   
2
= dU.
G M 2 (1 − e2 )
2 [(1 + e) + (1 − e)U 2 ]2 2 1 + U2

Expanding the numerator, we have

v3 r3 −e[(1 + e) + (1 + e)U 2 − (1 − e)U 2 − (1 − e)U 4 ]

Z  
2
= 2 20 0 2 dU .
G M (1 − e ) [(1 + e) + (1 − e)U 2 ]2 1 + U2

This is equivalent to

v03 r03 −e[(1 + e) + U 2 + 2eU 2 − U 2 − (1 − e)U 4 ]

Z  
2
= 2 2 dU .
G M (1 − e2 ) [(1 + e) + (1 − e)U 2 ]2 1 + U2

Rearranging terms in the numerator yields

v03 r03 −e[(1 − U 4 ) + e − 2eU 2 + eU 4 + 4eU 2 ]

Z  
2
= 2 2 dU .
G M (1 − e2 ) [(1 + e) + (1 − e)U 2 ]2 1 + U2

We now factor terms in the numerator and have

v03 r03 −e(1 − U 2 )(1 + U 2 ) − e2 (1 − U 2 )(1 − U 2 ) − e(2U )2
Z  
2
= 2 2 dU .
G M (1 − e2 ) [(1 + e) + (1 − e)U 2 ]2 1 + U2
 2
1+U 2
We now multiply by 1+U 2 such that we have

−e(1−U 2 )(1+U 2 )−e2 (1−U 2 )(1−U 2 )−e(2U )2

v3 r3
Z  
(1+U 2 )2 2
= 2 20 0 2 i2 dU .
G M (1 − e ) 1 + U2
h  
1+U 2 1−U 2
1+U 2 + e 1+U 2

We then rearrange terms, which yields

   2  2
1−U 2 2

v03 r03
Z −e
1+U 2 − e2 1−U1+U 2
2U
− e2 1+U 2

2

= i2 dU .
G M 2 (1 − e2 )
2 1 + U2
h 
2
1 + e 1−U
1+U 2

We can combine the first and second terms of the numerator such that
 h  i  2
Z −e 1−U 22 1 + e 1−U 22 − e2 2U 2
v03 r03
 
1+U 1+U 1+U 2
= 2 2 2 dU .
G M (1 − e2 ) 1 + U2
h  i
2
1 + e 1−U
1+U 2

29
We can now substitute for cos θ, sin θ, dθ. This yields

v3 r3 (1 + e cos θ)(−e cos θ) − (−e sin θ)(−e sin θ)

Z
= 2 20 0 2 dθ.
G M (1 − e ) (1 + e cos θ)2
Splitting the integral into two terms, we have

v3 r3 e2 sin2 θ
 Z Z 
e cos θ
= 2 20 0 2 − dθ − dθ .
G M (1 − e ) 1 + e cos θ (1 + e cos θ)2
We begin integration by parts on the second term such that

U = e sin θ dU = e cos θdθ

1 e sin θ
V = 1+e cos θ dV = (1+e cos θ)2 dθ

We now have
v03 r03
 Z Z 
e cos θ e sin θ e cos θ
= 2 2 − dθ − + dθ .
G M (1 − e2 ) 1 + e cos θ 1 + e cos θ (1 + e cos θ)
The first and third terms cancel and we are left with
v03 r03 −e sin θ
= .
G2 M 2 (1 − e2 ) 1 + e cos θ

We now integrate the second term,

v03 r03
Z
2
2 2 2
dU.
G M (1 − e ) (1 + e) + (1 − e)U 2
√
We multiply this by √1−e and 1+e
such that we have
1−e 1+e
√
v3 r3
Z 2 √1−e √ 1+e
√
1−e 1+e 1+e
= 2 20 0 2 dU.
G M (1 − e ) (1 + e) + (1 − e)U 2
Rearranging terms yields
√ 
√ 2√ √1−e
v03 r03
Z
1−e 1+e 1+e
= h i dU.
G M 2 (1 − e2 )
2 (1+e)+(1−e)U 2
1+e

This is equivalent to
√ 
√1−e
v03 r03
Z
2 1+e
= √ √ √ 2 dU.
G M 2 (1 − e2 )
2
1−e 1+e √1−e U
1+ 1+e

30
We now introduce another trigonometric substitution such that
r r
1+e 1+e
U= tan(V ) dU = sec2 (V )dV.
1−e 1−e
Substituting, we have
v3 r3 sec2 (V )
Z
2
= 2 20 0 2 √ √ dV.
G M (1 − e ) 1 − e 1 + e 1 + tan2 (V )
We know sec2 (V ) = 1 + tan2 (V ), and so we reduce and integrate yielding
v3 r3 v3 r3
Z
2 2
= 2 20 0 2 √ √ dV = 2 20 0 2 √ √ · V.
G M (1 − e ) 1 − e 1 + e G M (1 − e ) 1 − e 1 + e
Substituting for V yields
r !
v3 r3 2 1−e
= 2 20 0 2 √ tan−1 U .
G M (1 − e ) 1 − e2 1+e
Substituting for U , we have
"r  #
v3 r3 2 1−e θ
= 2 20 0 2 √ tan−1 tan .
G M (1 − e ) 1 − e2 1+e 2

Finally, we have
θ
v03 r03
Z
dθ
t= 2
G2 M 2 0 (1 + e cos(θ))
" "r  # #
v3 r3 −e sin θ 2 1−e θ
= 2 20 0 2 +√ tan−1 tan + nT ,
G M (1 − e ) 1 + e cos θ 1 − e2 1+e 2
where T is the period of the orbit, given by
2πa3/2
T = √
GM
and n is the number of complete revolutions the planet has made, or for φ > 2π,
φ = θ + 2πn.

6.3 Appendix C
Proving u1 + u2 + u3 = 1ε .
When we expand ε(u − u1 )(u2 − u)(u3 − u), the u2 term works out to be −(u1 +
u2 + u3 )u2 . Equating this with the u2 term in our energy equation, we have
−(u1 + u2 + u3 )u2 = −u2 , and thus
1
u1 + u2 + u3 = .
ε

31
7 Bibliography
Carroll, Sean M. “Spacetime and Geometry: An Introduction to General
Relativity.” Preposterous Universe. California Institute of Technology.
5 May 2009. http://preposterousuniverse.com/spacetimeandgeometry/

Foster, James and J.D. Nightingale. A Short Course in General Relativity. 1st.
London, Longman Mathematical Texts, 1979.

Kwok, Yue-Kuen. “Motion of an Artificial Satellite About the Earth.”

Undergraduate Mathematics and its Applications Project
695(1989): 130-147.

Linton, C.M. From Eudoxus to Einstein: A History of Mathematical Astronomy.

1st. Cambridge, UK: Cambridge UP, 2004.

Roseveare, N.T. Mercury’s Perihelion: From Le Verrier to Einstein. 1st.

Oxford, UK: Oxford UP, 1982.

“The Advance of the Perihelion of Mercury.” Cornell University Dept. of

Astronomy. 5 May 2009.
http://www.astro.cornell.edu/academics/courses/astro201/merc adv.htm

32