diff --git a/lib/plots.py b/lib/plots.py index d3029d5..29f793f 100755 --- a/lib/plots.py +++ b/lib/plots.py @@ -134,9 +134,8 @@ def display_parameters(E,L,sma,ecc,parameters,savename=""): fig3 = plt.figure(figsize=(15,7)) ax3 = fig3.add_subplot(111) - for i in range(len(L)): - ax3.plot(np.arange(sma[i].shape[0])*step[i]/yr, sma[i], label="a (step of {0:.2e}s)".format(step[i])) - ax3.plot(np.arange(ecc[i].shape[0])*step[i]/yr, ecc[i], label="e (step of {0:.2e}s)".format(step[i])) + ax3.plot(np.arange(sma[-1].shape[0])*step[-1]/yr, sma[i], label="a (step of {0:.2e}s)".format(step[-1])) + ax3.plot(np.arange(ecc[-1].shape[0])*step[-1]/yr, ecc[i], label="e (step of {0:.2e}s)".format(step[-1])) ax3.set(xlabel=r"$t \, [yr]$", ylabel=r"$a \, [au] \, or \, e$") ax3.legend() fig3.suptitle("Semi major axis and eccentricity "+title2) diff --git a/main_concurrent.py b/main_concurrent.py index 067673d..ea9ac29 100755 --- a/main_concurrent.py +++ b/main_concurrent.py @@ -14,9 +14,9 @@ from lib.units import * def main(): #initialisation m = np.array([1., 1., 1e-1],dtype=np.longdouble)*Ms/Ms # Masses in Solar mass - a = np.array([1., 1., 5.],dtype=np.longdouble)*au/au # Semi-major axis in astronomical units - e = np.array([0., 0., 0.],dtype=np.longdouble) # Eccentricity - psi = np.array([0., 0., 0.],dtype=np.longdouble)*np.pi/180. # Inclination of the orbital plane in degrees + a = np.array([1., 1., 5.],dtype=np.longdouble)/2.*au/au # Semi-major axis in astronomical units + e = np.array([0., 0., 1./4.],dtype=np.longdouble) # Eccentricity + psi = np.array([0., 0., 45.],dtype=np.longdouble)*np.pi/180. # Inclination of the orbital plane in degrees x1 = np.array([0., -1., 0.],dtype=np.longdouble)*a[0]*(1.+e[0]) x2 = np.array([0., 1., 0.],dtype=np.longdouble)*a[1]*(1.+e[1]) @@ -29,19 +29,20 @@ def main(): v = np.array([v1, v2, v3],dtype=np.longdouble) #integration parameters - duration, step = 100*yr, np.array([1./(365.25*24.), 12./(365.25*24.), 24./(365.25*24.)],dtype=np.longdouble)*yr #integration time and step in seconds + duration, step = 5000.*yr, np.array([10./(365.25)],dtype=np.longdouble)*yr #integration time and step in seconds step = np.sort(step)[::-1] integrator = "leapfrog" n_bodies = 3 display = False - savename = "{0:d}bodies_conc_{1:s}".format(n_bodies, integrator) + savename = "{0:d}bodies_psi45_{1:s}".format(n_bodies, integrator) bodies, bodysyst = [],[] for j in range(n_bodies): bodies.append(Body(m[j], q[j], v[j])) bin_syst = System(bodies[0:2]) dyn_syst = System(bodies, main=True) - bodysyst = [[deepcopy(bin_syst), deepcopy(dyn_syst)] for _ in range(n_bodies)] + bodysyst = [[deepcopy(bin_syst), deepcopy(dyn_syst)] for _ in range(len(step))] + #simulation start exe = ProcessPoolExecutor() future_ELae = [] diff --git a/plots/3bodies_psi45_leapfrog_a_e.png b/plots/3bodies_psi45_leapfrog_a_e.png new file mode 100644 index 0000000..ccb3998 Binary files /dev/null and b/plots/3bodies_psi45_leapfrog_a_e.png differ diff --git a/plots/3bodies_psi45_leapfrog_dEm.png b/plots/3bodies_psi45_leapfrog_dEm.png new file mode 100644 index 0000000..48bab20 Binary files /dev/null and b/plots/3bodies_psi45_leapfrog_dEm.png differ diff --git a/plots/3bodies_psi45_leapfrog_dL2.png b/plots/3bodies_psi45_leapfrog_dL2.png new file mode 100644 index 0000000..c561d07 Binary files /dev/null and b/plots/3bodies_psi45_leapfrog_dL2.png differ