COS 126

N-body Simulation
Programming Assignment 8

Due: Wednesday, 11:59pm

Write a program that simulates the motion of N bodies, mutually affected by gravitational forces, in a two-dimensional space.

This program will be a Java applet similar to the one you did for the ``warmup''. You may use the code /u/cs126/code/ as a starting point. Your task is to implement the methods move in class Body (which updates the position and velocity, using the current force) and moveall (which recalculates forces and uses move to update positions and velocities) and showall (which uses show to display the bodies) in class Nbody.

The method initall initializes a 3-body system for you to use in debugging. Each body is characterized by 6 variables: the x and y position coordinates; initial values for vx and vy (the x and y components of the velocity); the ``mass'' of the body; and the size (in pixels) to use in displaying it. The position and velocity are updated when we compute the effects of mutual gravitational interactions, as we will see.

If you have not done the warmup, do it; otherwise you will find it very difficult to complete this assignment. First, get showall implemented, and make sure that you understand the basic data structures and how main works. This is the easy part: don't even put in any ``erase'' code at this point; it's liable to get in the way of debugging. Note that initall sets a variable max that is the largest value of any x or y coordinate.

Next, do the code for updating the coordinates for each display. This involves a ``force'' calculation and a ``velocity and acceleration'' calculation, both from elementary physics. We consider these in the next two paragraphs.

The strength of the gravitational force between two bodies is given by the product of their masses divided by the square of the distance between them, scaled by a constant G. You need to calculate these forces for all possible pairs of bodies and add all the forces acting on each body. This is a vector addition: the pull of one body towards another is on the line between them. Therefore you will have to compute the x and y components of the forces and add them, to get Fx and Fy force components representing the total gravitational pull on each body of all the other bodies.

Next, use the basic formula F[] = ma[] to get the x and y components of the acceleration of each body that is implied by these forces over the given time quantum dt. From these, the velocity and position can be updated: the change in velocity is given by v[] = a[] dt and the change in position is given by p[] = v[] dt + a[]*dt*dt/2. The square brackets are intended to indicate that these are all vector formulae: for example p[] represents the x and y position coordinates, v[] represents the x and y coordinates of the velocity, and so forth.

Part of the challenge of this assignment is in coping with realistic scientific data values, as you will see when examining the data files. Use double floating point values at every opportunity, though this isn't even quite good enough (see below). The force, velocity, and acceleration update function only takes thirty or so lines of straightforward code, but the use of ``real data'' means that you will want to think carefully about what the code is doing, as debugging will be difficult.

To speed things up somewhat, you should fix things so that moveall does T steps, where T is set to, say, 100. Proper setting of this value will make the planets move smoothly, and may depend upon the machine that you are using.

To make debugging easier, borrow from your experience on the warmup to keep the bodies in the window by having them bounce off the edges. Though it is extremely dubious that this happens in the real-world situations, this makes for a ``virtual'' universe with fascinating properties.

Once you get the orbits plotting in ``track'' mode, put the erase code in, as you did for the warmup. As in the warmup, the display looks particularly good if you use a light color (or gray) for the background, white for the path traces, and dark colors for the bodies.

After you get the orbits for the bodies given in the initall procedure in working, try the system in the file SSfive.

Extra Credit: Add a slow moving comet to SSfive. Try to make it hit one of the planets.

Challenges for the bored: Most solutions to this problem will have serious roundoff errors in the computation. For example, this leads to all planets eventually spiraling into the sun in SSfive. Find a way to fix this. Make the simulation three-dimensional by doing calculations for x, y, and z coordinates, then using the z coordinate to vary the sizes of the planets.

Copyright © 1998 Robert Sedgewick