Python. Physics


#1

I am new to programming. Can anyone tell me what is RK ? how can i run this program? and what formula is used in this program? please explain with details, Thank you

import RK
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline

“”" Implements d\vec X/dt = \vec V, d\vec V/dt = - M G \vec X/|X|^3.
Takes as input a time and
a 6-d vector. The first three components are the components of \vec
X, the next three are the components of \vec V. The parameter “t” is
not used.
“”"

def Kepler(t, y):

M = 1.0
G = 1.0

position = y[0:3]
velocity = y[3:6]

# position_dot = velocity
pos_dot = velocity;

r = np.sqrt(np.dot(position,position))
r3 = r*r*r

# formula for velocity_dot
vel_dot = - position  * M * G / r3;

# y_dot 
return np.concatenate((pos_dot, vel_dot))

“”" Source term for a binary. It has components m1 and m2. Currently
assumes circular orbits. Can adjust binary separation and mass
ratio.
“”"

def Binary(t, y):

binary_separation  = 1.0

mass_ratio = 1.0
M = 1.0
G = 1.0

m1 = M * mass_ratio / (1.+mass_ratio)
m2 = M / (1.+mass_ratio)

binary_frequency = np.sqrt(M/binary_separation**3)
binary_angle = t * binary_frequency
X1 = np.array([np.cos(binary_angle), np.sign(binary_angle),0]) / (1.0 + mass_ratio)
X2 = -X1 * mass_ratio


position = y[0:3]
velocity = y[3:6]

# position_dot = velocity
pos_dot = velocity;

p1 = position - X1
p2 = position - X2

r1 = np.sqrt(np.dot(p1,p1))
r2 = np.sqrt(np.dot(p2,p2))

# formula for velocity_dot
vel_dot = - p1  * m1 * G / r1**3 - p2 * m2 * G / r2**3

# y_dot 
return np.concatenate((pos_dot, vel_dot))

Initial position

x0 = 2
y0 = 0
z0 = 0

#Initial velocity
Vx = 0
Vy = np.sqrt(1.0/2.0)
Vz = 0

tend = 10000
dt = 0.01
y = np.array((x0,y0,z0, Vx, Vy, Vz), dtype=np.float64)
t = 0

while t < tend:

#(t, y, dt) = RK.RK45_Step(t, y, dt, Kepler, tol=1.0e-12)

(t, y, dt) = RK.RK45_Step(t, y, dt, Binary, tol=1.0e-12)

print (y[0], y[1])


#2

The formula in the kepler function looks like his third law…

P**2 is proportional to R ** 3

Here is a link to a page on Kepler’s Third Law…

http://astro.physics.uiowa.edu/ITU/glossary/keplers-third-law/


#3

Thank you so much for answering my question

Sincerely, Niloo