Using physics to simulate moon in orbit around planet?

Quick question: is there a reason you are all avoiding vectors? I actually think that the equations are clearer and cleaner when expressed using vectors.

No reason, other than dave1707 had done it that way and I’m not too offay with vectors…

Below rewritten with vectors


displayMode(FULLSCREEN)
--supportedOrientations(PORTRAIT)

function setup()
    leaveTrails = false
    if leaveTrails then
        backingMode(RETAINED)
    end
    bodies = {}
    
    planetarySetup = 1
    
    if planetarySetup == 1 then
        --solar system
        --gravitation constant for scaling
        G = .005
        --initial planets
        --sun
        bodies[1] = { m= 40000, location = vec2(WIDTH/2, HEIGHT/2), velocity = vec2(-0.005, 0), r=200, g=255, b=0 }
        --planets
        bodies[2] = { m= 15, location = vec2(WIDTH/2, HEIGHT/2 - 65), velocity = vec2(1.8, 0), r=200, g=0, b=255 }
        bodies[3] = { m= 20, location = vec2(WIDTH/2, HEIGHT/2 - 100), velocity = vec2(1.5, 0) , r=0, g=255, b=0 }
        bodies[4] = { m= 30, location = vec2(WIDTH/2, HEIGHT/2 - 150), velocity = vec2(1.3, 0) , r=255, g=0, b=0}
        bodies[5] = { m= 40, location = vec2(WIDTH/2, HEIGHT/2 - 180), velocity = vec2(1.2, 0), r=0, g=200, b=255}
        bodies[6] = { m= 60, location = vec2(WIDTH/2, HEIGHT/2 - 300), velocity = vec2(0.85, 0), r= 255, g=255, b=255}
        --moons
        bodies[7] = { m= 0.1, location = vec2(WIDTH/2, HEIGHT/2 - 155), velocity = vec2(1.465, 0), r=150, g=150, b=150 }
    end
    
    if planetarySetup == 2 then
        --2 planets and a moon
        G = 1
        --initial planets
        bodies[1] = { m = 400, location = vec2(300, 500), velocity = vec2(1.4, .4) , r=0, g=200, b=255}
        bodies[2] = { m = 400, location = vec2(300, 600), velocity = vec2(-1.4, -.4), r=255, g=200, b=0 }
        bodies[3] = { m = 4, location = vec2(300, 480), velocity = vec2(-3.1, .4), r=130, g=130, b=130 }
    end
    
end

function draw()
    if leaveTrails == false then
        background(40,40,50)
    end
    strokeWidth(0)
    fill(255)
    
    applyGravity()
    
    moveBodies()
    
    -- plot bodies
    for k,v in pairs(bodies) do
        fill(v.r,v.g,v.b,255)
        ellipse(v.location.x, v.location.y, math.pow(v.m, .2)*2 + 2)
    end
end

function applyGravity()
    for k,v in pairs(bodies) do
        for k2,v2 in pairs(bodies) do
            --don't apply gravity to itself
            if k2 ~= k then
                --determine the distance between the planets (pythagorus theorem)
                distance = v.location:dist(v2.location)
                
                --determine force of gravity (newton F = G*(m1*m2)/r^2
                gravityForce = G*(v.m*v2.m)/(distance^2)
                
                --determine how much acceleration this will have on planet 1 (F = m*a -> a = F/m)
                acceleration = gravityForce/v.m
                
                --now apply the acceleration to the velocity multiply by the normalised vector between the locations
                gravityDirection = v.location - v2.location
                v.velocity = v.velocity - (acceleration * gravityDirection:normalize())
            end
        end
    end
end

function moveBodies() 
    for k,v in pairs(bodies) do
        v.location = v.location + v.velocity
    end
end

Wow. I didn’t think this was going to turn into a class for physics, I like all the formulas. I didn’t use vectors because for anyone not used to vectors, they turn into a black box. You put numbers in, you get numbers out (vector math), but you don’t know what happens to them. I just wanted to show a simple program that would display orbits under the influence of gravity that can’t be done with the way gravity works now in Codea. I like to keep things basic and let whoever uses my examples to build on them, @spacemonkey was a good example.

Hello, I’ve made this code when I tried to understand the physic engine. It uses applyForce to pull satellites toward the planet.


function setup()
    displayMode(FULLSCREEN)

    noSmooth()

    balls = {}
    touches = {}
    
    planet = {x=500,y=400,r=100,tspot=physics.body(CIRCLE, 50)}
    planet.tspot.position = vec2(planet.x, planet.y)
    planet.tspot.type = STATIC
    planet.tspot.restitution = 1
    
    physics.gravity(0,0)

    nextball = 0
end

function touched(touch)
    if touch.state == ENDED then
        touches[touch.id] = nil
    else
        --if touches[touch.id] == nil then
        touches[touch.id] = touch
        --end
    end    
end    

function touchActions()
    for k,v in pairs(touches) do 
        if CurrentTouch.state == ENDED then
            --if there are no current touches then we kill all current touches to avoid bugged ball producers
            touches[k] = nil
        else

            --add a new ball at the touch location
            size = math.random(10,20)
            tspot = physics.body(CIRCLE, size)
            tspot.position = vec2(v.x+math.random(-1,1), v.y+math.random(-1,1))
            tspot.restitution = 1
            tspot.linearVelocity = vec2(300,0)
            --tspot:applyForce(vec2(20,0))

            balls[nextball] = { tspot = tspot, size = size * 2, r = math.random(30,255), g = math.random(30,255), b = math.random(30,255) }

            nextball = nextball + 1
        end    
    end
end

function draw()
    background(0, 0, 0, 255)
    strokeWidth(0)
    
    fill(31, 91, 33, 255)
    ellipse(planet.x,planet.y,planet.r)

    for k,v in pairs(balls) do
        if v.tspot.x < -100 or v.tspot.x > WIDTH + 100 or v.tspot.y < -100 or v.tspot.y > HEIGHT + 100 then
            balls[k].tspot:destroy()
            balls[k] = nil  
        else
            v.tspot:applyForce(0.3*v.size*vec2(planet.x-v.tspot.x, planet.y-v.tspot.y)/math.sqrt((planet.x-v.tspot.x)^2 + (planet.y-v.tspot.y)^2))
            fill(v.r, v.g, v.b, 255)
            ellipse(v.tspot.x, v.tspot.y, v.size)
        end
    end

    touchActions()
end

Here is my original program using physics.body and applyForce. There isn’t any way to use gravityScale except to make it 0, and use my original gravity.


displayMode(FULLSCREEN)
supportedOrientations(PORTRAIT_ANY)

function setup()
    s1=physics.body(CIRCLE,5)
    s1.x=300
    s1.y=500
    s1.gravityScale=0
    s1:applyForce(vec2(0,6))
    s1g=400
    
    s2=physics.body(CIRCLE,5)
    s2.x=400
    s2.y=500
    s2.gravityScale=0
    s2:applyForce(vec2(0,-6))
    s2g=400
end

function draw()
    background(40,40,50)
    fill(255)
    
    d = s2g / ((s1.x - s2.x)^2 + (s1.y - s2.y)^2 )^1.5
    s1dx = ( s1.x - s2.x ) * d * -5
    s1dy = ( s1.y - s2.y ) * d * -5

    d = s1g / ((s2.x - s1.x)^2 + (s2.y - s1.y)^2 )^1.5
    s2dx = ( s2.x - s1.x ) * d * -5
    s2dy = ( s2.y - s1.y ) * d * -5
    
    s1:applyForce(vec2(s1dx,s1dy))
    s2:applyForce(vec2(s2dx,s2dy))
    
    ellipse(s1.x,s1.y,10)
    ellipse(s2.x,s2.y,10)
end