Computer Science Canada

Solar System Simulation

Author:  Cervantes [ Fri Oct 15, 2004 3:39 pm ]
Post subject:  Solar System Simulation

I've recently been thinking about making a solar system simulation, using Newtonian physics. Could anyone perhaps give me a quick rundown of what I'll need to know (in terms of the physics, not the coding Wink)

Author:  josh [ Fri Oct 15, 2004 3:51 pm ]
Post subject: 

by simulation do you mean a mathematical one or a 2D or a 3D one???

that is a sick idea by the way Very Happy Idea

Author:  Cervantes [ Fri Oct 15, 2004 3:58 pm ]
Post subject: 

2D, using math.
Just the basics, at first.

Author:  josh [ Fri Oct 15, 2004 4:12 pm ]
Post subject: 

I might have some books that mention it, I will try and find some info

Author:  Martin [ Sat Oct 16, 2004 2:15 am ]
Post subject: 

Newton's equations are wrong.

Author:  Cervantes [ Sat Oct 16, 2004 6:30 am ]
Post subject: 

I know, but I still want to use them.
Using Newton's equations, don't we end up with circular orbits? And with the General Theory of Relativity, we end up with elliptical oribts?

Author:  Andy [ Sat Oct 16, 2004 1:41 pm ]
Post subject: 

just use centripetal acceleration formulas... where MAc=Fg... with Ac=V^2/r or Ac=4pi^2rf^2

Author:  Martin [ Sun Oct 17, 2004 1:29 am ]
Post subject: 

Yeah, ignore friction (it's physics, real world forces don't apply), and just give them an initial velocity perpindicular to the object they're to be orbitting, and it should work. Remember, planets pull on each other Wink

And for god's sake, make it OO.

Author:  Cervantes [ Mon Oct 25, 2004 6:37 pm ]
Post subject: 

dodge_tomahawk wrote:
just use centripetal acceleration formulas... where MAc=Fg... with Ac=V^2/r or Ac=4pi^2rf^2

would you mind explaining that a little more please? Confused

Martin wrote:
Remember, planets pull on each other "wink:

Oh god, I know! I plan on making an initialization section of the program where you can edit the mass of the planets (along with various other things), in which case some funky things could happen.. perhaps with the right mass, the solar system will orbit earth! Laughing

Tony showed me a way of creating a gravity field a while back, and it worked well enough for its purpose at that time, but what would it do in a solar system simulation?

code:

View.Set ("graphics:600;600,nobuttonbar")
var ball :
    record
        x, y, dx, dy, GFdist, GFxdist, GFydist, GFangle, GFforce : real
    end record

ball.x := 300
ball.y := 500
ball.dx := 0.5
ball.dy := 0

var gravityfield :
    record
        x, y, strength : int
    end record

gravityfield.x := 300
gravityfield.y := 300
gravityfield.strength := 150

const highestspeed := 3
const ballradius := 6

function MathDistance (x1, y1, x2, y2 : real) : real
    result (((x2 - x1) ** 2) + ((y2 - y1) ** 2)) ** 0.5
end MathDistance

setscreen ("offscreenonly")
loop
    cls

    ball.x += ball.dx
    ball.y += ball.dy

    if ball.x > maxx - ballradius or ball.x < ballradius then
        ball.dx := -ball.dx
        ball.x += (2 * ball.dx)       %move it 2 steps ahead to help prevent sticking
    end if
    if ball.y > maxy - ballradius or ball.y < ballradius then
        ball.dy := -ball.dy
        ball.y += (2 * ball.dy)      %move it 2 steps ahead to help prevent sticking
    end if


    %%%%%%%%%GRAVITY FIELD EFFECTS%%%%%%
    ball.GFdist := MathDistance (gravityfield.x, gravityfield.y, ball.x, ball.y)
    ball.GFforce := gravityfield.strength * (1 / (ball.GFdist * ball.GFdist))
    if ball.x not= gravityfield.x then
        ball.GFxdist := gravityfield.x - ball.x
        ball.GFydist := gravityfield.y - ball.y
        ball.GFangle := arctand (ball.GFxdist / ball.GFydist)
    else
        if ball.y > gravityfield.y then
            ball.GFangle := 90
        else
            ball.GFangle := 270
        end if
    end if
    if ball.x < gravityfield.x then
        ball.GFangle := ball.GFangle + 180
    end if
    ball.dx -= ball.GFforce * cosd (ball.GFangle)
    ball.dy -= ball.GFforce * sind (ball.GFangle)
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    if ball.dx > highestspeed then
        ball.dx := highestspeed
    elsif ball.dx < -highestspeed then
        ball.dx := -highestspeed
    end if
    if ball.dy > highestspeed then
        ball.dy := highestspeed
    elsif ball.dy < -highestspeed then
        ball.dy := -highestspeed
    end if

    drawfilloval (round (ball.x), round (ball.y), ballradius, ballradius, red)


    drawfilloval (gravityfield.x, gravityfield.y, 11, 11, black)
    View.Update
    delay (10)
