19

MATLAB での基本的な算術演算の丸め誤差を理解しようとしていますが、次の興味深い例に出くわしました。

(0.3)^3 == (0.3)*(0.3)*(0.3)

ans = 0

左辺の計算方法を正確に知りたいです。MATLAB のドキュメントでは、整数のべき乗に対して「二乗によるべき乗」アルゴリズムが使用されることが示唆されています。

「行列累乗。p がスカラーの場合、X^p は X の p 乗です。p が整数の場合、累乗は二乗を繰り返すことによって計算されます。」

したがって、同じ値を返すと想定(0.3)^3しました。(0.3)*(0.3)^2しかし、そうではありません。丸め誤差の違いを説明するにはどうすればよいですか?

4

5 に答える 5

5

MATLAB については何も知りませんが、Ruby で試してみました。

irb> 0.3 ** 3
  => 0.026999999999999996
irb> 0.3 * 0.3 * 0.3
  => 0.027

Ruby のソース コードによると、べき乗演算子は、左側のオペランドが float の場合、右側のオペランドを float にキャストしてから、標準 C 関数を呼び出しますpow()。関数のfloatバリアントは、pow()整数以外の指数を処理するためのより複雑なアルゴリズムを実装する必要があります。これにより、丸め誤差が発生する演算が使用されます。おそらくMATLABも同様に機能します。

于 2013-01-23T00:50:14.777 に答える
4

興味深いことに、スカラー^は を使用して実装されているようですがpow、行列^は平方乗算を使用して実装されています。ウィット：

octave:13> format hex
octave:14> 0.3^3
ans = 3f9ba5e353f7ced8
octave:15> 0.3*0.3*0.3
ans = 3f9ba5e353f7ced9
octave:20> [0.3 0;0 0.3]^3
ans =

  3f9ba5e353f7ced9  0000000000000000
  0000000000000000  3f9ba5e353f7ced9

octave:21> [0.3 0;0 0.3] * [0.3 0;0 0.3] * [0.3 0;0 0.3]
ans =

  3f9ba5e353f7ced9  0000000000000000
  0000000000000000  3f9ba5e353f7ced9

これは、gdb で octave を実行し、 でブレークポイントを設定することで確認されますpow。

同じことがmatlabにも当てはまる可能性がありますが、実際には確認できません。

于 2013-01-23T04:12:08.100 に答える
3

これは、Apple の のシステムpow()がこの場合に行うことに従う小さなテスト プログラムです。Source/Intel/xmm_power.cLibm-2026

#include <stdio.h>
int main() {
    // basically lines 1130-1157 of xmm_power.c, modified a bit to remove
    // irrelevant things

    double x = .3;
    int i = 3;

    //calculate ix = f**i
    long double ix = 1.0, lx = (long double) x;

    //calculate x**i by doing lots of multiplication
    int mask = 1;

    //for each of the bits set in i, multiply ix by x**(2**bit_position)
    while(i != 0)
    {
        if( i & mask )
        {
            ix *= lx;
            i -= mask;
        }
        mask += mask;
        lx *= lx; // In double this might overflow spuriously, but not in long double
    }

    printf("%.40f\n", (double) ix);
}

これは を出力します。これは、Matlab とPython で0.0269999999999999962252417162744677625597取得した結果と一致します(後者は単にこのコードを呼び出していることがわかっています)。対照的に、私にとっては を取得します。これは、小数点以下の桁数まで出力するように要求した場合と同じものであり、おそらく最も近い double です。.3 ^ 3.3 ** 3.3 * .3 * .30.02699999999999999969468866822808195138350.027

そこでアルゴリズムです。各ステップで設定されている値を正確に追跡することはできますが、それを行うための別のアルゴリズムを考えると、非常にわずかに小さい数値に丸められることはそれほど驚くべきことではありません。

于 2013-01-23T02:05:16.800 に答える
3

@Dougalのおかげで、私はこれを見つけました：

#include <stdio.h>
int main() {
  double x = 0.3;
  printf("%.40f\n", (x*x*x));

  long double y = 0.3;
  printf("%.40f\n", (double)(y*y*y));
}

与える：

0.0269999999999999996946886682280819513835
0.0269999999999999962252417162744677625597

より多くの桁で計算すると最悪の結果が得られるため、ケースは奇妙です。これは、とにかく最初の数値 0.3 が数桁で概算されるため、比較的「大きな」誤差で開始されるためです。この特定のケースでは、数桁の計算で別の「大きな」エラーが発生しますが、符号が反対になります...したがって、最初のエラーを補償します。代わりに、より多くの桁で計算すると、2 番目に小さいエラーが発生しますが、最初のエラーは残ります。

于 2013-01-23T05:59:42.970 に答える
-7

Goldberg の「What Every Computer Scientist Should Know About Floating-Point Arithmetic」を読んでください(これは Oracle による転載です)。それを理解してください。浮動小数点数は微積分の実数ではありません。申し訳ありませんが、TL;DR バージョンはありません。

于 2013-01-23T00:36:10.583 に答える