## This sheet solves an initial value problem of the Lorenz Attractor
## Edward Lorenz explored this function while studying chaos and meteorological modeling
var('X Y Z t')
var('X0 Y0 Z0')
var('A B C')
A = 10
B = 30
C = 8/3
X0 = 10
Y0 = -3
Z0 = 0
## Calculation Parameters
end_points = 50
stepsize = 0.005
steps = end_points/stepsize
ics = [0,X0,Y0,Z0]
## Equations
r1 = (diff(X,t) == A*(Y-X))
r2 = (diff(Y,t) == X*(B-Z)-Y)
r3 = (diff(Z,t) == X*Y-C*Z)
des = [r1.rhs(), r2.rhs(), r3.rhs()]
sol = desolve_system_rk4(des,[X,Y,Z],ics,ivar=t,end_points=end_points,step=stepsize)
## Clean up to graph
sol_t=[]
sol_x=[]
sol_y=[]
sol_z=[]
for i in range(steps):
sol_t.append(sol[i][0])
sol_x.append(sol[i][1])
sol_y.append(sol[i][2])
sol_z.append(sol[i][3])
show(line(zip(sol_x,sol_y,sol_z)))