end loop


See how the ball goes right, then left, then back? Gravity shouldn't do that! Surprised or should it? One of those good ol' fashioned "Educational Films" showed that our solar system is not revolving around the centre of the galaxy in a 2D plane, it is revolving around it and going up and down as well. Picture the Milky Way as a disc, totally flat (close enough..) then the sun and our solar system would be, at one point, below that disc, and then, as time passes, it moves up, passing through the disc, continues upwards, reaches its peak, and then comes back down. Kinda like a sine wave. Razz
Why, oh, why, is that? Eh

Author:  Tony [ Mon Oct 25, 2004 6:48 pm ]
Post subject: 

well you need to give the ball a little push at first to make it spin around the center...

code:

ball.dx := 1.75
ball.dy := 0

gravityfield.strength := 200

that works out much nicer Wink though the code is still not perfect, as you can tell yourself from the resulting animation... Though it was enough back in the day Laughing

Author:  Cervantes [ Mon Oct 25, 2004 6:51 pm ]
Post subject: 

Oh, I know, in order to get the orbit it would need that push. But, what I was trying to show, is that, without the push, you'd think it would go straight towards the gravity field, but in fact it wavers back and forth on its way there.

Author:  Catalyst [ Mon Oct 25, 2004 7:42 pm ]
Post subject: 

this might be of some interest

Author:  Paul [ Mon Oct 25, 2004 10:20 pm ]
Post subject: 

Orbits as I understand them are the result of a warping of spacetime due to a large mass present. Since there isn't friction to slow an object down, and given an initial velocity, the object just follows the simplest path, in an orbit. And wasn't there the gravitational force formula, y'know with the mass, distance and the universal gravitational constant. Posted Image, might have been reduced in size. Click Image to view fullscreen.
It doesn't sound easy though, but you like math Wink

and if you look at catalysts program wouldn't the direction and magnitude of the force being exerted on the moving circles have to be somehow balanced with the circle's initial velocity to have an orbit?

Author:  Cervantes [ Tue Oct 26, 2004 5:53 pm ]
Post subject: 

What is the universal gravitational constant?
And, catalyst, is what you posted this only open source instead of .exe? My winRAR hates me, and gives me annoying error messages whenever I try to extract something. Could you, or anyone willing, perhaps send it to me as a winzip?
I don't want to get offtopic here, but does anyone know where to get a free full version of winRAR (not a trail version)?

Author:  Andy [ Tue Oct 26, 2004 6:17 pm ]
Post subject: 

6.67300 × 10-11 m3 kg-1 s-2 kepler discovered that

Author:  Catalyst [ Tue Oct 26, 2004 7:32 pm ]
Post subject: 

yes it is open source that
heres a zip Smile

Author:  zylum [ Tue Oct 26, 2004 11:17 pm ]
Post subject: 

hey thats pretty cool... but over 400 lines? here's a simpler gravitation class i made:

