matlab - ODE45和Runge-Kutta方法的绝对误差与解析解的比较

如果有人能帮忙解决以下问题，我将不胜感激。
我有以下颂歌：

dr/dt = 4*exp(0.8*t) - 0.5*r   ,r(0)=2, t[0,1]       (1)

我用两种不同的方法解决了（1）。
利用Runge-Kutta方法（四阶）和Matlab中的ode45我将这两个结果与解析解进行了比较，解析解如下：

r(t) = 4/1.3 (exp(0.8*t) - exp(-0.5*t)) + 2*exp(-0.5*t)

当我绘制每种方法相对于精确解的绝对误差时，我得到以下结果：
对于RK方法，我的代码是：

h=1/50;
x = 0:h:1;
y = zeros(1,length(x));
y(1) = 2;
F_xy = @(t,r) 4.*exp(0.8*t) - 0.5*r;
for i=1:(length(x)-1)
    k_1 = F_xy(x(i),y(i));
    k_2 = F_xy(x(i)+0.5*h,y(i)+0.5*h*k_1);
    k_3 = F_xy((x(i)+0.5*h),(y(i)+0.5*h*k_2));
    k_4 = F_xy((x(i)+h),(y(i)+k_3*h));
    y(i+1) = y(i) + (1/6)*(k_1+2*k_2+2*k_3+k_4)*h;  % main equation
end

对于ode45：

tspan = 0:1/50:1;
x0 = 2;
f = @(t,r) 4.*exp(0.8*t) - 0.5*r;
[tid, y_ode45] = ode45(f,tspan,x0);

我的问题是，当我使用ode45时，为什么会有振荡（我指的是绝对误差）两种解决方案都是准确的（1e-9），但在这种情况下ode45会发生什么？
当我计算RK方法的绝对误差时，为什么它看起来更好？

最佳答案

RK4函数正在执行的固定步骤比ode45正在执行的步骤小得多您真正看到的是由于polynomial interpolation导致的错误，该错误用于生成ode45所采取的真正步骤之间的点这通常被称为“密集输出”（参见Hairer & Ostermann 1990）。
当您使用两个以上的元素指定TSPAN向量时，Matlab's ODE suite solvers会产生固定的步长输出但这并不意味着他们实际使用了固定的步长或使用了TSPAN中指定的步长您可以看到实际使用的步长，并且仍然可以通过ode45输出结构和使用deval获得所需的固定步长输出：

sol = ode45(f,tspan,x0);
diff(sol.x) % Actual step sizes used
y_ode45 = deval(sol,tspan);

您将看到，在0.02的初始步骤之后，因为您的ODE很简单，所以在随后的步骤中会收敛到0.1与默认最大步长限制（十分之一积分区间）相结合的默认容差确定了这一点。让我们以正确的步骤绘制错误：

exactsol = @(t)(4/1.3)*(exp(0.8*t)-exp(-0.5*t))+2*exp(-0.5*t);
abs_err_ode45 = abs(exactsol(tspan)-y_ode45);
abs_err_ode45_true = abs(exactsol(sol.x)-sol.y);
abs_err_rk4 = abs(exactsol(tspan)-y);
figure;
plot(tspan,abs_err_ode45,'b',sol.x,abs_err_ode45_true,'k.',tspan,abs_err_rk4,'r--')
legend('ODE45','ODE45 (True Steps)','RK4',2)

如您所见，真实步骤的错误比RK4的错误增长得慢（ode45实际上是一个比RK4更高阶的方法，所以您会预料到这一点）由于插值，积分点之间的误差增大如果要限制此范围，则应通过odeset调整公差或其他选项。
如果你想强制ode45使用1/50的步骤，你可以这样做（因为你的ODE很简单）：

opts = odeset('MaxStep',1/50,'InitialStep',1/50);
sol = ode45(f,tspan,x0,opts);
diff(sol.x)
y_ode45 = deval(sol,tspan);

对于另一个实验，尝试扩大积分间隔，使之积分到t = 10可能您将在错误中看到许多有趣的行为（绘制相对错误在这里很有用）你能解释一下吗你能用ode45和odeset产生表现良好的结果吗将大区间指数函数与自适应步长方法相结合具有挑战性，ode45不一定是这项工作的最佳工具不过，也有alternatives但可能需要一些编程。