浮動小数点の罠:「a==b」を「fabs(a-b)<E」に書き換えるは間違い。</pp>

a==bで比較してはいけない理由は計算精度都合で「1/10*10 != 1/2*2」となる場合があるからです。

その対策として、微小値E(doulbeだと約1E-15,float 1E-7ぐらい）を使って、

「fabs(a-b)<E」と実装しなおしている人が多いかと思います。

（サラリーマンとしてコード書いていると、a==bってコードは機械的に禁止されてる。）

多くの場合は上のやり方で動くんです。でも以下の例を見て下さい。

fabs( (a-b)*1.0E+100) <E

実はこないだ、ぼーっとコードを書いていてやってしまったんですがこいつは正しく動かない可能性が高いです。だって1.E-15*1.0E+100=1.0E+85ですから。

私の場合は、数学的にはf(x)<g(x)をやりたくてf(x)-g(x)<Eと書いていました。そして同じ問題にぶつかりました。

実際にはf(x)もしくはg(x)のオーダに合わせてEを定数倍する必要があるんですね。

さらに、浮動小数点を「==」だけじゃなく「<」でも考えなきゃいけないんです。

ただ、厄介なのはオーダってどないして調べようという話があります。例えば「a-b == c-d」においてEの定数倍を考える時必要なオーダはaやbのオーダであってa-bのオーダではない。という問題があります。

この問題は、大きな値同士を引き算(or足し算）した結果が元のオーダより小さくなる場合（桁落ち）の時に起きます。

かけ算しかないような場所なら桁落ちの心配はいらないので、単純に以下の様にしてOKです。

fabs(a-b)<E*max(fabs(a),fabs(b))

ここで、再び桁落ちについて考えます。計算途中に桁落ちでおちる桁数をdとすると以下の数式のようにも解釈できます。

fabs(a-b)<E*1.0E+[桁落ちした桁数] *max(fabs(a),fabs(b))

桁落ちした分だけ誤差範囲Eを定数倍して、さらに比較対象のオーダにあわせてEを定数倍する。

上の様に浮動小数点の比較（==だけじゃなく<も)ってのは計算を追わないと解らない根が深い問題なんです。世の中のコーディング規約(MISRA-C) なんかではさくっと無視されたりする問題だけど、すっごい希に見つかるので謎のバグとして頭なやます事になります。気をつけてね。

書いていて自信がなくなったぞ、、、いつもは深い事考えずにEは許される精度の範囲でできるだけ大きくとって逃げてるので。。。orz。だれか詳しい人まともな解説お願い。