MATLAB 显式欧拉法求解常微分方程初值问题及稳定区间分析

使用 MATLAB 编程利用显式 Euler 法求解常微分方程初值问题

本文将使用 MATLAB 编程实现显式欧拉法，并用它来求解如下常微分方程组的初值问题：

du/dt = -4u + 6v, u(0) = 1,dv/dt = 3u - 7v, v(0) = 2

同时，我们还会分析时间步长的稳定区间。

MATLAB 代码matlab% 设置初始条件u(1) = 1;v(1) = 2;t(1) = 0;

% 设置步长和总时间dt = 0.01; % 此处取了一个较小的步长T = 10;

% 利用显式 Euler 法求解for i = 1:T/dt u(i+1) = u(i) + dt*(-4u(i) + 6v(i)); v(i+1) = v(i) + dt*(3u(i) - 7v(i)); t(i+1) = t(i) + dt;end

% 绘制图像figure;subplot(2,1,1);plot(t, u, 'b');xlabel('t');ylabel('u');title('Solution of du/dt = -4u + 6v');subplot(2,1,2);plot(t, v, 'r');xlabel('t');ylabel('v');title('Solution of dv/dt = 3u - 7v');

% 计算稳定区间lambda1 = -4;lambda2 = -7;stable_dt = 2/abs(lambda1+lambda2);fprintf('The stable time step is %f. ', stable_dt);

运行结果

其中，蓝色线表示 $u$ 的解，红色线表示 $v$ 的解。可以看出，两个方程的解都是指数衰减的。

稳定区间分析

通过计算，可以得到时间步长的稳定区间为：

$$\Delta t \leq \frac{2}{|\lambda_1 + \lambda_2|} = \frac{2}{11} \approx 0.1818$$

也就是说，只有当步长小于等于 $0.1818$ 时，显式 Euler 法才能稳定地求解这个问题。如果步长过大，可能会导致数值解发散或者出现非物理解。

总结

本文使用 MATLAB 编程实现显式欧拉法求解常微分方程组的初值问题，并分析时间步长的稳定区间。代码示例和结果可视化，可以帮助读者更直观地理解该方法的应用和注意事项。