|
The argument can be made that object-oriented language provide the most efficient way of modeling the real world. Todd states his case, as he uses C++ to build a simulation system that mimics the force of one planet on another.
September 01, 1989
URL:http://www.drdobbs.com/cpp/force-based-simulations/184408204
Copyright © 1989, Dr. Dobb's Journal
Copyright © 1989, Dr. Dobb's Journal
Todd is a programmer/analyst with the Institute of Geophysics and Planetary Physics at UCLA. He is also associated with the NASA/JPL Planetary Data Systems project. Todd can be reached at 1104 N Orchard, Burbank, CA 91506.
You've probably heard it said that the languages that best model the real world are object-oriented. When you try to build a simulation system you realize just how true this is. This article presents a force-based simulation system written in Zortech's C++ that demonstrates how and why an object-oriented language like C++ is a natural choice for implementing simulation systems.
The most common type of simulation systems are probably circuit type simulations, where components have inputs and perform transformation of the input into some form of output. The output is then fed to another component. With force-based systems, you have a collection of bodies that forms a universe. Each body exerts a force on every other body in the universe. This means that just the presence of a body has an influence on all other bodies. This is significantly different from circuit-based simulation, where influences are channeled from component to component.
Probably the most difficult and critical element in a simulation system is managing the objects that exist in the simulation. The function of the object manager is to maintain the communication between bodies and to enforce the physical laws. Listing One contains the source code for the administrator for the system presented in this article. In this listing, the administrator method is called big_bang( ), of the object class UNIVERSE.
It is crucial to maintain the integrity of each body. In C++, this is best accomplished by making the physical parameters private and providing methods to return the physical parameters. Because the physical parameters for any body in the system are the same, the class for bodies is called BODY (also in Listing One). The method UNIVERSE::service( ) registers a body with the UNIVERSE. The collection of bodies to be serviced is maintained as an array of objects of class BODY, and is repeatedly stepped through when you call UNIVERSE::big_bang( ).
In the class of BODY, the most notable method is apply_ force( ), which takes an outside influence on the body and converts it into a change in physical parameters. The physical law that determines how bodies influence each other is encoded in this method. In the present system, we are using the derivative of the force law of gravity, which results in an instantaneous velocity vector. This velocity vector is then applied to the body's current velocity vector to determine its new vector. The other member functions of the BODY class are simple and I'll refer you to the source for details.
Viewing the results of the simulations can be done in a variety of ways. In some cases, a table of numbers is sufficient, but in most cases the results are best seen in some animated form. Listing Two contains those methods and functions that are directly related to the display of a simulation. The method BODY::update( ) displays a body on the screen. It also updates the position of the body based on its velocity vector. Even though I choose to display things as characters (so everyone can watch the simulation), it should be relatively easy to alter the code to work with graphic images.
Let's look at some examples. Example 1 is a simulation of the Earth and Moon. The source is straight forward. We create the bodies, define their physical parameters (mass, location, and velocity), and then let the UNIVERSE take care of the rest. The numbers used in all the methods are real. So if the physical laws on which the system is built match reality, the Moon should orbit the earth every 27.3 days. In the example, one day is equal to one tick.
/*------------------------------------------------
Simulates the orbital dynamics of the Earth
and Moon.
Todd King
--------------------------------------------------*/
#include "simul.hpp"
#include "simulscr.hpp"
main() {
BODY earth;
BODY moon;
UNIVERSE universe;
earth.set_mass(5.98e24);
earth.set_position(5.0e8, 5.0e8);
earth.set_icon('E');
moon.set_mass(7.36e22);
moon.set_position (5.0e8, 8.8e8);
moon.set_icon('M');
moon.set_velocity(-1020.0, 0.0);
universe.service(&earth);
universe.service(&moon);
universe.big_bang();
Another thing you will observe when you run this example is that the Earth/Moon pair moves gradually to the left. While trying to understand this result, I realized that the simulation system had revealed the notion of "center of mass" for systems, even though this notion was not bred into the system. The center of mass for a system is the point that two (or more) orbiting bodies revolve around. It is this point that moves at the velocity of the system as a whole. So when the Moon was given its initial velocity, the center of mass was given the same velocity. Hence the entire system moves to the left. Figure 1 illustrates this system.
Example 1 is a verification-type simulation that demonstrates that the simulation system produces substantiated results. This allows us to perform a hypothetical simulation and trust the results. Let's look at what would happen if Planet X passed through our system. Assume that the planet is twice as big as the Moon. Listing Three is the source for this simulation. The system is the same as that presented in Example 1, with the addition of a third body that moves through the system from the lower-left corner to the upper-right corner. When running this example you'll see just how disruptive such an occurrence would be, even if Planet X doesn't collide with any other body. Figure 2, shows the path of Planet X through the Earth/Moon system.
Simulation systems have broader applications than the model presented here. They are useful in developing and testing new theories. With simulation systems you can change the physical laws and constants and observe the results. One such constant is G, the gravitational constant, which can be found in Listing Four. Simulations also allow us to duplicate events in nature at a faster rate. For example, it takes the moon a lot less than 27.3 days of wall clock time to revolve about the earth in Example 1.
Games are a form of simulation systems, and with a little work you can create a game with this simulation system by adding a space ship body to one of the examples and allowing the user to control the rockets on it. The goal would be to obtain an orbit about one of the other bodies. Such a skill could be useful as we begin to move off the surface of our planet.
The system presented here has a few limitations. First, the maximum distance allowed between any two objects is the square root of the largest possible value for a double. On a PC, this value is 1.3e154. Next, the number of bodies that can exist in a UNIVERSE are limited. This is set with the parameter MAX_BODIES and has an upper boundary that is the same as for an int, and is limited because we use an array to track the bodies. Converting this to a link list would eliminate this limitation.
Rusnick, Robert, and Halliday, David. Physics, Part 1. New York: John Wiley & Sons, 1977.
Negoita, Constantin Virgil, and Ralescu, Dan. Simulation, Knowledge-Based Computing and Fuzzy Statistics. New York: Van Nostrand Reinhold Company, Inc., 1987.
All source code for articles in this issue is available on a single disk. To order, send $14.95 (Calif. residents add sales tax) to Dr. Dobb's Journal, 501 Galveston Dr., Redwood City, CA 94063, or call 800-356-2002 (from inside Calif.) or 800-533-4372 (from outside Calif.). Please specify the issue number and format (MS-DOS, Macintosh, Kaypro).
_FORCED-BASED SIMULATIONS_ by Todd King
[LISTING ONE]
Copyright © 1989, Dr. Dobb's Journal
/*-------------------------------------
Basic object classes and method
definitions for simulation software
File: simul.hpp
Todd King
----------------------------------------*/
#include <stream.hpp>
#include <disp.h>
#include <conio.h>
#include <math.h>
#include <time.h>
#include <math.h>
#include "simconst.h"
#define ESC 27
#define MAX_BODY_POOL 100
typedef struct VECTOR_2D {
double x, y;
};
class BODY {
VECTOR_2D world;
VECTOR_2D velocity;
double mass;
double gmass;
char icon;
public:
BODY();
set_mass(double m);
set_velocity(double x, double y);
set_position(double x, double y);
apply_force(VECTOR_2D from, double amount);
VECTOR_2D position();
double get_gmass();
update();
set_icon(char c);
};
/* Distance and time units are converted to screen units and ticks */
BODY::set_mass(double m) {
double pow();
mass = m;
gmass = G * m;
};
BODY::set_velocity(double x, double y) {
velocity.x = x;
velocity.y = y;
};
BODY::set_position(double x, double y) {
world.x = x;
world.y = y;
};
BODY::apply_force(VECTOR_2D from, double gmass) {
VECTOR_2D d;
double rs;
double v;
double r;
d.x = world.x - from.x;
d.y = world.y - from.y;
rs = (d.x * d.x) + (d.y * d.y);
if(rs != 0.0) { // there's a seperation
r = sqrt(rs);
v = (gmass / rs) * SECS_PER_TIC;
velocity.x += v * d.x / r;
velocity.y += v * d.y / r;
}
};
BODY::BODY() {
world.x = 0;
world.y = 0;
velocity.x = 0;
velocity.y = 0;
icon = '*';
};
VECTOR_2D BODY::position() {
VECTOR_2D vec;
vec.x = world.x;
vec.y = world.y;
return(vec);
};
double BODY::get_gmass() {
return(gmass);
}
BODY::set_icon(char c) { icon = c; }
class UNIVERSE {
unsigned int body_cnt;
BODY *body_pool[MAX_BODY_POOL];
public:
UNIVERSE();
service(BODY *bptr);
big_bang();
};
UNIVERSE::UNIVERSE() {
body_cnt = 0;
};
UNIVERSE::service(BODY *bptr) {
if(body_cnt >= MAX_BODY_POOL) return(0);
body_pool[body_cnt] = bptr;
body_cnt++;
};
UNIVERSE::big_bang() {
int i, j;
init_screen();
print_message(" Press ESC to stop.");
for(;;) {
print_tick();
if(kbhit()) {
switch(getch())
{
case ESC:
return(0);
}
}
sleep(0.1);
/* Let each body influence all others */
for(i = 0; i < body_cnt; i++) {
for(j = 0; j < body_cnt; j++) {
if(j == i) continue; // don't apply force to self
body_pool[i]->apply_force(body_pool[j]->position(),
body_pool[j]->get_gmass());
}
}
/* Display all bodies */
for(i = 0; i < body_cnt; i++) {
body_pool[i]->update();
}
}
deinit_screen();
};
sleep(double seconds) {
time_t ltime1, ltime2;
time(<ime1);
time(<ime2);
while(difftime(ltime1, ltime2) < seconds) {
time(<ime2);
}
}
[LISTING TWO]
/*--------------------------------------------
Methods which are related to screen I/O.
File: simulscr.hpp
Todd King
----------------------------------------------*/
#include "simconst.h"
#define DISPLAY_Y 24
#define DISPLAY_X 80
BODY::update() {
extern VECTOR_2D Extent;
VECTOR_2D screen;
screen.x = DISPLAY_X * (world.x / EXTENT_X);
screen.y = DISPLAY_Y * (world.y / EXTENT_Y);
if( (screen.x < DISPLAY_X && screen.x >= 0.0) &&
(screen.y < DISPLAY_Y && screen.y >= 0.0) ) {
disp_move(DISPLAY_Y - (int) screen.y, (int) screen.x);
disp_printf(" ");
}
world.x += velocity.x * SECS_PER_TIC;
world.y += velocity.y * SECS_PER_TIC;
screen.x = DISPLAY_X * (world.x / EXTENT_X);
screen.y = DISPLAY_Y * (world.y / EXTENT_Y);
if( (screen.x < DISPLAY_X && screen.x >= 0.0) &&
(screen.y < DISPLAY_Y && screen.y >= 0.0) ) {
disp_move(DISPLAY_Y - (int) screen.y, (int) screen.x);
disp_printf("%c", icon);
}
disp_move(0,0);
};
init_screen() {
disp_open();
disp_move(0, 0);
disp_eeop();
}
deinit_screen() {
disp_close();
}
print_tick() {
static unsigned int Tick = 0;
disp_move(0, 0);
disp_printf("Tick: %u", Tick);
Tick++;
}
print_message(char mesg[]) {
disp_move(24,0);
disp_printf(mesg);
}
[LISTING THREE]
/*------------------------------------------------
Simulation of what would happen if a planet
about twice the size of the Moon passed
close to the earth within the Moon's orbit.
Todd King
--------------------------------------------------*/
#include "simul.hpp"
#include "simulscr.hpp"
main() {
BODY earth;
BODY moon;
BODY planet_x;
UNIVERSE universe;
earth.set_mass(5.98e24);
earth.set_position(5.0e8, 5.0e8);
earth.set_icon('E');
moon.set_mass(7.36e22);
moon.set_position(5.0e8, 8.8e8);
moon.set_icon('M');
moon.set_velocity(-1020.0, 0.0);
planet_x.set_mass(14.8e22);
planet_x.set_position(1.0, 1.0e4);
planet_x.set_icon('X');
planet_x.set_velocity(1800, 2000);
universe.service(&earth);
universe.service(&moon);
universe.service(&planet_x);
universe.big_bang();
}
[LISTING FOUR]
/*-----------------------------------------------------
Various constants which affect the simulation
system. Units are in meters, seconds and Kilograms.
File: simconst.h
Todd King
-------------------------------------------------------*/
#define G -6.67e-11 /* Gravitational constant */
#define SECS_PER_TIC 86400 /* One day */
#define EXTENT_X 10e8 /* Width of displayed universe */
#define EXTENT_Y 10e8 /* Hieght of displayed universe */
[EXAMPLE 1]
/*-----------------------------------------------
Simulates the orbital dynamics of the Earth
and Moon.
Todd King
-------------------------------------------------*/
#include "simul.hpp"
#include "simulscr.hpp"
main() {
BODY earth;
BODY moon;
UNIVERSE universe;
earth.set_mass(5.98e24);
earth.set_position(5.0e8, 5.0e8);
earth.set_icon('E');
moon.set_mass(7.36e22);
moon.set_position(5.0e8, 8.8e8);
moon.set_icon('M');
moon.set_velocity(-1020.0, 0.0);
universe.service(&earth);
universe.service(&moon);
universe.big_bang();
}