------------------------------------------------------- Orbit7: A Simple Gravity Simulator ------------------------------------------------------- by Mike Sussman <plexus23@aol.com> and Ben Collins-Sussman <sussman@red-bean.com> and others: <kfogel@red-bean.com>, <cmpilato@xnet.com> (C) 1999 under the GNU GPL version 2 ------------------------------------------------------- HISTORY This program was originally written by Mike Sussman, a young astrophysicist, as a very first C program. (Originally using DOS, Turbo C, and the Borland Graphics Library.) His brother Ben ported the code to Xaw, Lesstif, and then finally found an excuse to learn GTK+. REQUIREMENTS The ability to compile programs written for GTK+ 1.2 or higher. (The program was developed on a Linux 2.2.x system, but should build anywhere GTK is available.) BUILDING Just run `make`. GUI NOTES This program uses the "gtkdial" widget from the GTK Tutorial; it was originally meant as a simple example of writing a custom widget, but we find it provides a nice "Star Trek" feel to the user interface. PHYSICS NOTES Mike Sussman's description of the code is below. CONTACT & UPDATES Please feel free to hack on this program. It's fun. We welcome any changes or improvements! Email orbit7@red-bean.com. You can get the latest version from red-bean.com's anonymous CVS repository: export CVSROOT=:pserver:anonymous@cvs.red-bean.com:/usr/local/cvs cvs login cvs co orbit7 (The password for anonymous access is "the key".) ------------------------------------------------ here, incidentally, is a little explanation for that orbit engine...it's not totally thorough, but it should give a bit of an overview for someone trying to disect my trainwreck of code.... it all starts with a 2-dimensional array of doubles (far better than floats, as sensitivity to precision is very important in non-linear dynamics), where the horizontal axis is the hamiltonian expressed as 5 separate variables - (mass, x-position, y-position, x-velocity, y-velocity)... (in case you were wondering, a Hamiltonian is what is known as the state descriptor - if you know the hamiltonian, then you know you have necessary and sufficient information to predict the outcome of the system...in classical mechanics, the hamiltonian is the position and the momentum....however, in quantum mechanics, the hamiltonian must contain imaginary numbers - thus why you can only know the position or the momentum - because the other variable contains uninterpretable data - after all, what does it mean to say a neutron be traveling at 5i miles per second?) in my array, the 2nd and 3rd variables are the position, and since momentum is just mass times velocity, the 1st, 4th, and 5th variables are the momentum descriptor...so that's how i set up the columns...the first row is then the sun, while each subsequent row represents different planets.... something like: Mass x-positions y-position x-velocity y-velocity sun 5000 0 0 0 0 planet1 100 10 0 0 -10 planet2 100 5 5 5 -5 planet3 100 15 -5 -3 2 ...and so on, for as many planets as your little heart desires....i believe the code involves some kind of general random function (thanks again, karl) to generate some positions and velocities for the planets...actually what i do is give everything a standard radius from the central star, with a complementary tangential motion for that radial position, and then offset each planet's position and momentum by some random amount...currently that randomness is unweighted, but i suppose i could work in a function to give it a standard bellcurve randomness.... once this is all said and done, i calculate, for each planet and the sun, what its change in velocity must be based on the position of every other planet (since gravity is a function of the distance...) - this is where the real meat of the proceesor time occurs - it must calculate every interaction between every two bodies...(note that if you have 2 objects, there's only 1 interaction, but 3 objects means 3 interactions, 4 bodies is 6 interactions, 5 bodies is 10, etc...in other words, if the number of bodies is N, the number of interactions is (N^2 - N)/2...this eats up processor time very quickly as you increase the objects...i'll have to check the code again, but i'm afraid i might actually calculate each interaction twice...once from object A to object B, then again later from object B to object A...thanks to newton's law about every-action--means-an-opposite-and-equal-reaction, tho', we should only have to calculate each interaction once, and then just switch the direction of that vector for the opposite action... the acceleration due to gravity is exactly as you'll find in the textbooks - it's proportional to (mass of gravitating body)/(distance^2)...of course - this is just a proportion - the actual number is found by multiplying by a gauge constant...as the units of velocity, position, and mass are all arbitrary, this is the only part which requires some tweaking to figure out what will produce the desired affect - you don't want gravity too strong, or everything will fall into the sun and produce ugly anamolies when the distance between objects enters the same order of magnitude as the imprecision of a double variable...but, you don't want it too weak, or else the planets will go flying off in hyperbolic orbits, never to return - this makes for a very boring demo after the first few seconds. (Interestingly, gauge constants are also the only values in physics which cannot be derived from theory - they must be empirically determined...my sloppiness/arbitrariness is thus a consequence of our own slightly arbitrary universe...) anyway, each planets/sun now has a new velocity based on the strength of that gravitational interaction, and based on the new velocity, i calculate what its new position is...then we redraw all the planets/sun, and the process begins anew....thus, it goes something like this: acceleration = gauge constant * (SUM (from element 1 to N) of ((Mass of N)/distance^2)) (note that this is all vector addition) velocity += gravitational acceleration position += velocity redraw new positions. as you can see, it takes place as a series of steps - the size of the steps (just how much each planet moves between each subsequent calculation) is equivalent to the size of a step one would take to solve a numerical integration....it's technically arbitrary, but there are guidelines - you want them large enough to show visible swirling motion, but if you make them too large, you'll start introducing inaccuracies and the orbits will start to precess.... Big Fun. confused? good. write me email, and I'll explain more.
