I have converted an equation to the first order ODE and now i would like to solve the motion equations for multiple periods with given conditions.
The equations shall be solved with the following 0 values: x(0), y(0), vx(0), vy(0) = 3, 1, 2, 1.3 × π
This is the first order ODE, derived from motion equation of Newton's gravitational law for cellestial bodies:
How can I solve this in python? Could I use Runge Kutta or Keplers methods? I feel like I'm doing something wrong
import math
import matplotlib.pyplot as plt
import numpy as np
g = 6.674 * 10**(-11)
M = 4*(math.pi**2)/g
x0 = 3 #x-position of the center or h
y0 = 1 #y-position of the center or k
vx0 = 2
vy0 = 1.3* math.pi
#Trying second order
def RK2(y0, f, tlist):
t = [tlist[0]]
tf = tlist[1]
dt = tlist[2]
y = [y0]
while t[-1] < tf:
k1 = dt * f(y[-1],t[-1])
k2 = dt * f(y[-1]+0.5*k1 , t[-1]+0.5*dt) ##Take half step
y.append( y[-1] + k2 )
t.append(t[-1] + dt)
return(np.array(y), np.array(t) )
