1997-98
NEW MEXICO
HIGH SCHOOL
SUPERCOMPUTING
CHALLENGE


Team Number: 89
School Name: Portales High School
Area of Science: Earth and Space Science
Project Title: Planetary Dynamics
Project Abstract: http://mode.lanl.k12.nm.us/97.98/abstracts/089.html
Interim Report: http://mode.lanl.k12.nm.us/97.98/interims/089.html
Final Report: http://mode.lanl.k12.nm.us/97.98/finalreports/089/finalreport.html
Program Code: http://mode.lanl.k12.nm.us/~ch089lxw/project.c
Note: There are images in this HTML document and would require the reader of this
document to have access to the internet to fully understand what is described in the report.

Please click on one of the following to go to that section of our Interim Report
I. Abstract
II. Newtonian Mechanics
III. Particle Dynamics Simulation
IV. Super Computer Implementation
V. Computer Graphical Animation
VI. Progress and Result to Date
VII. Apendix (Our Code)

I. Abstract

This year's project for the Portales challenge team will be based upon the theories of particle dynamics as applied to calculations of the dynamic behavior of real or hypothetical solar systems comprised of planets, a sun (or multiple orbiting suns), asteroids, comets, moons, and spaceships. The "particles" we will be using may include any object that might affect or be affected by the other objects in the given solar system. In our programming, a user will specify initial conditions for each particle (planets, sun, asteroid, etc.) such as particle mass starting positions (xi, yi, zi), and initial velocities (Vxi, Vyi, Vzi) in a given coordinate system. Other starting conditions could be mass and initial velocity of particles. The program will calculate these data in the context of Newton's Universal Law of Gravitation to find the inter-particle forces at word and then calculate new coordinates and new velocities of each particle for a given time interval(DeltaT). The new time interval will represent a new movie frame in our particle dynamics "motion picture". The motion picture of solar system particle dynamics will be produced via computer graphics animation by hooking together many frames representing many successive time increments DeltaT. This way we will find out how the objects within a given solar system act and interact and how changes (both real and hypothetical) can act to perturb the solar system.

II. Newtonian Mechanics

Newton’s kinematic equations of motion govern the motions of large bodies at low speeds and modest accelerations. These equations are listed below in vector form and in vector component form.



These equations give the position and the velocity of particle i in terms of initial position, initial velocity, and acceleration. These equations will be adapted to form the basis for a Newtonian particle dynamics simulation.

Newton’s Universal Law of Gravitation which expresses the potential energy of gravitational interaction between particles i and j with masses mi and mj is given by:



where the constant G is Newton’s gravitational constant and the inter-particle distance rij is expressed as:



For a multi-particle system, the gravitational potential energy affecting particle i is:



The force that particle I experiences as a result of interactions with the all other particles is given as the negative gradient of the gravitational potential energy:



The individual components of the force vectors are:

III. Particle Dynamics Simulation


The particle dynamics simulation is carried out by calculating Newton’s kinematic equations of motion iteratively through n time-steps, each of duration Dt. At each time-step the acceleration of particle i is defined as the force on particle i divided by its mass .



Particle positions at each time step are noted with a time-step superscript. These particle positions are stored in an array for subsequent plotting and animation programs. The forces and positions at each time-step are computed as:



IV. Supercomputer Implementation


The program we have constructed is in the C++ language and consists of 3 main loops. User inputs into the program include: all particle masses, initial positions, and initial velocities. In its present form, the program can handle particles that are either suns and/or planets. In the future, our program will have the provision to handle comets, asteroids, moons, and spaceship launches. The first loop is the time-step loop which goes from 1 to n time steps of duration Dt, n being the total number of time-steps in the entire program: for (int t=1;t <= n;t++) (n is a user input). Within the time-step loop are the particle i and j loops. In the particle loop i we initialize the forces Fi(x, y, z). The third and final loop is within loop i and is the particle loop j. In this loop are the main equations. Here we compute all the forces, velocities, and positions of each particle but only if i does not equal j. All the equations which you have seen in section II. (Newtonian Mechanics) are used in this loop to figure out what will happen to each planet. Overall in this code there are 10 equations to compute for each time one body interacts with another body. Apply that to 10 bodies: ; then multiplying that number by n time-steps, there are literally millions of calculations which would crash a normal PC.

V. Computer Graphical Animation


The visual presentation of our project will be shown to the judges by a program called TechPlot. The preliminary plotting of information will be done in Microsoft Excel. We will run our program which will output all the information into a file-format which TechPlot or Excel could recognize. Animation effects will be achieved by having the graphics package step through the sequence of particle positions as calculated in each particle dynamics time-step.

