simple solar system simulation - four objects moving with gravitation forces between them

gravorbit.py

import sys, time, math
from pyglet.gl import *
from euclid import *

GRAV_CONSTANT = 1

window = pyglet.window.Window(512,512)

def makeCircle(numPoints):
    verts = [0,0]
    for i in range(numPoints+1):
        angle = math.radians(float(i)/numPoints * 360.0)
        x = math.cos(angle)
        y = math.sin(angle)
        verts = verts + [x,y]
    vlist = pyglet.graphics.vertex_list(len(verts)/2, ('v2f', verts))
    return vlist

circle = makeCircle(16)

class Planet:
    def __init__(self, pos, vel, mass, color):
        self.pos = pos
        self.vel = vel
        self.mass = mass
        self.color = color
    def draw(self):
        glColor3f(self.color[0], self.color[1], self.color[2])
        glPushMatrix()
        glTranslatef(self.pos[0], self.pos[1], self.pos[2])
        circle.draw(GL_TRIANGLE_FAN)
        glPopMatrix()
    def update(self, dt, planets):
        force = Vector3(0, 0, 0)
        for p in planets:
            if p != self:
                d = p.pos - self.pos
                f = GRAV_CONSTANT * self.mass * p.mass / d.magnitude_squared()
                dir = d.normalized()
                force = force + dir * f
        acc = force / self.mass
        self.vel = self.vel + acc * dt
        self.pos = self.pos + self.vel * dt

planets = []
planets.append(Planet(Vector3(0,0,0), Vector3(0,0,0), 10000, [1,1,0]))
planets.append(Planet(Vector3(10,0,0), Vector3(0,28,0), 100.0, [0,1,0]))
planets.append(Planet(Vector3(-14,0,0), Vector3(0,-24,0), 50.0, [0,1,1]))
planets.append(Planet(Vector3(0,-16,0), Vector3(-30,0,0), 2.0, [1,1,1]))


@window.event
def on_draw():
    glClear(GL_COLOR_BUFFER_BIT)
    glMatrixMode(GL_PROJECTION)
    glLoadIdentity()
    glOrtho(-25,25, -25,25, -25,25)
    glMatrixMode(GL_MODELVIEW)
    for p in planets: p.draw()


def update(dt):
    for p in planets:
        p.update(dt, planets)


pyglet.clock.schedule_interval(update,1/60.0)
pyglet.app.run()