fused-multiply-add浮動小数点の不正確さを処理する一般的な方法

Question

昨日私は数時間後に多かれ少なかれ次のようなことをしていたコードに絞り込んだプロジェクトのバグを追跡していましたこの：コンパイルと実行後

#include <iostream> 
#include <cmath> 
#include <cassert> 

volatile float r = -0.979541123; 
volatile float alpha = 0.375402451; 

int main() 
{ 
    float sx = r * cosf(alpha); // -0.911326 
    float sy = r * sinf(alpha); // -0.359146 
    float ex = r * cosf(alpha); // -0.911326 
    float ey = r * sinf(alpha); // -0.359146 
    float mx = ex - sx;  // should be 0 
    float my = ey - sy;  // should be 0 
    float distance = sqrtf(mx * mx + my * my) * 57.2958f; // should be 0, gives 1.34925e-06 

// std::cout << "sv: {" << sx << ", " << sy << "}" << std::endl; 
// std::cout << "ev: {" << ex << ", " << ey << "}" << std::endl; 
// std::cout << "mv: {" << mx << ", " << my << "}" << std::endl; 
    std::cout << "distance: " << distance << std::endl; 

    assert(distance == 0.f); 
// assert(sx == ex && sy == ey); 
// assert(mx == 0.f && my == 0.f); 
} 

：

$ g++ -Wall -Wextra -Wshadow -march=native -O2 vfma.cpp && ./a.out 
distance: 1.34925e-06 
a.out: vfma.cpp:23: int main(): Assertion `distance == 0.f' failed. 
Aborted (core dumped) 

ビュー何かの私の視点から、私は（私は2つを取得することが予想2ビット単位、同一ペアの2つの減算を求めてきたように、間違っていますそれらを二乗（2つのゼロ）し、それらを一緒に（ゼロ）加算することによって、

問題の根本的な原因は、fused-multiply-add演算の使用であることが判明しました。これは、行のどこかで結果が不正確になります（私の観点から）。一般的に私はこの最適化に対して何も持っていません。より正確な結果を出すことが約束されていますので、と正確ですが、この場合1.34925e-06は私が期待していた0から本当に遠いです。

テストケースは非常に「壊れやすい」 - それ以上のプリントを有効にすると、それ以上アサートすると、コンパイラはfused-multiply-addを使用しないためアサートを停止します。

私は、コンパイラにかかわる問題であることを考えられてきたように

$ g++ -Wall -Wextra -Wshadow -march=native -O2 vfma.cpp && ./a.out 
sv: {-0.911326, -0.359146} 
ev: {-0.911326, -0.359146} 
mv: {0, 0} 
distance: 0 

は、私がいることを報告してきたが、これは正しい動作であることを説明して閉じてしまった：たとえば、私は、すべての行のコメントを解除した場合。

https://gcc.gnu.org/bugzilla/show_bug.cgi?id=79436

だから私は思ったんだけど - べき一つのコードこのような計算は、問題を回避する方法？私は、一般的な解決策を考えていたが、より良い何か：

mx = ex != sx ? ex - sx : 0.f; 

は私が修正したり、自分のコードを改善したいと思います - fused-として、代わりに私のプロジェクト全体のために-ffp-contract=offを設定する - 改善/修正するために何があるかどうコンパイラライブラリの内部でmultiply-addが使用されています（sinf（）とcosf（）で多く参照されています）、解決策ではなく「部分的な回避策」になります。 （「浮動小数点を使用しないでください」のようなソリューションを避けるために、

Answer 1

をなし、一般的には：これはまさにあなたが偶然（-ffp-contract=fastを使用するために支払う価格で、それはまさにこの例でWilliam Kahan notes in the problems with automatic contraction）

その

理論的には、C（C++ではなく）を使用していて、コンパイラがC-1999プラグマ（つまり、ないGCC）、あなたは

#pragma STDC FP_CONTRACT OFF 
// non-contracted code 
#pragma STDC FP_CONTRACT ON 

Answer 2

を使用することができ興味深いことに、FMAのおかげで、フロートは、MXと私のあなたにrとCOSを掛ける際に行われた丸めエラーが発生します。

fma(r,cos, -r*cos) = theoretical(r*cos) - float(r*cos) 

だから何とか取得結果は、理論的（SX、SY）から計算された（SX、SY）はフロートの乗算によるものであった（ただし、COSの計算における丸め誤差から会計ではないか遠くを示し、罪）。

浮動小数点の丸めに関連する不確実性の範囲内の差異（ex-sx、ey-sy）に、プログラムがどのように依存するのかという疑問がありますか？