Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
x,y,z=var('x,y,z')
# Next we define the parameters
sigma=10
rho=28
beta=8/3
# The Lorenz equations
lorenz=[sigma*(y-x),x*(rho-z)-y,x*y-beta*z]
# Time and initial conditions
N=250000
tmax=250
h=tmax/N
t=srange(0,tmax+h,h)
ics=[0,1,1]
sol=desolve_odeint(lorenz,ics,t,[x,y,z],rtol=1e-13,atol=1e-14)
X=sol[:,0]
Y=sol[:,1]
Z=sol[:,2]
# Plot the result
from mpl_toolkits.mplot3d import axes3d
from matplotlib import pyplot as plt
# Call the plot function if you want to plot the data
fig = plt.figure()
ax = fig.gca(projection='3d')
ax.plot(X, Y, Z, lw=0.5)
ax.set_xlabel("X Axis")
ax.set_ylabel("Y Axis")
ax.set_zlabel("Z Axis")
ax.set_title("Lorenz Attractor")
plt.show()
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.