VI. Progress and Result to Date


As of this date (1/8/98) all particle dynamics equations appropriate for this project have been coded in C++. The code compiles and is now being tested for computational and scientific accuracy. A data set of planetary masses, positions, and orbital velocities has been obtained for our solar system:



We now lack only minor C++ syntax required for our code.

VII. Appendixes (code)

//This code is produced by Liangfei Wang and Dr. John Kenney, III
//Contributions have been made to this code by Gina Weber (lanl) and
// Ron Obenhaus (ENMU).

//Begin Program!

//This file defines the parameters for all C++ commands
#include < iostream.h >

//This file defines the mathematical functions in C++
#include < math.h >

//This file defines the function to output information into files
#include < fstream.h >

main()
{
//Create and open a file named stuff.txt
ofstream newfile ("stuff.txt");
//Initializing number of time intervals, DT, suns, and planets
int t;
float delt;
int i, j;

//Input information
cout << "How many time intervals?" << endl;
cin >> t;
cout << "How many seconds per time interval?" << endl;
cin >> delt;
cout << "How many planets in solar system?" << endl;
cin >> i;
cout << "How many suns in solar system?" << endl;
cin >> j;

//Initializing mass, starting positions, and starting velocities
float mass[i+j];
float xi[i+j][t], yi[i+j][t], zi[i+j][t];
float vix[i+j], viy[i+j], viz[i+j];

//Integer c is the sum of suns and planets
int c = i + j;

//This loop is for entering data for each mass, starting position,
//and starting velocities
for (int a = 1; a <= c; a++)
{
cout << "Enter mass for a planetary body (Sun or Planet) in mass * 10^23" << endl;
cin >> mass[a];
cout << "Enter starting positions (x, y, z):" << endl;
cin >> xi[a][1] >> yi[a][1] >> zi[a][1];
cout << "Enter starting velocities in m/sec (Vx, Vy, Vz)" << endl;
cin >> vix[a] >> viy[a] >> viz[a];
}

//First of the triple loop is the time interval loop
for (int time=1; time <= t; time++)
{
//Initializing forces to 0
//G is Newton's universal constant of gravitation
float rij;
float fix = 0;
float fiy = 0;
float fiz = 0;
double G = (0.00000000006772);

//Second loop is for k=1 to c particles
for (int k = 1; k <= c; k++)
{
//Third loop is for l=1 to c particles
for (int l = 1; l <= c; l++)
{
//Only run equations if l does not equal k
if (l != k)
{

//Compute interparticle positions and force vectors for
//particle i at time step(t-1)
rij = sqrt(pow(xi[k] - xi[l],2) + pow(yi[k] - yi[l],2) + pow(zi[k] - zi[l],2));

fix = (-(G * mass[k] * mass[l]) / sqrt(pow(rij,5))) * (xi[k][time] - xi[l][time]) + fix;
fiy = (-(G * mass[k] * mass[l]) / sqrt(pow(rij,5))) * (yi[k][time] - yi[l][time]) + fiy;
fiz = (-(G * mass[k] * mass[l]) / sqrt(pow(rij,5))) * (zi[k][time] - zi[l][time]) + fiz;

vix[k] = (delt * (fix / mass[k])) + vix[k];
viy[k] = (delt * (fiy / mass[k])) + viy[k];
viz[k] = (delt * (fiz / mass[k])) + viz[k];

xi[k][time] = xi[k][time] + (delt * vix[k]) + (1 / 2) * pow(delt,2) * (fix / mass[k]);
yi[k][time] = xi[k][time] + (delt * viy[k]) + (1 / 2) * pow(delt,2) * (fiy / mass[k]);
zi[k][time] = xi[k][time] + (delt * viz[k]) + (1 / 2) * pow(delt,2) * (fiz / mass[k]);

//Output planetary positions to the screen
cout << xi[k][time] << ", " << yi[k][time] << ", " << zi[k][time] << endl;

//Output planetary positions to the file: stuff.txt
newfile << xi[k][time] << " is the x coordinate on the map\n";
newfile << yi[k][time] << " is the y coordinate on the map\n";
newfile << zi[k][time] << " is the z coordinate on the map\n";
newfile << "\n";

//close file: stuff.txt
newfile.close();
}
}
}
}
return 0;
}
//End Program!

Team Members:

Sponsoring Teacher(s):

Project Advisor(s):

New Mexico High School Supercomputing Challenge