atan2f対fmodf対公正な減算

Question

私は統合中に角度を囲むように書いたコードに問題があり、私が取り組んでいる小さなシミュレーションの一部です。ですから、基本的には、角度が大きくなるのを防ぐには、それが常に正当な価値を持っていることを確認することです。私は3つの異なるアプローチを試みましたが、私は同じ結果を期待しています。そして、ほとんどの時間彼らはします。しかし、最初の2つは、アングル値がラップアラウンドするポイントの周りにアーチファクトを与えます。私が角度値から波形を生成すると、これらの精度誤差のために望ましくない結果が得られます。

だから、最初のアプローチは、この（-8PI + 8PI範囲に限界角度）のようなものです：

：の

第二のアプローチ：これは、このように見えるアーティファクトを作成

self->state.angle = atan2f(sinf(angle/8), cosf(angle/8)) * 8; 

self->state.angle = fmodf(angle, (float)(2.f * M_PI * 8)) 

同じ結果を作成します。

は、しかし、私はちょうどこのようにそれを行う場合：

float limit = (8 * 2 * M_PI); 
if(angle > limit) angle -= limit;   
if(angle < 0) angle += limit;    
self->state.angle = a; 

予想通り、それはどのアーティファクトなしで動作します。

だから私はここで何をしないのですか？なぜ他の2つの方法で精度誤差が生じるのでしょうか？私はそれらのすべてが同じ結果を生み出すことを期待するでしょう（私は角度の範囲が異なっていることを知っていますが、角度がさらにsin関数に渡されると、結果は同じになると思います）。

編集：小さなテスト

// g++ -o test test.cc -lm && ./test 
#include <stdio.h> 
#include <stdlib.h> 
#include <math.h> 
#include <stdint.h> 

int main(int argc, char **argv){ 
    float a1 = 0; 
    float a2 = 0; 
    float a3 = 0; 
    float dt = 1.f/7500.f; 

    for(float t = -4.f * M_PI; t < (4.f * M_PI); t+=dt){ 
     a1 += dt; 
     a2 += dt; 
     a3 += dt; 

     float b1 = a1; 
     if(b1 > 2.f * M_PI) b1 -= 2.f * M_PI; 
     if(b1 < 0.f) b1 += 2.f * M_PI; 
     float b2 = atan2f(sinf(a2), cosf(a2)); 
     float b3 = fmodf(a3, 2 * M_PI); 

     float x1 = sinf(b1); 
     float x2 = sinf(b2); 
     float x3 = sinf(b3); 

     if((x1 * x2 * x3) > 1e-9){ 
      printf("%f: x[%f %f %f],\tx1-x2:%f x1-x3:%f x2-x3:%f]\n", t, x1, x2, x3, (x1 - x2) * 1e9, (x1 - x3) * 1e9, (x2 - x3) * 1e9); 
     } 
    } 

    return 0; 
} 

出力：

-9.421306: x[0.001565 0.001565 0.001565],  x1-x2:0.000000 x1-x3:0.000000 x2-x3:0.000000] 
-9.421172: x[0.001431 0.001431 0.001431],  x1-x2:0.000000 x1-x3:0.000000 x2-x3:0.000000] 
-9.421039: x[0.001298 0.001298 0.001298],  x1-x2:0.000000 x1-x3:0.000000 x2-x3:0.000000] 
-9.420905: x[0.001165 0.001165 0.001165],  x1-x2:0.000000 x1-x3:0.000000 x2-x3:0.000000] 
-9.420772: x[0.001032 0.001032 0.001032],  x1-x2:0.000000 x1-x3:0.000000 x2-x3:0.000000] 
-6.275573: x[0.001037 0.001037 0.001037],  x1-x2:0.000000 x1-x3:174.855813 x2-x3:174.855813] 
-6.275439: x[0.001171 0.001171 0.001171],  x1-x2:0.000000 x1-x3:174.855813 x2-x3:174.855813] 
-6.275306: x[0.001304 0.001304 0.001304],  x1-x2:0.000000 x1-x3:174.855813 x2-x3:174.855813] 
-6.275172: x[0.001438 0.001438 0.001438],  x1-x2:0.000000 x1-x3:174.855813 x2-x3:174.855813] 
-6.275039: x[0.001571 0.001571 0.001571],  x1-x2:0.000000 x1-x3:174.855813 x2-x3:174.855813] 
-6.274905: x[0.001705 0.001705 0.001705],  x1-x2:0.000000 x1-x3:174.855813 x2-x3:174.855813] 
-6.274772: x[0.001838 0.001838 0.001838],  x1-x2:0.116415 x1-x3:174.855813 x2-x3:174.739398] 

