12

ある浮動小数点数が別の浮動小数点数の逆数であるかどうかを（C ++で）判別したいと思います。問題は、それを行うために3番目の変数を使用する必要があることです。たとえば、このコード：

float x=5,y=0.2;
if(x==(1/y)) cout<<"They are the multiplicative inverse of eachother"<<endl;
else cout<<"They are NOT the multiplicative inverse of eachother"<<endl;

次のように出力されます：「彼らはそうではありません...」これは間違っており、このコードは次のとおりです。

float x=5,y=0.2,z;
z=1/y;
if(x==z) cout<<"They are the multiplicative inverse of eachother"<<endl;
else cout<<"They are NOT the multiplicative inverse of eachother"<<endl;

「彼らは...」と出力されます。これは正しいです。
なぜこうなった？

4

5 に答える 5

36

浮動小数点精度の問題

    ここには2つの問題がありますが、どちらも同じルートから来ています

フロートを正確に比較することはできません。それらを正確に減算または除算することはできません。あなたは彼らのために何も正確に数えることはできません。それらを使用した操作は、結果にエラーをもたらす可能性があります（ほとんどの場合はそうなります）。a=0.2f正確な操作でもありません。そのより深い理由は、ここにある他の回答の著者によって非常によく説明されています。（私の感謝と彼らへの投票。）

ここに、最初のより単純なエラーがあります。決して、決して、決して、決して、決して、それらに==または同等の言語を使用しないでください。

の代わりにa==b、代わりに使用Abs(a-b)<HighestPossibleErrorしてください。

    しかし、これはあなたの仕事の唯一の問題ではありません。

Abs(1/y-x)<HighestPossibleError どちらも機能しません。少なくとも、それは十分な頻度で機能しません。なんで？

x=1000とy=0.001のペアを考えてみましょう。yの「開始」相対誤差を10-6としましょう。

（相対誤差=誤差/値）。

値の相対誤差は、乗算と除算で増加しています。

1/yは約1000です。その相対誤差は同じ10-6です。（「1」にはエラーはありません）

これにより、絶対誤差= 1000 * 10 -6 =0.001になります。後でxを引くと、そのエラーだけが残ります。（絶対誤差は加算と減算で加算され、xの誤差は無視できるほど小さいです。）確かに、それほど大きな誤差は期待していません。HighestPossibleErrorは確実に低く設定され、プログラムはxの適切なペアをスローします。 y

したがって、フロート操作の次の2つのルールは、大きい方の値を小さい方の値で除算しないようにすることです。そうすれば、神はその後の近い値を減算することからあなたを救います。

この問題を回避する簡単な方法は2つあります。

  • x、yの絶対値が大きいものを見つけ、1を大きい方の値で除算し、後で小さい方の絶対値を減算します。

  • 1/y against x値ではなく文字を使用して作業しているときに比較したい場合で、操作でエラーが発生しない場合は、比較の両側にyを掛けると、が得られます1 against x*y。（通常、その操作で符号を確認する必要がありますが、ここではabs値を使用しているため、クリーンです。）結果の比較にはまったく除算がありません。

簡単に言うと：

1/y V x   <=>   y*(1/y) V x*y   <=>   1 V x*y

1 against x*yそのような比較を行う必要があることはすでにわかっています。

