Странное поведение numpy.sum при добавлении нулей

Я понимаю, как математически эквивалентные арифметические операции могут привести к разным результатам из-за числовых ошибок (например, суммирование чисел с плавающей запятой в разных порядках).

Однако меня удивляет, что добавление нулей к 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

Кажется, этого не происходит при меньших значениях 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

Я использую numpy V1.9.2.


person shx2    schedule 23.06.2015    source источник
comment
Я могу воспроизвести с 1.9.2, но не с 1.6.1. Я предполагаю, что более длинный массив каким-то образом заставляет элементы добавляться в другом порядке, например, для облегчения SMID.   -  person amaurea    schedule 23.06.2015
comment
С math.fsum(), который отличается от sum(), похоже, этого не происходит (по крайней мере, в 2.7.x).   -  person Armin Rigo    schedule 23.06.2015
comment
Рискну предположить, что это попарное суммирование, github.com/numpy/numpy/pull/3685   -  person ev-br    schedule 23.06.2015
comment
Невозможно воспроизвести с python 2.7.6 и numpy 1.8.0   -  person Noel Segura Meraz    schedule 24.06.2015
comment
Также не воспроизводится на python 2.7.5 и numpy 1.6.2.   -  person skippy    schedule 26.06.2015


Ответы (1)


Короткий ответ: вы видите разницу между

a + b + c + d

а также

(a + b) + (c + d)

что из-за неточностей с плавающей запятой не то же самое.

Длинный ответ: Numpy реализует попарное суммирование как оптимизацию как скорости (это упрощает векторизацию), так и ошибки округления.

Реализацию суммы numpy можно найти здесь ( функция pairwise_sum_@TYPE@). По сути, он делает следующее:

  1. Если длина массива меньше 8, выполняется обычное суммирование цикла for. Вот почему странный результат не наблюдается, если W < 4 в вашем случае - в обоих случаях будет использоваться одно и то же суммирование цикла for.
  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.

person jornb87    schedule 26.06.2015