Skip to content

Instantly share code, notes, and snippets.

@BlogBlocks
Created December 2, 2017 19:58
Show Gist options
  • Save BlogBlocks/24e7e4599081cf04c19a0484f1f733d7 to your computer and use it in GitHub Desktop.
Save BlogBlocks/24e7e4599081cf04c19a0484f1f733d7 to your computer and use it in GitHub Desktop.
orbits
#!/usr/bin/env python3
import math from turtle import *
# The gravitational constant G
G = 6.67428e-11
# Assumed scale: 100 pixels = 1AU. AU = (149.6e6 * 1000)
# 149.6 million km, in meters. SCALE = 250 / AU class Body(Turtle): """Subclass of Turtle representing a gravitationally-acting body. Extra attributes: mass : mass in kg vx, vy: x, y velocities in m/s px, py: x, y positions in m """ name = 'Body' mass = None vx = vy = 0.0 px = py = 0.0 def attraction(self, other): """(Body): (fx, fy) Returns the force exerted upon this body by the other body. """ # Report an error if the other object is the same as this one. if self is other: raise ValueError("Attraction of object %r to itself requested" % self.name) # Compute the distance of the other body. sx, sy = self.px, self.py ox, oy = other.px, other.py dx = (ox-sx) dy = (oy-sy) d = math.sqrt(dx**2 + dy**2) # Report an error if the distance is zero; otherwise we'll # get a ZeroDivisionError exception further down. if d == 0: raise ValueError("Collision between objects %r and %r" % (self.name, other.name)) # Compute the force of attraction f = G * self.mass * other.mass / (d**2) # Compute the direction of the force. theta = math.atan2(dy, dx) fx = math.cos(theta) * f fy = math.sin(theta) * f return fx, fy def update_info(step, bodies): """(int, [Body]) Displays information about the status of the simulation. """ print('Step #{}'.format(step)) for body in bodies: s = '{:<8} Pos.={:>6.2f} {:>6.2f} Vel.={:>10.3f} {:>10.3f}'.format( body.name, body.px/AU, body.py/AU, body.vx, body.vy) print(s) print() def loop(bodies): """([Body]) Never returns; loops through the simulation, updating the positions of all the provided bodies. """ timestep = 24*3600 # One day for body in bodies: body.penup() body.hideturtle() step = 1 while True: update_info(step, bodies) step += 1 force = {} for body in bodies: # Add up all of the forces exerted on 'body'. total_fx = total_fy = 0.0 for other in bodies: # Don't calculate the body's attraction to itself if body is other: continue fx, fy = body.attraction(other) total_fx += fx total_fy += fy # Record the total force exerted. force[body] = (total_fx, total_fy) # Update velocities based upon on the force. for body in bodies: fx, fy = force[body] body.vx += fx / body.mass * timestep body.vy += fy / body.mass * timestep # Update positions body.px += body.vx * timestep body.py += body.vy * timestep body.goto(body.px*SCALE, body.py*SCALE) body.dot(3) def main(): sun = Body() sun.name = 'Sun' sun.mass = 1.98892 * 10**30 sun.pencolor('yellow') earth = Body() earth.name = 'Earth' earth.mass = 5.9742 * 10**24 earth.px = -1*AU earth.vy = 29.783 * 1000 # 29.783 km/sec earth.pencolor('blue') # Venus parameters taken from # http://nssdc.gsfc.nasa.gov/planetary/factsheet/venusfact.html venus = Body() venus.name = 'Venus' venus.mass = 4.8685 * 10**24 venus.px = 0.723 * AU venus.vy = -35.02 * 1000 venus.pencolor('red') loop([sun, earth, venus])
if __name__ == '__main__': main()
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment