Adams 预报校正法
一、直观理解:用“历史斜率”外推,再回炉修正
Euler 与 Runge–Kutta 都是单步法:走下一步只依赖当前点。Adams 方法属于多步法:它记得前几步已经算过的斜率,用多项式去拟合 的历史,再积分到下一步。
预报—校正的两拍节奏是:
- Adams–Bashforth(AB)预报:只用已知历史斜率,显式得到 ;
- Adams–Moulton(AM)校正:把预报点处的新斜率纳入,用隐式公式的显式化版本更新为 。
直观比喻:预测明天天气时,不只看今天,还参考近几天趋势(预报);得到初步预报后,再结合“如果明天是这样,大气方程是否更自洽”做一次修正(校正)。
多步法的优势是:阶数提高后,平均每步调用 的次数可以少于同阶 RK,适合右端昂贵或需要长时间光滑积分的情形。MATLAB 的 ode113 即是变步长、变阶 Adams–Bashforth–Moulton 实现。
二、适用场景与优缺点
| 类型 | 说明 |
|---|---|
| 适用场景 | 解足够光滑、需要长时间积分; 计算成本高,希望减少调用;作为 ode113 的原理说明;与 RK 交叉验证非刚性轨迹。 |
| 核心优点 | 复用历史信息,渐近效率高;预报校正可提供简单误差指示 ;高阶格式对光滑问题精度优秀。 |
| 局限与缺点 | 需要单步法提供启动值;变步长时公式系数不能照搬等步长表;稳定域通常窄于同阶 RK,对突变/刚性更敏感;实现复杂度高于 RK4。 |
使用前应确认:
- 问题非刚性,解在启动段后较光滑;
- 有 RK4/
ode45生成前若干步; - 等步长假设成立,或明确使用变步长算法;
- 准备监控校正量与守恒量。
三、核心原理:多项式积分得到线性多步公式
精确关系:
1. 四步 Adams–Bashforth 预报(显式)
用 做牛顿插值再积分,得
2. 四步 Adams–Moulton 校正(隐式,用预报显式化)
其中 。这是一次校正的 PECE 模式(Predict–Evaluate–Correct–Evaluate 的简化实现中常再算一次最终 供下一步使用)。
3. 阶与误差
上述四步 AB 为 4 阶,配套 AM 校正后整体仍保持高阶(具体实现细节影响误差常数)。局部误差可用
作粗估,供定步长是否过大参考。
4. 启动
步方法需要 。常用同阶或更高阶 RK 生成,不要用粗 Euler 随便填,否则启动误差会污染后续多步过程。
四、关键决策
-
步数(阶)选择
建模手写常用 2–4 步;再高阶对光滑性要求更强,稳定域也可能更苛刻。 -
步长
等步长公式仅在固定 下成立。若必须变步长,优先ode113,避免错误套用系数。 -
校正次数
一次校正通常足够;迭代校正到收敛接近隐式 AM,成本上升。 -
与 RK 分工
短区间、强非线性、频繁事件:RK/ode45;长区间光滑:Adams/ode113。 -
失败信号
校正量持续增大、出现锯齿振荡、与ode45分道扬镳 → 减小 或换方法。
五、MATLAB 示例
%% 四步 Adams 预报校正:y'=-y, y(0)=1
f = @(t, y) -y;
h = 0.05;
t = 0:h:2;
y = zeros(size(t));
% 用精确启动(实践中可换 RK4 启动)
y(1:4) = exp(-t(1:4));
fHist = zeros(1,4);
for j = 1:4
fHist(j) = f(t(j), y(j)); % fHist(1) 对应最旧或需小心索引
end
% 重新按时间顺序存储:fn, fnm1, fnm2, fnm3
% 这里用数组 F(n)=f(t_n,y_n)
F = zeros(size(t));
for n = 1:4
F(n) = f(t(n), y(n));
end
for n = 4:numel(t)-1
% AB4 预报
yP = y(n) + h/24 * (55*F(n) - 59*F(n-1) + 37*F(n-2) - 9*F(n-3));
fP = f(t(n+1), yP);
% AM3/四系数校正(用到 fP 与前三历史)
y(n+1) = y(n) + h/24 * (9*fP + 19*F(n) - 5*F(n-1) + F(n-2));
F(n+1) = f(t(n+1), y(n+1)); % PECE:用校正后的值再评估
end
yTrue = exp(-t);
fprintf("Adams PC 最大误差 = %.3e\n", max(abs(y - yTrue)));
figure;
plot(t, y, "o", t, yTrue, "-"); grid on;
legend("Adams 预报校正", "解析解");
title("四步 Adams 预报校正");
%% 与 ode45 / ode113 对比:轻度振荡问题
fOsc = @(t, z) [z(2); -z(1) - 0.05*z(2)];
z0 = [1; 0];
opts = odeset("RelTol", 1e-7, "AbsTol", 1e-9);
[t1, z1] = ode45(fOsc, [0 40], z0, opts);
[t2, z2] = ode113(fOsc, [0 40], z0, opts);
figure;
plot(t1, z1(:,1), "-", t2, z2(:,1), "--"); grid on;
legend("ode45", "ode113");
xlabel("t"); ylabel("z_1");
title("Adams 型求解器 ode113 与 RK 型 ode45");
%% 手写:RK4 启动 + ABM4
f = @(t, y) -y + sin(t); % y' + y = sin t
tEnd = 10; h = 0.02;
t = 0:h:tEnd;
y = zeros(size(t));
y(1) = 0; % 初值
% RK4 启动前 3 步(得到 y2,y3,y4,共 4 个点)
for n = 1:3
tn = t(n); yn = y(n);
k1 = f(tn, yn);
k2 = f(tn+h/2, yn+h/2*k1);
k3 = f(tn+h/2, yn+h/2*k2);
k4 = f(tn+h, yn+h*k3);
y(n+1) = yn + h/6*(k1+2*k2+2*k3+k4);
end
F = arrayfun(@(n) f(t(n), y(n)), 1:4);
F = F(:).';
Fstore = zeros(size(t));
Fstore(1:4) = F;
for n = 4:numel(t)-1
yP = y(n) + h/24*(55*Fstore(n)-59*Fstore(n-1)+37*Fstore(n-2)-9*Fstore(n-3));
fP = f(t(n+1), yP);
y(n+1) = y(n) + h/24*(9*fP + 19*Fstore(n) - 5*Fstore(n-1) + Fstore(n-2));
Fstore(n+1) = f(t(n+1), y(n+1));
end
% 真解:y = 0.5*(sin t - cos t) + 0.5 e^{-t}
yTrue = 0.5*(sin(t) - cos(t)) + 0.5*exp(-t);
fprintf("RK4 启动 ABM 最大误差 = %.3e\n", max(abs(y - yTrue)));
六、如何解释结果
报告 Adams 预报校正时应写清:
- 多步阶数、等步长 、启动方法(RK4 等);
- 预报公式与校正公式(或说明使用
ode113); - 与解析解或
ode45的误差; - 典型校正量水平,以及是否在某时段异常增大;
- 为何选择多步法(长时间、效率、对照);
- 方法假设:光滑性、非刚性、步长固定(若手写)。
若 ode113 与 ode45 结果高度一致,可增强对非刚性轨迹的信心;若显著不一致,应检查容差、刚性与事件。
七、常见误区
1. 没有可靠启动值
用随意猜测或大步 Euler 填历史,等于从源头注入 误差。
2. 变步长仍套用经典系数表
等步长 AB/AM 系数在 变化时不再正确。
3. 在刚性化学反应上强行高阶 Adams
稳定域不足会导致振荡。应换 BDF/ode15s。
4. 混淆 AB 与 AM 的历史点索引
实现时 off-by-one 非常常见,必须用已知解析例做回归测试。
5. 忽视解的光滑性
解有折点、切换或随机跳跃时,高阶多步可能还不如稳健的低阶单步。
八、小结
Adams 预报校正法的核心思想可以概括为:
用过去若干步的斜率构造多项式,积分得到下一步的显式预报,再把新斜率纳入隐式型公式做一次校正。
它在光滑非刚性的长时间积分中效率出色,也是理解 ode113 的钥匙。建模实践中,手写 ABM 重在公式正确与启动可靠;工程上更常见的是直接调用 ode113,并与 ode45 交叉验证。遇到刚性、强间断或复杂事件,应果断回到单步刚性/非刚性求解器。