洛伦兹方程求解与可视化分析

1. 洛伦兹方程为： $$ \begin{cases} \frac{dx}{dt} = \sigma(y-x) \\ \frac{dy}{dt} = x(\rho-z)-y \\ \frac{dz}{dt} = xy-\beta z \end{cases} $$ 其中，$\sigma, \rho, \beta$ 为常数。 雅可比矩阵为： $$ J = \begin{bmatrix} \frac{\partial f_1}{\partial x} & \frac{\partial f_1}{\partial y} & \frac{\partial f_1}{\partial z} \\ \frac{\partial f_2}{\partial x} & \frac{\partial f_2}{\partial y} & \frac{\partial f_2}{\partial z} \\ \frac{\partial f_3}{\partial x} & \frac{\partial f_3}{\partial y} & \frac{\partial f_3}{\partial z} \end{bmatrix} $$ 其中，$f_1 = \sigma(y-x), f_2 = x(\rho-z)-y, f_3 = xy-\beta z$。 计算过程如下： $$ J = \begin{bmatrix} -\sigma & \sigma & 0 \\ \rho-z & -1 & -x \\ y & x & -\beta \end{bmatrix} $$ 2. 使用Matlab求解洛伦兹方程的代码如下： ```matlab sigma = 10; rho = 28; beta = 8/3; dt = 0.01; t = 0:dt:50; x = zeros(size(t)); y = zeros(size(t)); z = zeros(size(t)); x(1) = 1; y(1) = 1; z(1) = 1; for i = 2:length(t) x(i) = x(i-1) + sigma*(y(i-1)-x(i-1))*dt; y(i) = y(i-1) + (x(i-1)*(rho-z(i-1))-y(i-1))*dt; z(i) = z(i-1) + (x(i-1)*y(i-1)-beta*z(i-1))*dt; end plot(x,y) xlabel('x') ylabel('y') title('Lorenz Attractor') ``` 选取相对瑞利数为28，初始状态为$(1,1,1)$，得到的轨迹如下图所示：  3. 使用Matlab求解洛伦兹方程的代码如下： ```matlab sigma = 10; rho = 28; beta = 8/3; dt = 0.01; t = 0:dt:50; x1 = zeros(size(t)); y1 = zeros(size(t)); z1 = zeros(size(t)); x2 = zeros(size(t)); y2 = zeros(size(t)); z2 = zeros(size(t)); x1(1) = 1; y1(1) = 1; z1(1) = 1; x2(1) = 1.01; y2(1) = 1.01; z2(1) = 1.01; for i = 2:length(t) x1(i) = x1(i-1) + sigma*(y1(i-1)-x1(i-1))*dt; y1(i) = y1(i-1) + (x1(i-1)*(rho-z1(i-1))-y1(i-1))*dt; z1(i) = z1(i-1) + (x1(i-1)*y1(i-1)-beta*z1(i-1))*dt; x2(i) = x2(i-1) + sigma*(y2(i-1)-x2(i-1))*dt; y2(i) = y2(i-1) + (x2(i-1)*(rho-z2(i-1))-y2(i-1))*dt; z2(i) = z2(i-1) + (x2(i-1)*y2(i-1)-beta*z2(i-1))*dt; end plot(t,x1,'b',t,x2,'r') xlabel('t') ylabel('x') title('Lorenz Attractor') legend('Initial state: (1,1,1)','Initial state: (1.01,1.01,1.01)') ``` 选取相对瑞利数为28，初始状态分别为$(1,1,1)$和$(1.01,1.01,1.01)$，得到的轨迹如下图所示：