code:
class SYSTEM
    export ADDBODY, DRAWSYSTEM, UPDATEBODY, UPDATEBODIES, SIMULATE, CHANGEG

    type BODYDATA :
        record
            _X : real
            _Y : real
            _VX : real
            _VY : real
            _AX : real
            _AY : real
            _MASS : real
            _COLOR : int
        end record

    var BODY : flexible array 1 .. 0 of BODYDATA
    var _G := 3.0


    fcn GETANGLE (x1, y1, x2, y2 : real) : real
        var dx, dy, ratio, angle : real
        dx := x1 - x2
        dy := y1 - y2
        if dx not= 0 then
            ratio := dy / dx
        else
            ratio := dy / 0.00000001
        end if
        angle := arctand (abs (ratio))
        if dx < 0 then
            angle := 180 - angle
        end if
        if dy < 0 then
            angle := 360 - angle
        end if
        result angle
    end GETANGLE

    proc CHANGEG (g : real)
        _G := g
    end CHANGEG

    proc ADDBODY (x, y, vx, vy, mass : real, clr : int)
        var i := upper (BODY) + 1
        new BODY, i
        BODY (i)._X := x
        BODY (i)._Y := y
        BODY (i)._VX := vx
        BODY (i)._VY := vy
        BODY (i)._MASS := mass
        BODY (i)._COLOR := clr
        BODY (i)._AX := 0
        BODY (i)._AY := 0
    end ADDBODY

    proc UPDATEBODY (i : int)
        BODY (i)._AX := 0
        BODY (i)._AY := 0
        for j : 1 .. upper (BODY)
            if i ~= j then
                BODY (i)._AX += ((_G * BODY (i)._MASS * BODY (j)._MASS) / ((BODY (i)._X - BODY (j)._X) ** 2 + (BODY (i)._Y - BODY (j)._Y) ** 2))
                    * (cosd (GETANGLE (BODY (j)._X, BODY (j)._Y, BODY (i)._X, BODY (i)._Y))) / BODY (i)._MASS
                BODY (i)._AY += ((_G * BODY (i)._MASS * BODY (j)._MASS) / ((BODY (i)._X - BODY (j)._X) ** 2 + (BODY (i)._Y - BODY (j)._Y) ** 2))
                    * (sind (GETANGLE (BODY (j)._X, BODY (j)._Y, BODY (i)._X, BODY (i)._Y))) / BODY (i)._MASS
            end if
        end for
        BODY (i)._VX += BODY (i)._AX
        BODY (i)._VY += BODY (i)._AY
        BODY (i)._X += BODY (i)._VX
        BODY (i)._Y += BODY (i)._VY
    end UPDATEBODY

    proc UPDATEBODIES
        for i : 1 .. upper (BODY)
            UPDATEBODY (i)
        end for
    end UPDATEBODIES

    proc DRAWSYSTEM
        for i : 1 .. upper (BODY)
            drawfilloval (round (BODY (i)._X), round (BODY (i)._Y), ceil (sqrt (BODY (i)._MASS)) + 1, ceil (sqrt (BODY (i)._MASS)) + 1, BODY (i)._COLOR)
        end for
    end DRAWSYSTEM

    proc SIMULATE
        var x, y, b : int
        loop
            mousewhere (x, y, b)
            UPDATEBODIES
            DRAWSYSTEM
            View.Update
            exit when b = 1
            cls
        end loop
    end SIMULATE
end SYSTEM



setscreen ("offscreenonly, graphics:max;max, nobuttonbar")

var SS : ^SYSTEM
new SYSTEM, SS

SS -> ADDBODY (maxx / 2, maxy / 2, 0, 0, 5000, 7)
SS -> ADDBODY (maxx / 2 + 400, maxy / 2, 0, 5, 20, 7)
SS -> ADDBODY (maxx / 2 + 420, maxy / 2, 0, 6.5, 1, 7)

SS -> SIMULATE

Author:  Tony [ Wed Oct 27, 2004 12:28 am ]
Post subject: 

your object representation is quite a bit off...

based on mass, radius does not equate to ceil (sqrt (BODY (i)._MASS)) + 1 Confused density of a 3D object should be considered.

V = (4/3)(pi)(R^3) (for the sphere) where V could be replaced by _MASS (assuming same density for all objects)... thus we get

R = (3*pi*_MASS/4)**(1/3) I'd recommend calculating that value just once and keeping it in object's record as well Wink

Author:  Cervantes [ Wed Oct 27, 2004 7:49 am ]
Post subject: 

Ahh, but the densities are not all the same. Earth has the highest density, Saturn has the lowest. I was thinking of hardcoding the values of each planets density and radius and mass.

Thanks very much for the zip catalyst, and also thanks to zylum. This will indeed be a complex program for me to make, but I'll try my best Smile

Author:  wtd [ Wed Oct 27, 2004 10:21 pm ]
Post subject: 

I like your use of data structures, but if I may suggest... there are lots of places where you have x, y pairs. Create a Point or Coordinate, or Vector record type and group them together.


: