0

numpyやや大きな値で外積を使用すると、奇妙な動作と思われる現象が発生しています。

たとえば、次は正しいようです。

r = 1e15
a = array([1, 2, 3]) * r;
b = array([-1, 2, 1]) * r;

c = cross(a / norm(a), b / norm(b));

print(dot(c, a))  # outputs 0.0

しかし、指数を 1 上げると、次のようになります。

r = 1e16
a = array([1, 2, 3]) * r;
b = array([-1, 2, 1]) * r;

c = cross(a / norm(a), b / norm(b));

print(dot(c, a))  # outputs 2.0

指数の値が大きいほど、数値はさらに奇妙になります。ここで何が起こっているか知っている人はいますか？ありがとう！

4

1 に答える 1

2

丸め誤差が発生しています。デフォルトでarray()は、 を含むオブジェクトを返しますdtype=float64。rどんどん大きくしていくと、配列積を正確に表すための仮数スペースが不足しています。これをテストする方法は次のとおりです。

def testcross(r, dt):
    a = array([1, 2, 3], dtype=dt)*r
    b = array([-1, 2, 1], dtype=dt)*r
    c = cross(a/norm(a), b/norm(b))
    return dot(c, a)

for rr in logspace(4, 15, 10):
    print "%10.2f %10.2f %g" % (testcross(rr, float32), testcross(rr, float64)

結果:

     -0.00       0.00 10000
      0.00      -0.00 166810
      0.00       0.00 2.78256e+06
     -4.00       0.00 4.64159e+07
    -64.00       0.00 7.74264e+08
   1024.00       0.00 1.29155e+10
      0.00       0.00 2.15443e+11
-524288.00       0.00 3.59381e+12
      0.00      -0.02 5.99484e+13
-134217728.00       0.00 1e+15

でさえ「完璧」ではないことに注意しfloat64てr=5.99484e13ください。これは、 に到達するずっと前に精度が低下し始めていることを示していr=1e15ますfloat64。予想通り、精度が低いと事態はさらに悪化しfloat32ます。

OPの提案に従ってください。32ビットと64ビットの浮動小数点表現の仮数フィールドは、それぞれ24ビットと53ビットです（暗黙のビットを含む）。を取るlog10([2**24, 2**53])と、これはそれぞれ約 7 桁と 16 桁に相当することがわかります。これは、最初に指摘されたように、およびの周りr=4.6e7に表示される表形式のエラーに対応しています。丸めは、内積によって基礎となる行列計算が大きな数を減算する場合に発生し、その差はいずれかの大きな数の仮数で表すことができません。float32r=1e16

于 2016-06-03T01:13:27.127 に答える