基于matlab的欧拉法，改进欧拉法，显示欧拉法，隐式你欧拉法，四阶龙格库塔法求解微分方程
微分方程这玩意儿，碰到手算真是要命。但有了MATLAB，咱们就能用数值方法暴力破解各种难解的方程。今天主要聊聊几个经典方法：从最基础的欧拉法到四阶龙格库塔，手把手教你用代码把它们玩出花。
先拿显式欧拉法开刀。这方法简单粗暴，核心思想就是用当前点的斜率推下一步。比如解方程y'=f(x,y)，代码骨架大概长这样：
function [x, y] = explicit_euler(f, x0, y0, h, N) x = zeros(1,N+1); y = zeros(1,N+1); x(1) = x0; y(1) = y0; for i=1:N y(i+1) = y(i) + h*f(x(i), y(i)); % 关键在这行，直接拿当前斜率往前怼 x(i+1) = x(i) + h; end end比如解y'=y这种简单方程，你会发现步长h选大了结果直接起飞。这就是显式欧拉的毛病——稳定性差，得小心控制步长。
那隐式欧拉咋整？它得用未来斜率算当前值，方程得反着解。举个y'= -y的例子：
function [x, y] = implicit_euler(f, x0, y0, h, N) x = x0 + h*(0:N); y = zeros(1,N+1); y(1) = y0; for n=1:N % 牛顿迭代解方程 y_{n+1} = y_n + h*f(x_{n+1}, y_{n+1}) yn_guess = y(n); for k=1:10 % 最多迭代10次 residual = yn_guess - y(n) - h*f(x(n+1), yn_guess); if abs(residual) < 1e-6 break; end yn_guess = yn_guess - residual/(1 - h*dfdy(x(n+1), yn_guess)); end y(n+1) = yn_guess; end end这里得注意，隐式方法每步都要解方程，要么迭代要么用符号计算。上面的牛顿迭代就是个典型操作，不过实际用的时候可能要调收敛条件。
基于matlab的欧拉法，改进欧拉法，显示欧拉法，隐式你欧拉法，四阶龙格库塔法求解微分方程
改进欧拉法算是个折中方案，先显式预测再隐式校正。代码结构贼有意思：
function [x, y] = improved_euler(f, x0, y0, h, N) x = x0 + h*(0:N); y = zeros(1,N+1); y(1) = y0; for n=1:N k1 = h*f(x(n), y(n)); % 预测步 k2 = h*f(x(n)+h, y(n)+k1); % 校正步 y(n+1) = y(n) + (k1 + k2)/2; % 取平均 end end这个双斜率策略让精度从O(h)提升到O(h²)，代价只是多算一次函数值。比如解震荡方程时，改进欧法的轨迹明显比显式欧拉更接近真实解。
重头戏还是四阶龙格库塔(RK4)，经典中的经典。直接上硬货：
function [x, y] = RK4(f, x0, y0, h, N) x = x0 + h*(0:N); y = zeros(1,N+1); y(1) = y0; for n=1:N k1 = h*f(x(n), y(n)); k2 = h*f(x(n)+h/2, y(n)+k1/2); k3 = h*f(x(n)+h/2, y(n)+k2/2); k4 = h*f(x(n)+h, y(n)+k3); y(n+1) = y(n) + (k1 + 2*k2 + 2*k3 + k4)/6; end end四个斜率计算看着麻烦，但MATLAB向量化处理起来很爽。比如解刚性问题时，虽然RK4计算量大了四倍，但步长可以放得更开，总体效率可能反而更高。
实际用的时候得注意：对于震荡剧烈或刚度大的方程，可能需要自适应步长。不过这些基础方法已经能解决大部分教学案例了。下次遇到微分方程别慌，挑个合适的算法，代码一怼，图形一画，啥妖魔鬼怪都现形。
)
