Matlab 编程实现数值方法误差分析：步长与精度关系

% 设置初始条件和参数 t0 = 0; tf = 10; y0 = [1; 0]; mu = 1; h = 0.05;

% 定义ODE右侧函数 f = @(t, y) [y(2); -mu*(y(1)^2-1)*y(2)-y(1)];

% 使用不同步长进行计算 for i = 1:3 h = h/2; [t, y] = euler(f, t0, tf, y0, h); [t_rk4, y_rk4] = rk4(f, t0, tf, y0, h); [t_im, y_im] = imp_euler(f, t0, tf, y0, h); [t_bdf, y_bdf] = bdf(f, t0, tf, y0, h); % 计算整体误差 err_euler(i) = norm(y(:, end) - y_rk4(:, end)); err_rk4(i) = norm(y_rk4(:, end) - y_rk4(:, end)); err_im(i) = norm(y_im(:, end) - y_rk4(:, end)); err_bdf(i) = norm(y_bdf(:, end) - y_rk4(:, end)); end

% 画图进行比较 log_h = log([0.05, 0.025, 0.0125]); log_err_euler = log(err_euler); log_err_rk4 = log(err_rk4); log_err_im = log(err_im); log_err_bdf = log(err_bdf);

figure hold on plot(log_h, log_err_euler, '-o'); plot(log_h, log_err_rk4, '-o'); plot(log_h, log_err_im, '-o'); plot(log_h, log_err_bdf, '-o'); legend('Euler', 'RK4', 'Implicit Euler', 'BDF'); xlabel('log(h)'); ylabel('log(enum)');