添加零时奇怪的numpy.sum行为

shx*_*hx2 22 python numpy sum numerical-stability

我理解,由于数值误差(例如,以不同的顺序求和浮点数),数学上等效的算术运算会如何导致不同的结果.

但是,添加零sum可以改变结果让我感到惊讶.无论如何,我认为这总是适用于花车:x + 0. == x.

这是一个例子.我预计所有的线都是零.任何人都可以解释为什么会这样吗?

M = 4  # number of random values
Z = 4  # number of additional zeros
for i in range(20):
    a = np.random.rand(M)
    b = np.zeros(M+Z)
    b[:M] = a
    print a.sum() - b.sum()

-4.4408920985e-16
0.0
0.0
0.0
4.4408920985e-16
0.0
-4.4408920985e-16
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
2.22044604925e-16
0.0
4.4408920985e-16
4.4408920985e-16
0.0
Run Code Online (Sandbox Code Playgroud)

对于较小的M和,似乎不会发生Z.

我也确定了a.dtype==b.dtype.

这是另一个例子,它还演示了python的内置sum行为符合预期:

a = np.array([0.1,      1.0/3,      1.0/7,      1.0/13, 1.0/23])
b = np.array([0.1, 0.0, 1.0/3, 0.0, 1.0/7, 0.0, 1.0/13, 1.0/23])
print a.sum() - b.sum()
=> -1.11022302463e-16
print sum(a) - sum(b)
=> 0.0
Run Code Online (Sandbox Code Playgroud)

我正在使用numpy V1.9.2.

jor*_*b87 8

简短的回答:你看到的区别

a + b + c + d
Run Code Online (Sandbox Code Playgroud)

(a + b) + (c + d)
Run Code Online (Sandbox Code Playgroud)

因浮点不准确而不一样.

答案很长: Numpy实现成对求和作为速度优化(它允许更容易的矢量化)和舍入误差.

numpy sum-implementation可以在这里找到(函数pairwise_sum_@TYPE@).它基本上做了以下事情:

  1. 如果阵列的长度小于8,则执行常规的循环求和.这就是为什么W < 4在你的情况下没有观察到奇怪的结果- 在这两种情况下将使用相同的for-loop求和.
  2. 如果长度在8到128之间,它会将总和累计在8个区间中,r[0]-r[7]然后将它们相加((r[0] + r[1]) + (r[2] + r[3])) + ((r[4] + r[5]) + (r[6] + r[7])).
  3. 否则,它递归地将数组的两半相加.

因此,在第一种情况下,你得到a.sum() = a[0] + a[1] + a[2] + a[3]和在第二种情况下b.sum() = (a[0] + a[1]) + (a[2] + a[3])导致a.sum() - b.sum() != 0.