差分スキームの誤差と情報落ち [件名訂正]

差分スキームの誤差は O(Δx) のようにステップの冪であらわされる。つまり、アルゴリズムをなにも改良しなくても、ステップを小さくすれば誤差を小さくできる。しかしそれは無限精度の実数で計算ができる場合で、現実の浮動小数点演算においては別種の誤差もあるので無限高精度の積分を達成することはできない。ステップが小さいのに漫然と加算をしていれば桁落ち[訂正：情報落ち]がまず思いつく脅威である。

あまりいい例ではないかもしれないが、わかりやすいので、前進差分で dy/dx = 2x
を x: [0, 1] で積分してみた。もちろん解析解は ｙ = x^2 で答は 1 である。

出力をみると、ステップ 1e-8 程度で精度は頭打ちになっている。

import java.text.DecimalFormat;

public strictfp class Z {
public static void main(String[] args) {
DecimalFormat df = new DecimalFormat(" 0.0;-0.0");
df.setMinimumIntegerDigits(1);
df.setMaximumIntegerDigits(1);
df.setMinimumFractionDigits(16);
df.setMaximumFractionDigits(16);
double dxs[] = {1.0e-4, 3.0e-5, 1.0e-5, 3.0e-6, 1.0e-6, 3.0e-7,
1.0e-7, 3.0e-8, 1.0e-8, 3.0e-9, 1.0e-9};
double xmin = 0.0;
double xmax = 1.0;
double ymin = 0.0;
double x, y, dydx;
boolean continuep;
for (double step: dxs) {
double dx = step;
double xstop = xmax - dx;
for (x = xmin, y = ymin, continuep = true; continuep; x += dx) {
dydx = x * 2.0;
if (x > xstop) {
dx = xmax - x;
continuep = false;
}
y += dx * dydx;
}
System.out.println(
df.format(step)
+ " " + df.format(y)
+ " " + df.format(y - (xmax * xmax - xmin * xmin))
);
}
}
}

/* output

0.0001000000000000 0.9999000000001081 -0.0000999999998919
0.0000300000000000 0.9999700002007929 -0.0000299997992071
0.0000100000000000 0.9999900000031827 -0.0000099999968173
0.0000030000000000 0.9999970000096747 -0.0000029999903253
0.0000010000000000 0.9999989999833878 -0.0000010000166122
0.0000003000000000 0.9999997001106760 -0.0000002998893240
0.0000001000000000 0.9999999003876101 -0.0000000996123899
0.0000000300000000 0.9999999690878816 -0.0000000309121184
0.0000000100000000 0.9999999863603514 -0.0000000136396486
0.0000000030000000 0.9999999888525976 -0.0000000111474024
0.0000000010000000 1.0000000151388213 0.0000000151388213

*/