const float HighestPossibleError=1e-10;
if(Abs(x*y-1.0)<HighestPossibleError){...

以上です。

PS本当にすべてを一行で必要とする場合は、以下を使用してください。

if(Abs(x*y-1.0)<1e-10){...

しかし、それは悪いスタイルです。私はそれをお勧めしません。

PPS 2番目の例では、コンパイラーがコードを最適化して、コードを実行する前にzを5に設定します。したがって、5を5に対してチェックすることは、フロートに対しても機能します。

于 2012-02-03T23:33:39.980 に答える
13

問題は0.2、バイナリ展開の桁数が無限であるため、バイナリで正確に表現できないことです。

 1/5: 0.0011001100110011001100110011001100110011...

1/3これは、10進数で正確に表すことができない方法に似ています。xは有限ビット数のに格納されているためfloat、これらの数字はある時点で切り捨てられます。次に例を示します。

   x: 0.0011001100110011001100110011001

この問題は、CPUが内部でより高い精度を使用することが多いために発生します。したがって、計算したばかりの結果はより多くの桁になり、それらを比較するために1/yロードすると、CPUの内部精度に一致するように拡張されます。xx

 1/y: 0.0011001100110011001100110011001100110011001100110011
   x: 0.0011001100110011001100110011001000000000000000000000

したがって、ビットごとに直接比較する場合、それらは異なります。

ただし、2番目の例では、結果を変数に格納すると、比較を行う前に結果が切り捨てられるため、この精度で比較すると、次のようになります。

   x: 0.0011001100110011001100110011001
   z: 0.0011001100110011001100110011001

多くのコンパイラには、一貫性を保つためにすべてのステップで中間値を強制的に切り捨てることができるスイッチがありますが、通常のアドバイスは、浮動小数点値を直接比較することを避け、代わりにイプシロン値よりも小さい値の違いがあるかどうかを確認することです。 Gangnusが示唆していること。

于 2012-02-03T23:52:13.207 に答える
5

2つの近似が逆数になることの意味を正確に定義する必要があります。そうしないと、何をテストすることになっているのかわかりません。

0.2正確なバイナリ表現はありません。精度が制限された正確な表現がない数値を格納すると、正確に正しい答えが得られません。

同じことが10進数でも起こります。たとえば、1/3正確な10進表現はありません。として保存できます.333333。しかし、あなたは問題を抱えています。3と.333333乗法逆数ですか？それらを掛けると、が得られます.999999。答えを「はい」にしたい場合は、乗算して1に等しいかどうかをテストするほど単純ではない、乗法逆数のテストを作成する必要があります。

同じことがバイナリでも起こります。

于 2012-02-03T23:42:14.143 に答える
2

他の返信での議論は素晴らしいので、私はそれらのどれも繰り返さないでしょうが、コードはありません。これは、floatのペアが乗算されたときに正確に1.0になるかどうかを実際にチェックするためのコードです。

コードはいくつかの仮定/アサーションを行います（通常はx86プラットフォームで満たされます）：
- float'sは32ビットバイナリ（AKA single precision）IEEE-754
- int'sまたはlong'sは32ビットです（可用性に依存しないことにしました） of uint32_t）
-8873283.0fmemcpy()が0x4B076543になるようにfloatをint / longにコピーします（つまり、特定の「エンディアン」が期待されます）

追加の前提条件の1つは、次のとおりです。-乗算
される実際のfloatを受け取ります*（つまり、floatの乗算では、数学ハードウェア/ライブラリが内部で使用できるより高い精度の値は使用されません）。 

#include <stdio.h>
#include <string.h>
#include <limits.h>
#include <assert.h>

#define C_ASSERT(expr) extern char CAssertExtern[(expr)?1:-1]

#if UINT_MAX >= 0xFFFFFFFF
typedef unsigned int uint32;
#else
typedef unsigned long uint32;
#endif
typedef unsigned long long uint64;

C_ASSERT(CHAR_BIT == 8);
C_ASSERT(sizeof(uint32) == 4);
C_ASSERT(sizeof(float) == 4);

int ProductIsOne(float f1, float f2)
{
  uint32 m1, m2;
  int e1, e2, s1, s2;
  int e;
  uint64 m;

  // Make sure floats are 32-bit IEE754 and
  // reinterpreted as integers as we expect
  {
    static const float testf = 8873283.0f;
    uint32 testi;
    memcpy(&testi, &testf, sizeof(testf));
    assert(testi == 0x4B076543);
  }

  memcpy(&m1, &f1, sizeof(f1));
  s1 = m1 >= 0x80000000;
  m1 &= 0x7FFFFFFF;
  e1 = m1 >> 23;
  m1 &= 0x7FFFFF;
  if (e1 > 0) m1 |= 0x800000;

  memcpy(&m2, &f2, sizeof(f2));
  s2 = m2 >= 0x80000000;
  m2 &= 0x7FFFFFFF;
  e2 = m2 >> 23;
  m2 &= 0x7FFFFF;
  if (e2 > 0) m2 |= 0x800000;

  if (e1 == 0xFF || e2 == 0xFF || s1 != s2) // Inf, NaN, different signs
    return 0;

  m = (uint64)m1 * m2;

  if (!m || (m & (m - 1))) // not a power of 2
    return 0;

  e = e1 + !e1 - 0x7F - 23 + e2 + !e2 - 0x7F - 23;
  while (m > 1) m >>= 1, e++;

  return e == 0;
}

const float testData[][2] =
{
  { .1f, 10.0f },
  { 0.5f, 2.0f },
  { 0.25f, 2.0f },
  { 4.0f, 0.25f },
  { 0.33333333f, 3.0f },
  { 0.00000762939453125f, 131072.0f }, // 2^-17 * 2^17
  { 1.26765060022822940E30f, 7.88860905221011805E-31f }, // 2^100 * 2^-100
  { 5.87747175411143754E-39f, 1.70141183460469232E38f }, // 2^-127 (denormalized) * 2^127
};

int main(void)
{
  int i;
  for (i = 0; i < sizeof(testData) / sizeof(testData[0]); i++)
    printf("%g * %g %c= 1\n",
           testData[i][0], testData[i][1],
           "!="[ProductIsOne(testData[i][0], testData[i][1])]);
  return 0;
}

出力（ideone.comを参照）：

0.1 * 10 != 1
0.5 * 2 == 1
0.25 * 2 != 1
4 * 0.25 == 1
0.333333 * 3 != 1
7.62939e-06 * 131072 == 1
1.26765e+30 * 7.88861e-31 == 1
5.87747e-39 * 1.70141e+38 == 1
于 2012-02-16T15:14:41.647 に答える
0

驚くべきことは、丸め規則が何であれ、2つのバージョンの結果が同じになることを期待していることです（2回間違っているか2回正しい）！

ほとんどの場合、最初のケースでは、x == 1 / yを評価するときに、FPUレジスタの高精度への昇格が行われますが、z = 1/yは実際には単精度の結果を格納します。

他の寄稿者は、なぜ5 == 1 / 0.2が失敗する可能性があるのか​​を説明していますが、それを繰り返す必要はありません。

于 2012-02-14T11:49:07.627 に答える