教学相长。舍入误差对计算精度的影响,一般计算方法的教材里都有讲,但是里面的例子总感觉不够生动。
万没想到的是,前几天在课上,博主遇到了一个意想不到的情况,生动地展示了舍入误差的影响。
考虑方程 x = 4 sinx。这个方程在x=2.5附近有一个非平庸的零点。这是函数4sinx的一个不稳定的不动点。我们用aitken方法来逼近它。
我们有两个解析上等价的迭代公式。但是数值上,这两个表达式会导致截然不同的结果。
clear all; close all; clc;
format long;
a = 4;
N = 6;
disp("1st method")
xlist = zeros(1, N);
x = 2.5;
for s = 1 : N
x1 = a * sin(x);
x2 = a * sin(x1);
x = x2 - (x2-x1)^2/(x2 + x - 2*x1);
xlist(s) = x;
disp(x)
end
disp("2nd method")
xlist = zeros(1, N);
x = 2.5;
for s = 1 : N
x1 = a * sin(x);
x2 = a * sin(x1);
x = (x1^2 - x*x2)/(2*x1 - x - x2);
xlist(s) = x;
disp(x)
end
在第二个办法中,分子里的两项都是1的量级,而我们需要把他们的精度算到epsilon的平方的量级,才能实现aitken方法的精度。在一般的双精度计算机上,这意味着当误差epsilon达到10的负8次方的量级后,误差就无法再减小了。这时候,对计算机而言,(1+epsilon)^2 = 1 + 2 epsilon,平方项被舍去了。
第一个方法则不存在这个问题,误差可以一直减小到机器精度,即10的负16次方。
输出结果:
1st method
2.473939114711634
2.474576406173987
2.474576787369692
2.474576787369829
2.474576787369829
2.474576787369829
2nd method
2.473939114711633
2.474576406174009
2.474576787257246
2.474576839952956
2.474576787532312
2.474576956802725