3

質問

正確な IEEE 754 算術演算を実装する C99 コンパイラの場合、型の の値がf存在divisorするfloatf / divisor != (float)(f * (1.0 / divisor))?

編集: 「正確な IEEE 754 演算を実装する」とは、FLT_EVAL_METHOD を 0 として正当に定義するコンパイラを意味します。

環境

IEEE 754 準拠の浮動小数点を提供する AC コンパイラは、逆数自体が として正確に表現できる場合にのみ、単精度除算を定数で置き換え、単精度乗算を逆数で置き換えることができfloatます。

実際には、これは 2 のべき乗でのみ発生します。f / 2.0fしたがって、プログラマーのアレックスは、 が のようにコンパイルされると確信しているかもしれませんが、アレックスが 10 で割る代わりにf * 0.5fを掛けることが許容される場合0.10f、アレックスは、プログラムに掛け算を書くか、 GCC の-ffast-math.

この質問は、単精度除算を倍精度乗算に変換することに関するものです。常に正しく丸められた結果を生成しますか? それが安くなる可能性はあり-ffast-mathますか?

反例を見つけることなく、1 から 2 までのすべての単精度値に対してとを比較(float)(f * 0.10)しました。これにより、法線のすべての分割がカバーされ、法線の結果が生成されます。f / 10.0fffloat

次に、以下のプログラムを使用して、テストをすべての除数に一般化しました。

#include <float.h>
#include <math.h>
#include <stdio.h>

int main(void){
  for (float divisor = 1.0; divisor != 2.0; divisor = nextafterf(divisor, 2.0))
    {
      double factor = 1.0 / divisor; // double-precision inverse
      for (float f = 1.0; f != 2.0; f = nextafterf(f, 2.0))
        {
          float cr = f / divisor;
          float opt = f * factor; // double-precision multiplication
          if (cr != opt)
            printf("For divisor=%a, f=%a, f/divisor=%a but (float)(f*factor)=%a\n",
                   divisor, f, cr, opt);
        }
    }
}

検索空間は、これを興味深いものにするのに十分な大きさです (2 46 )。プログラムは現在実行中です。終了する前に、おそらく理由または理由を説明して、何かを印刷するかどうかを誰かに教えてもらえますか?

4

3 に答える 3

2

デフォルト以外の丸めモードが可能な場合、それは確かに不可能です。たとえば、 に置き換える3.0f / 3.0f3.0f * C、正確な逆数よりも小さい値はC、下方またはゼロ方向の丸めモードでC間違った結果をもたらしますが、正確な逆数よりも大きい値は、上方への丸めモードで間違った結果をもたらします。

デフォルトの丸めモードに制限すると、探していることが可能かどうかはわかりません。私はそれについて考えて、何か思いついたらこの答えを修正します。

于 2013-10-25T12:55:13.623 に答える
1

ランダム検索の結果、例が得られました。

結果が「非正規/非正規」の数値の場合、不等式が発生する可能性があります。しかし、私のプラットフォームは IEEE 754 に準拠していないのでしょうか?

f        0x1.7cbff8p-25 
divisor -0x1.839p+116
q       -0x1.f8p-142
q2      -0x1.f6p-142

int MyIsFinite(float f) {
  union {
    float f;
    unsigned char uc[sizeof (float)];
    unsigned long ul;
  } x;
  x.f = f;
  return (x.ul & 0x7F800000L) != 0x7F800000L;
}

float floatRandom() {
  union {
    float f;
    unsigned char uc[sizeof (float)];
  } x;
  do {
    size_t i;
    for (i=0; i<sizeof(x.uc); i++) x.uc[i] = rand();
  } while (!MyIsFinite(x.f));
  return x.f;
}

void testPC() {
  for (;;) {
    volatile float f, divisor, q, qd;
    do {
      f = floatRandom();
      divisor = floatRandom();
      q = f / divisor;
    } while (!MyIsFinite(q));
    qd = (float) (f * (1.0 / divisor));
    if (qd != q) {
      printf("%a %a %a %a\n", f, divisor, q, qd);
      return;
    }
  }

}

Eclipse PC バージョン: Juno Service Release 2 ビルド ID: 20130225-0426

于 2013-10-25T14:49:41.007 に答える