これはバグだとしか思えません。最初のアサートは成功し、2 番目のアサートは失敗します。
double sum_1 = 4.0 + 6.3;
assert(sum_1 == 4.0 + 6.3);
double t1 = 4.0, t2 = 6.3;
double sum_2 = t1 + t2;
assert(sum_2 == t1 + t2);
バグでない場合、なぜですか?
浮動小数点数を比較しています。そうしないでください。浮動小数点数には、状況によっては固有の精度エラーがあります。代わりに、2 つの値の差の絶対値を取得し、その値が小さい数 (イプシロン) より小さいことをアサートします。
void CompareFloats( double d1, double d2, double epsilon )
{
assert( abs( d1 - d2 ) < epsilon );
}
これはコンパイラとは関係なく、すべて浮動小数点数の実装方法に関係しています。ここにIEEE仕様があります:
http://www.eecs.berkeley.edu/~wkahan/ieee754status/IEEE754.PDF
これは私も噛んだものです。
はい、丸め誤差のため、浮動小数点数は等しいかどうかを比較するべきではありません。おそらくそれはご存知でしょう。
しかし、この場合、 をt1+t2
計算してから、もう一度計算します。確かにそれは同じ結果を生み出さなければなりませんか?
これがおそらく起こっていることです。これを x86 CPU で実行していると思いますよね?x86 FPU は内部レジスタに 80 ビットを使用しますが、メモリ内の値は 64 ビット double として格納されます。
したがってt1+t2
、最初に 80 ビットの精度で計算され、次に (おそらく) 64 ビットの精度でメモリに格納されsum_2
、いくつかの丸めが行われます。アサートの場合は、浮動小数点レジスタに再度ロードされ、t1+t2
80 ビットの精度で再度計算されます。sum_2
ここで、以前に 64 ビット浮動小数点値に丸められた を、より高い精度 (80 ビット) で計算されたと比較しています。これt1+t2
が、値が完全に同一でない理由です。
編集では、なぜ最初のテストがパスするのでしょうか? この場合、コンパイラはおそらく4.0+6.3
コンパイル時に評価し、代入とアサートの両方で 64 ビット量として格納します。したがって、同一の値が比較され、アサートはパスします。
2 番目の編集コードの 2 番目の部分 (gcc、x86) 用に生成されたアセンブリ コードとコメントを次に示します。上記のシナリオにほぼ従っています。
// t1 = 4.0
fldl LC3
fstpl -16(%ebp)
// t2 = 6.3
fldl LC4
fstpl -24(%ebp)
// sum_2 = t1+t2
fldl -16(%ebp)
faddl -24(%ebp)
fstpl -32(%ebp)
// Compute t1+t2 again
fldl -16(%ebp)
faddl -24(%ebp)
// Load sum_2 from memory and compare
fldl -32(%ebp)
fxch %st(1)
fucompp
興味深い補足: これは最適化なしでコンパイルされました。でコンパイルすると、コンパイラはすべてのコードを-O3
最適化します。
浮動小数点数の近さを比較するときは、通常、それらの相対的な差を測定する必要があります。これは、次のように定義されます。
if (abs(x) != 0 || abs(y) != 0)
rel_diff (x, y) = abs((x - y) / max(abs(x),abs(y))
else
rel_diff(x,y) = max(abs(x),abs(y))
例えば、
rel_diff(1.12345, 1.12367) = 0.000195787019
rel_diff(112345.0, 112367.0) = 0.000195787019
rel_diff(112345E100, 112367E100) = 0.000195787019
アイデアは、数字に共通する有効数字の先頭の数を測定することです。0.000195787019の-log10を取得すると、3.70821611が得られます。これは、すべての例に共通する10進数の先頭の数です。
2つの浮動小数点数が等しいかどうかを判断する必要がある場合は、次のようにする必要があります。
if (rel_diff(x,y) < error_factor * machine_epsilon()) then
print "equal\n";
ここで、マシンイプシロンは、使用されている浮動小数点ハードウェアの仮数に保持できる最小の数値です。ほとんどのコンピューター言語には、この値を取得するための関数呼び出しがあります。error_factorは、数値xおよびyの計算で、丸め誤差(およびその他)によって消費されると思われる有効桁数に基づいている必要があります。たとえば、xとyが約1000の合計の結果であり、合計される数値の範囲がわからない場合、error_factorを約100に設定します。
これらをリンクとして追加しようとしましたが、これが私の最初の投稿であるため、追加できませんでした。
Intel Core 2 Duo で問題を再現し、アセンブリ コードを確認しました。これが何が起こっているかです: コンパイラが を評価するときt1 + t2
、それは行います
load t1 into an 80-bit register
load t2 into an 80-bit register
compute the 80-bit sum
中に収納するsum_2
と
round the 80-bit sum to a 64-bit number and store it
次に、==
比較では 80 ビットの合計と 64 ビットの合計が比較されます。主に、小数部 0.3 は 2 進浮動小数点数を使用して正確に表すことができないため、これらは異なります。 2 つの異なる長さに切り捨てられたバイナリ)。
本当に腹立たしいのは、gcc -O1
またはgcc -O2
を指定してコンパイルすると、コンパイル時に gcc が間違った計算を行い、問題が解消されることです。標準ではこれで問題ないのかもしれませんが、gcc が私のお気に入りのコンパイラではないもう 1 つの理由です。
PS 80 ビットの合計と 64 ビットの合計を比較すると言うとき==
、もちろん、64 ビットの合計の拡張バージョンを比較することを意味します。考えてもいいかもしれません
sum_2 == t1 + t2
に解決
extend(sum_2) == extend(t1) + extend(t2)
と
sum_2 = t1 + t2
に解決
sum_2 = round(extend(t1) + extend(t2))
素晴らしい浮動小数点の世界へようこそ!
場合によっては、64 ビットの double を 80 ビットの内部レジスタと比較することになるかもしれません。GCCが2つのケースに対して発行するアセンブリ命令を見ると啓発されるかもしれません...
倍精度数の比較は本質的に不正確です。たとえば、false0.0 == 0.0
を返すことがよくあります。これは、FPU が数値を保存および追跡する方法によるものです。
ウィキペディアは次のように述べています。
等しいかどうかのテストには問題があります。数学的に等しい 2 つの計算シーケンスは、異なる浮動小数点値を生成する可能性があります。
正確な値ではなく、デルタを使用して比較の許容範囲を指定する必要があります。
この「問題」は、次のオプションを使用して「修正」できます。
-msse2 -mfpmath=sse
このページで説明されているように:
http://www.network-theory.co.uk/docs/gccintro/gccintro_70.html
これらのオプションを使用すると、両方のアサートが成功しました。