1) In this problem we are going to compute the orbits in the Kepler and Navarro, Frenk, and White
(NFW) potentials. We will begin with a case in which we already have information analytically,
the Kepler potential with 0 = -GM/r. We will be solving the equation for the evolution of the
azimuthal coordinate 0 for a star with angular momentum L moving within this potential
= and the equation for the evolution of the radial coordinate of this star
where E is the total energy, G is the gravitational constant, and M is
the mass at the center of the potential.
a. Before we begin programming, our first step will be to nondimensionalize these
equations, by rescaling each of the units. To do this we'll define the following four
dimensionless variables i and = E is
an arbitrary scale radius. Replace r, t, L,
= these variables and write the new versions of L
and SEE = G and ro will not appear in your new equations.)
b. Now we will study orbits that have their apocenter at r = 1. Recall that this means that at
this value the radial orbit is just about to move inward, and ar
show that you rewrite your nondimensional form of eq. ar =
part (a) as =
To solve this equation numerically, along with the nondimension form
we will be integrating forward in time, choosing very small-time steps. The values of F and 0
during these steps will then trace out the motion of the star along the orbit. This means that
we start out with 00 and ro=(1-e) where € is a small number, say 10-6. Here the
subscripts indicate that these are the values at the 0th time step. We then use the equations
of motion to advance $ and r to the next time step as
Pi+1 = Qi + AT
Ti+1 = ri + AT = 10-3
Carry this out for orbits with L = 1/4, 1/2, 1/1.5, 1/1.2, 1/1.1 and 1/1.05 and plot them all
on a single plot with on the x-axis and r on the y-axis. The code will have problems as the
orbits move within the pericenter, as the argument of the square-root will go negative
there. Don't worry about this, and just plot up to the time that you reach this point. Add a
vertical line on the plot at $ = TT, and show that all the orbits reach this value at pericenter.
This means your code is working. Perhaps you may have small deviations, which are largest
for the orbits that pass closest to the center. Why do you think these are the ones that
might give the most trouble? On a new plot set x = r cos 0 and y = sin and plot all your
orbits in real space with X going from -1 to 1 and y going from 0 to 1. Which are the most
eccentric orbits? Those with the lowest or the greatest value of L?
C. Now consider the NFW halo. In this case the equation for is again
using equation where = 3.24, the radial equation
= + To nondimensionalize these equations, we'll define
the dimensionlessvariables= nd E=
Replace r, t, L,
E in equation with these
variables and write the newversionsot=*=and =
L² (Hint: an and IN will not appear in your new equations).
Note that now r is not just r scaled by an arbitrary radius, but instead r = 1 has a special
meaning. It is the r value at which the NFW density distribution starts to flatten out.
d. Now you'll study orbits that have their apocenter at r = 1. Show that in this case we can
rewrite the nondimensional
form of that you
- Solve this equation numerically,
exactly as you did above. Now make a plot for orbits with L = Lc
where LC = /In(2) - 1/2 is the angular momentum for a circular orbit at this radius. Again,
plot these all on a single figure with 0 on the x-axis and r on the y-axis. Add vertical lines on
the plot at = TT, the answer for the Kepler potential, and = T/2, the answer for the
constant density distribution, and show that all the orbits lie between these values. On a
new plot set x = r cos 0 and y = r sin and plot your orbits in real space with X going from -1
to 1 and y going from 0 to 1. Which are the most "eccentric" orbits?
e. Finally, consider NFW orbits with their apocenter at r = 2. Rewrite eq.
= - with the correct values of Ê for this distribution.
Then make D-r, and x-y plots for for orbits with Lc
2/1 In(3)/2 - 1/3. Don't forget that because we are starting further out, the total time to
reach pericenter is longer than in part (d). For the same value of L/Lc, do the orbits get
further in o in one radial orbit at this distance or at r = 1? Why do you think this is the
These solutions may offer step-by-step problem-solving explanations or good writing examples that include modern styles of formatting and construction
of bibliographies out of text citations and references. Students may use these solutions for personal skill-building and practice.
Unethical use is strictly forbidden.