Author Topic: [Solved] Two body simulation (Read 3034 times)

ElementCoder · « **on:** May 12, 2014, 09:25:35 am »

I'm trying to make a two body simulation, but it doesn't seem to work

Plotting x vs y should be circular or elliptical right? Am I choosing the wrong starting conditions or is my way of doing it just wrong?

Code: [Select]

#!/usr/bin/env python
'''
Two body simulation.

'''
from __future__ import division
import math
import matplotlib.pyplot as pyp
import numpy as np

# Gravitational constant (N*m^2*kg^-2)
G = 6.674e-11
# A giga year in seconds.
gyrts = math.pi * 1e7 * 1e9

def Fgrav(M, m, r):
    ''' Gravitational force of two bodies attracting each other.
    '''
    return -G * M * m / r**2
    
# Solar mass (kg)
M = 1.989e30
# Earth mass (kg)
m = 5.9721e24

# Initial conditions.
# Position (m).
x, y = 150e9, 0
X, Y = 0, 0

# Velocity (m/s).
vx, vy = 0, 30e3
Vx, Vy = 0, 0

time = 5*math.pi*1e7
t = 0
dt = 1000
xarr, yarr, rarr, tarr = [], [], [], []
while t < time:
    # Update the positions.
    x += vx * dt
    y += vy * dt
    # Calculate the separation vector.
    dx = x - X
    dy = y - Y
    r2 = dx*dx + dy*dy + 1e-4
    r = np.sqrt(r2)
    #print 'dx, dy: ', dx, dy
    # Calculate components of the gravitational force.
    fgrav_x = Fgrav(M, m, r) * dx / r
    fgrav_y = Fgrav(M, m, r) * dy / r
    #print 'Fgrav: ', fgrav_x, fgrav_y
    # Update the velocities.
    vx += (fgrav_x / m) * dt
    vy += (fgrav_y / m) * dt
    Vx += (fgrav_x / M) * dt
    Vy += (fgrav_y / M) * dt
    #print 'vx, vy: ', vx, vy
    
    X += Vx * dt
    Y += Vy * dt
    
    # Store all data.
    xarr.append(dx); yarr.append(dy); rarr.append(r)
    tarr.append(t)
    
    t += dt
    


fig, (ax1, ax2) = pyp.subplots(2, 'nom')
ax1.plot(xarr, yarr)
ax1.set_xlabel('$x$ (m)'); ax1.set_ylabel('$y$ (m)')
ax2.plot(tarr, rarr)
ax2.set_xlabel('$t$ (s)'); ax2.set_ylabel('$r$ (m)')
fig.tight_layout()

pyp.show()

Edit: Updated with non-broken code.

lkj · « **Reply #1 on:** May 12, 2014, 01:32:21 pm »

The starting conditions are wrong (if you're simulating the earth and the sun):
- The initial positions: the real distance is ~150e9, but placing earth at (100e9, 50e9) gives us only a distance of 110e9. I'd choose (0, 0) for the sun and (150e9, 0) for earth.
- The initial speeds: why (3e3, 5e3)? If you choose (150e9, 0) as the starting position of earth, it's about (0, 30e3).

I don't have mathplotlib and numpy installed, so I can't test if it really works with my values. The results seem more reasonable, but I may be wrong

ElementCoder · « **Reply #2 on:** May 12, 2014, 02:17:01 pm »

Those were just some test conditions.

I forgot a minus sign in the gravitational force. Everything works fine now

lkj · « **Reply #3 on:** May 12, 2014, 02:53:30 pm »

Ah lol, that would actually have been quite obvious from the results. Good to hear

Author Topic: [Solved] Two body simulation (Read 3034 times)

ElementCoder

[Solved] Two body simulation

lkj

Re: Two body simulation

ElementCoder

Re: Two body simulation

lkj

Re: [Solved] Two body simulation