The Lorentz system: dx/dt= sigma (y-x) dy/dt = x(rho-z)-y dz/dt =xy-beta z Use R

Question

The Lorentz system: dx/dt= sigma (y-x) dy/dt = x(rho-z)-y dz/dt =xy-beta z Use RK4 to solve the Lorenz equations for (x, y, z) with = 10, = 8/3, = 28 and initial conditions x(0) = 1, y(0) = 1, z(0) = 1, and separately with initial conditions x(0) = 1.01, y(0) = 1.01, z(0) = 1.01, for t [0, 20] with h = 0.01. Use MATLAB’s plot3 function to plot both solution trajectories (in three dimensions), carefully marking the start and end point of each path with appropriate markers.

Thanks The Lorentz system: dx/dt= sigma (y-x) dy/dt = x(rho-z)-y dz/dt =xy-beta z Use RK4 to solve the Lorenz equations for (x, y, z) with = 10, = 8/3, = 28 and initial conditions x(0) = 1, y(0) = 1, z(0) = 1, and separately with initial conditions x(0) = 1.01, y(0) = 1.01, z(0) = 1.01, for t [0, 20] with h = 0.01. Use MATLAB’s plot3 function to plot both solution trajectories (in three dimensions), carefully marking the start and end point of each path with appropriate markers.

Thanks dx/dt= sigma (y-x) dy/dt = x(rho-z)-y dz/dt =xy-beta z Use RK4 to solve the Lorenz equations for (x, y, z) with = 10, = 8/3, = 28 and initial conditions x(0) = 1, y(0) = 1, z(0) = 1, and separately with initial conditions x(0) = 1.01, y(0) = 1.01, z(0) = 1.01, for t [0, 20] with h = 0.01. Use MATLAB’s plot3 function to plot both solution trajectories (in three dimensions), carefully marking the start and end point of each path with appropriate markers.

Thanks

Explanation / Answer

Here is the fully functioning python code. I don't know matlab. Hope this helps. It was hard work!

import numpy as np
from scipy.integrate import odeint
import matplotlib.pyplot as plt
from numpy import pi
import math
import matplotlib.patches as mpatches

#function for odeint to call:
def func(var,t):
sigma, rho, beta = 10.,28.,8./3
x = var[0]
y = var[1]
z = var[2]
x1 = sigma*(y-x)
y1 = x*(rho -z) - y
z1 = x*y - beta*z
return [x1,y1,z1]

#time steps:
t = np.linspace(0.,20.,10000)

#initial conditions:, can be given directly also. That is done for the second

solution.
init = [1.,1.,1.]

#using odeint to solve
sol = odeint(func, init,t)
sol_2 = odeint(func,[1.01,1.01,1.01],t)

#plotting the solution
plt.plot(t,sol[:,0],color='r',linewidth=0.9)

plt.xlabel('Time')
plt.ylabel('x(t)')
plt.show()
plt.plot(t,sol[:,1],color='b',linewidth=0.9)

plt.xlabel('Time')
plt.ylabel('y(t)')
plt.show()
plt.plot(t,sol[:,2],color='g',linewidth=0.9)

plt.xlabel('Time')
plt.ylabel('z(t)')
plt.show()
plt.plot(t,sol[:,0],color='r',linewidth=0.9)
plt.plot(t,sol_2[:,0],color='c',linewidth=0.9)

red_patch = mpatches.Patch(color='red', label='x')
cyan_patch = mpatches.Patch(color='cyan', label='x_new')

plt.legend(handles=[red_patch,cyan_patch])
plt.xlabel('time')
plt.ylabel('x(t)')
plt.show()

plt.plot(t,sol[:,1],color='r',linewidth=0.9)
plt.plot(t,sol_2[:,1],color='c',linewidth=0.9)

red_patch = mpatches.Patch(color='red', label='y')
cyan_patch = mpatches.Patch(color='cyan', label='y_new')

plt.legend(handles=[red_patch,cyan_patch])

plt.xlabel('time')
plt.ylabel('y(t)')
plt.show()

plt.plot(t,sol[:,2],color='r',linewidth=0.9)
plt.plot(t,sol_2[:,2],color='c',linewidth=0.9)
plt.xlabel('time')
plt.ylabel('z(t)')

Related Questions

Navigate

• Browse (All)
• Browse T
• Subjects

---

---

Integrity-first tutoring: explanations and feedback only — we do not complete graded work. Learn more.