The Restricted Three-Body Problem: Sun-Jupiter & "X"

AP1601
Columbia University

Introduction

Some definitions for Jupiter's orbit. For simplicity, we'll take Jupiter's orbit to be a circle, and we'll study only motion on the plane.

aJ = 5.203; (* Jupiter's Radius in AU *)
ωJ = 2π/aJ^(3/2);(* Jupiter's rotation in AU/year *)
massJS= 1/1048.;

rJ = {aJ Cos[ωJ t],aJ Sin[ωJ t],0}

Equations of Motion: Sun & Jupiter

The equations of motion for a planet describe the time evolution of the planet's phase-space location and its initial position.

p = {x[t], y[t], vx[t], vy[t]};

fJ = - 4 (π^2)mJ  dr/(dr . dr)^(3/2) /. dr -> rJ - {x[t],y[t],0}; (* Force from Jupiter *)

dpdt = {D[x[t],t] == vx[t],
D[y[t],t] == vy[t],
D[vx[t],t] == -4 (π^2) x[t]/r^3 (* Sun *) + fJ[[1]],
D[vy[t],t] == -4 (π^2) y[t]/r^3 (* Sun *) + fJ[[2]]} /. r -> Sqrt[x[t]^2 + y[t]^2]

p0 = {x[0] == x0, y[0] == y0,vx[0] == vx0, vy[0] == vy0}

eqs = dpdt ~ Join ~p0

Earth's Orbit

In astronomical units (AUand year), the trajectory of the Earth is relatively easy to calculate. The key is finding the appropriate initial condition. 
Let's set the position of Earth to be {x, y} = {1, 0}. The initial velocity is approximately the circumference, 2πr = 2π, dividedby 1 year, or {vx, vy} = {0, 2π}.

earth = NDSolve[eqs /.x0 -> 1 /. y0 -> 0 /. vy0->2 π/. vx0 -> 0/. mJ -> massJS, p, {t, 0, 1.2}]

ParametricPlot[Evaluate[{x[t],y[t]}/.earth],{t,0,1},AspectRatio->Automatic]

The orbit is a circle because of our choice of the initial condition.

Longer Integration Times

earth = NDSolve[eqs /.x0 -> 1 /. y0 -> 0 /. vy0->2 π/. vx0 -> 0/. mJ -> massJS, p, {t, 0, 200}, MaxSteps -> 20000]

Plot[Evaluate[(Sqrt[((x[t])^2)+((y[t])^2)])-1./.earth],{t,0,50},PlotLabel->"Jovian-Earth Resonance"]

Asteroid Resonances

The asteroids lie within a range from 2 < a < 3.5. Orbits within this range have low-order resonances with Jupiter. The ratio of the Earth's orbital frequency to Jupiter's is ωE/ωJ ≈ 11.8681, and this is not near a low-order resonance.

Kirkwood Gaps

Kirkwood gaps shown in a histogram of asteroid orbits in the main asteroid belts. The gaps labeled "3:1", "5:2", "7:3", "2:1" are caused by mean-motion resonances between an asteroid and Jupiter. For example, at the 3:1 resonance, the asteroid completes three orbits for every one orbit of Jupiter.

Sun-Jupiter Gravitational Resonances

ω[a_] := 2 π /a^(3/2)

Plot[((ω[a])/(ω[aJ])),{a,2,5},GridLines->Automatic,PlotRange->{0,4},PlotLabel->"Orbital Rotation Ratios"]

FindRoot[ω[a]/ω[aJ]==3, {a, 2.5}]

FindRoot[ω[a]/ω[aJ]==5/2, {a, 3.0}]

FindRoot[ω[a]/ω[aJ]==7/3, {a, 3.0}]

FindRoot[ω[a]/ω[aJ]==2, {a, 3.2}]

Let's examine some orbits near some low-order resonances…

Example Asteroid Orbit

aStart = 3.278;

asteroid = First[NDSolve[eqs /.x0 -> aA /. y0 -> 0 /. vy0->aA ω[aA]/. vx0 -> 0/. mJ -> massJS /. aA -> aStart, p, {t, 0, 1000}, MaxSteps -> 20000]]

radius[t_] = (Sqrt[((x[t])^2)+((y[t])^2)])/.asteroid

Plot[radius[t]- aStart,{t,0,1000},PlotLabel->"Jovian-Asteroid Resonance"]

p1=ParametricPlot[Evaluate[{x[t],y[t]}/.asteroid],{t,0,10},AspectRatio->Automatic,DisplayFunction->Identity];
p2=ParametricPlot[Evaluate[{x[t],y[t]}/.asteroid],{t,150,160},AspectRatio->Automatic,PlotStyle->Red,DisplayFunction->Identity];
p3=ParametricPlot[Evaluate[{x[t],y[t]}/.asteroid],{t,380,390},AspectRatio->Automatic,PlotStyle->Blue,DisplayFunction->Identity];
Show[p1,p2,p3,DisplayFunction->$DisplayFunction]

radList = Table[radius[t], {t, 0, 1000, 0.5}];

Histogram[radList,{2.0,3.8,0.05}, PlotLabel -> "Asteroid Start Radius = "<>ToString[aStart], PlotRange -> {0,500}]

Histograms

asteroidLoc[as_] := Module[{sol,radList},
sol = First[NDSolve[eqs /.x0 -> aA /. y0 -> 0 /. vy0->aA ω[aA]/. vx0 -> 0/. mJ -> massJS /. aA -> as, p, {t, 0, 1000}, MaxSteps -> 30000]];
radList = Table[(Sqrt[((x[t])^2)+((y[t])^2)])/.sol, {t, 0, 1000, 0.5}];
Histogram[radList,{2.0,3.8,0.05}, PlotLabel -> "Start Radius = "<>ToString[as]<>", (Subscript[ω, a])/(Subscript[ω, J]) = "<>ToString[Rationalize[ω[as]/ω[aJ],0.1],InputForm], PlotRange -> {0,2000}]]

asteroidLoc[2.5]
asteroidLoc[3.2]

Manipulate[asteroidLoc[a], {a, 2.1, 3.4}]

Summary

The "restricted" 3-body problem is an relatively simple example of complexity in many-body dynamics. They illustrate the important concept of resonance in dynamics. When an oscillating object (like an asteroid) interacts with another oscillating object (like Jupiter), large orbital variations are possible. The orbital dynamics of the "restricted 3-body problem" can become complicated, but it is not chaotic.

Of course, our calculations are idealized. They do not take into account the other planets, ellipticity of Jupiter, and gravitational collisions among asteroids. The worst ommission is time. We only integrated orbits for a few hundred years!