Answer 1

より多くの情報がなければ、それは説明を提供するのは困難だが、私はとにかく頑張ります。

fmodと「プレーン・サブトラクション」（または追加）を使用することの違いは、値がすでに範囲外である場合（たとえば、800000 * M_PIなど）、加算/減算方法ではアーチファクトが見られないので、値を大きく変更する（ほとんど効果がありません）と非常に大きな（絶対値で）アングルが問題なく計算機能に当たっています。

fmod（またはatan2）を使用すると、値が定義した範囲内にあることが保証されますが、これは同じではありません。やっていること

注：

float limit = (8 * 2 * M_PI); 
while(angle > limit) angle -= limit;   
while(angle < 0) angle += limit;    
self->state.angle = a; 

がfmodと同等（約）だろう（しかし、それが原因で繰り返し加算または減算の浮動小数点累積誤差を導入するので、大きな値のためのよりも悪いfmodになり）。

計算に非常に大きな値を入力すると正しい結果が得られる場合は、数値を数学ライブラリに残すのではなく、角度を正規化することをお勧めします。

EDIT：...

他の違いこの答えの最初の部分は、この超範囲外のケースが発生すると仮定し、さらに、問題の編集が、これはそうではなかったことを示したので、 fmod及び2の試験の間に実装点の不正確が元の値に小さな値をsubstractできるフローティング、value - int(value/modulus)*modulus;ようなものであれば例えばfmod

を呼び出すときに値が範囲内に既にいる場合、同じであることを保証はありませんということです。

atan2fとsinを組み合わせて使用​​すると、すでに範囲内にある場合も結果が変更されます。

（値が範囲外のわずかであっても、そして、あなたがやっているように下塗り/追加することが切り捨て/乗算/除算含まない）

あなただけ追加することにより、範囲内の値を調整することができますので、あなたのケースでは、fmodfまたはatan2fを使用して過剰な場合は、あなたの単純なサブ/追加に固執することができます（elseを追加するとテストが保存されます：あなたはあまりにも低い値をreajustedすれば、値が大きすぎる）

Answer 2

float対double数学。

もちろん、3番目の方法が最適です。それはdouble数学を使用しています。

---

b1, b3をご覧ください。 b3はfmodf()コールのために確かにfloatの精度で計算されます。

M_PIは、通常のでb1 -= 2.f * M_PI;がありそうdouble精度の数学で行われdoubleで、より正確な答えを提供して注意してください。 2.fのfは、製品を2.f * M_PIにfloatに強制しません - 製品はdoubleであり、したがって-=です。

b1 -= 2.f * M_PI; 
// same as 
b1 = (float)((double)b1 - (2.f * M_PI)); 

はさらに：最適化とFLT_EVAL_METHOD > 0で、Cは、型の精度よりも高いでFPコードの実行が許可されています。コードがfloatと表示されても、b1はdoubleで計算されます。より精度の高い、とM_PI（有理数）は、正確にπ（無理数）ではないという事実、公正を行うためにvolatile float b1 = a1;を使用し、float結果を保証するために、より正確なb1

fmodf(a3, 2 * M_PI);より

float b1 = a1; 
if(b1 > 2.f * M_PI) b1 -= 2.f * M_PI; // double math 
if(b1 < 0.f) b1 += 2.f * M_PI;   // double math 
float b3 = fmodf(a3, 2 * M_PI); 

につながります比較して#define M_PIf ((float) M_PI)のような定数float

さらに、公正な比較を行うと、より使いやすくなりますif(b1 < -2.f * M_PIf) b1 += 2.f * M_PIf;

OPプリントFLT_EVAL_METHODを推奨します。

1. 敏感ラジアン削減のためdoubleのように広い数学を使用します。

   #include <float.h> 
   printf("%d\n", FLT_EVAL_METHOD); 
   

   ---

   OPは2つのソリューションを提供しています。

   float b3 = fmod(a3, 2 * M_PI); // not fmodf 
   

2. ラジアンを使用していますが、角度などの測定やBAMと正確な範囲の削減を行わないでください。角度は、三角呼び出しの前にラジアンの変換に角度を必要とします。

   float b3 = fmodf(a3, 360.0f); // use fmodf, a3, b3 are in degrees 
   

---

注：float b2 = atan2f(sinf(a2), cosf(a2));方法が合理的候補ではありません。