Generalize the integration package you wrote for PS5 so that it can
solve the -body problem in 3-D for particles of arbitrary mass
using the PP method with fixed timesteps. The code should read the
masses and initial positions and velocities of all the particles from
a text file (use the following input format, one particle per line:
). Use
`#define`s or command-line arguments to specify the integrator
(either leapfrog or Runge-Kutta), the softening parameter, the
timestep, the number of steps, the output frequency, and the name of
the initial conditions file. The output format can be the same as the
input format, with one file per output (use `sprintf()` to
imbed the current step number to create unique output
filenames)--this has the advantage that you can restart from any
output step. However, for problem 1 you may want to have the option
of appending outputs to a single file, with one line per output step
and each line listing data for all particles in order.

- (Worth 85% of the grade) Check your code by solving the
two-body problem for equal masses () without
softening: plot phase diagrams of
*vs*. , where is the separation and is the radial component of the relative velocity, for 100 loops around orbits of eccentricity and , using both leapfrog and Runge-Kutta. Start at unit apoapse: (you will need to derive the starting for yourself). By symmetry your initial conditions are then , , , , where we have chosen to start on the -axis with the orbital angular momentum aligned with the positive -axis. Use a timestep of 0.05 for and 0.003 for (plot approximately 1000 points for each case). Also plot the energy*vs*. time and comment.You may find the following relations helpful for this problem:

Here is the semi-major axis and is the total mass. Assume . You will also need Newton's version of Kepler's Third Law to get the orbital period.- BONUS: extra 2% Compare your results with the analytical solution.
You will need to solve Kepler's Equation.

- BONUS: extra 2% Compare your results with the analytical solution.
You will need to solve Kepler's Equation.
- (Worth 15% of the grade) Generate random initial conditions of
your choosing for between 100 and 1000 particles, e.g. particles in
a spherical region, particles in a disk, whatever you want. Use a
range of masses. Integrate this system long enough for something
``interesting'' to happen. What should the timestep be? Do you
need softening? Is energy conserved? Does the code scale with
as you expect? Summarize the evolution of your system in a plot of
some kind to hand in.
- Make a movie of the results.