Initial commit

Exercise 1/exercise1.py

@@ -0,0 +1,68 @@
+import numpy as np
+from scipy.integrate import solve_ivp
+import matplotlib.pyplot as plt
+
+
+
+#Part one Diffeq
+def f_part1(t, state, Me, Mm, G):
+    xm, ym, vx, vy = state
+
+    dvxdt = -(Me*G*xm)/((xm**2+ym**2)**(3/2))
+    dxmdt = vx
+    dvydt = -(Me*G*ym)/((xm**2+ym**2)**(3/2))
+    dymdt = vy
+
+    return (dxmdt, dymdt, dvxdt, dvydt)
+
+
+# Initial Conditions
+t_min = 0
+t_max = 86400*20
+numpoints = 200
+t = np.linspace(t_min, t_max, numpoints)
+
+xm0 = 384400000
+ym0 = 0
+vx0 = 0
+vy0 = 1017.79
+
+Me = 5.97*(10**24)
+Mm = 7.35*(10**22)
+G = 6.67*(10**-11)
+
+
+#Solver
+results = solve_ivp(f_part1, (t_min,t_max), (xm0, ym0, vx0, vy0), args=(Me, Mm, G), t_eval=t)
+
+
+#Graph plotting
+ax=plt.axes()    # This creates some axes, so that we
+ax.set_aspect(1) # can set the aspect ratio to 1 i.e. 
+                 # x and y axes are scaled equally.                                
+ax.set_xlabel("x coordinate (m)") # Must label axes (with
+ax.set_ylabel("y coordinate (m)") # units) and give
+ax.set_title("Orbit of moon around earth") # plot title.
+#
+ax.plot(results.y[0],results.y[1]) # Make the plot
+ax.legend(['Moon'])                # and add a key.
+plt.show()
+
+
+# MyInput = '0'
+# while MyInput != 'q':
+#     MyInput = input('Enter a choice, "1", "2" or "q" to quit: ')
+#     print('You entered the choice: ',MyInput)
+#     if MyInput == '1':
+#         print('You have chosen part (1): simulation of a lunar orbit')
+#         #
+#         # put your code for part (1) here
+#         #
+#     elif MyInput == '2':
+#         print('You have chosen part (2): earth-moon-probe system')
+#         #
+#         # put your code for part (2) here
+#         #
+#     elif MyInput != 'q':
+#         print('This is not a valid choice')
+# print('You have chosen to finish - goodbye.')