from scipy.integrate import odeint
import numpy as np
import pylab as plt
mu = 1/82.45; lamda= 1-mu
dz = lambda z, t: [z[1], 2z[3]+z[0]-lamda(z[0]+mu)/
((z[0]+mu)2+z[2]2)(3/2)-mu*(z[0]-lamda)/
((z[0]+lamda)2+z[2]2)(3/2), z[3],
-2z[1]+z[2]-lamdaz[2]/((z[0]+mu)2+z[2]2)(3/2)-
mu*z[2]/((z[0]+lamda)2+z[2]2)(3/2)]
t = np.linspace(0, 100, 1001);
s = odeint(dz, [1.2, 0, 0, -1.0494], t)
plt.plot(s[:,0], s[:,2]); plt.show()
2023310143007
