Merge branch 'revert-6e26de08' into 'main'

Revert "Added test binary star system, currently is unstable"

See merge request !3

core_functions.py

@@ -2,28 +2,6 @@ import numpy as np
 
 # True Acceleration Function
 def f_true(u_ki, index, u_k, masses):
-    """Returns acceleration for the i-th body in a two-body system.  You can find the overall acceleration for a body by summing the acceleration from all bodies in a n-body system. 
-    
-    Parameters
-    ----------
-    r_i : array
-        Position vector of body that we are finding acceleration for (i-th body)
-    r_j : array
-        Position vector of other body that is affecting the main body (j-th body)
-    m_j : float_like
-        Mass of other j-th body
-        
-    Returns
-    -------
-    a : array
-        Acceleration vector of output dynamics for i-th body
-    
-    """
-
-    #so apparently we need to calc r_dot (velocity), since the input is position
-    #vf = a*delta_t+vi
-    #take in v_i
-
     G = 6.6743e-11
     r1 = np.array([u_ki[0], u_ki[1], u_ki[2]])
     

test.py

@@ -22,7 +22,7 @@ for i in range(len(times)):
     rz.append(u[i][1][2])
 
 import matplotlib.pyplot as plt
-plt.axes(projection='3d')
-plt.plot(rx,ry,rz)
+
+plt.plot(rx,ry)
 
 plt.show()

test_TOI1338.py

-import numpy as np
-from core_functions import *
-
-#sun mass: 1.9885*10^30 https://nssdc.gsfc.nasa.gov/planetary/factsheet/sunfact.html
-#earth mass: 5.9724e+24
-#TOI-1338 data: https://arxiv.org/abs/2004.07783
-#orbital speed calc from https://archive.org/details/handbookofmathem00harr/page/386/mode/2up
-masses = [1.127*1.9885e+30, .3128*1.9885e+30, 33.0*5.9724e+30] #median mass in kg
-
-AU = 1.496e+11 #m
-
-G = 6.6743e-11
-
-i_bin = np.radians(89.696)
-i_pl = np.radians(89.37)
-long_ascend_pl = np.radians(0.91)
-
-e_bin = 0.15603
-a_bin = 0.1321*AU
-P_bin = 14.608559*3600*24 #s
-v_bin = 2*np.pi*a_bin/P_bin*(1 - 1/4*e_bin**2 - 3/64*e_bin**4 - 5/256*e_bin**6) #average orbital speed, m/s 
-
-e_pl = 0.0880
-a_pl = 0.4607*AU
-P_pl = 95.174*3600*24 #s
-v_pl = 2*np.pi*a_pl/P_pl*(1 - 1/4*e_pl**2 - 3/64*e_pl**4 - 5/256*e_pl**6) #average orbital speed, m/s
-
-u_0 = np.array([[0,a_bin*np.cos(i_bin),a_bin*np.sin(i_bin),0,v_bin*np.cos(i_bin),v_bin*np.sin(i_bin)],
-                [0,-a_bin*np.cos(i_bin),-a_bin*np.sin(i_bin),0,-v_bin*np.cos(i_bin),-v_bin*np.sin(i_bin)],
-                [a_pl*np.sin(long_ascend_pl), a_pl*np.cos(long_ascend_pl), 0, 0, 7823.6, 0]])
-
-T = 1000000
-
-delta_t = 10
-
-u, times = ivp_RK4(u_0, T, delta_t, masses)
-
-#binary star 0
-r0x = []
-r0y = []
-r0z = []
-#binary star 1
-r1x = []
-r1y = []
-r1z = []
-#planet
-rpx = []
-rpy = []
-rpz = []
-
-for i in range(len(times)):
-    r0x.append(u[i][0][0])
-    r0y.append(u[i][0][1])
-    r0z.append(u[i][0][2])
-    r1x.append(u[i][1][0])
-    r1y.append(u[i][1][1])
-    r1z.append(u[i][1][2])
-    rpx.append(u[i][2][0])
-    rpy.append(u[i][2][1])
-    rpz.append(u[i][2][2])
-
-import matplotlib.pyplot as plt
-plt.axes(projection='3d')
-# plt.plot(r0x, r0y, r0z)
-plt.plot(r1x, r1y, r1z)
-# plt.plot(rpx, rpy, rpz)
-
-plt.show()