具有重叠切片的 NumPy 就地操作

Question

考虑这个不合适的操作：

>>> b = numpy.asarray([1, 2, 3])
>>> b[1:] = b[1:] - b[:-1]
>>> b
array([1, 1, 1])

现在考虑就地操作：

>>> a = numpy.asarray([1, 2, 3])
>>> a[1:] -= a[:-1]
>>> a
array([1, 1, 2])

他们给出了不同的结果，这是我没想到的。

更新：这个行为似乎改变了；现在他们给出了相同的结果。

我本以为 NumPy 会以正确的顺序（向后）进行减法，这样它们就会给出与异地减法相同的结果。

我的问题是：这是 NumPy 的预期行为，还是错误，或者结果是否未定义？

Answer 1

未定义，或至少难以理解，可能是最佳答案。另一个 SO 回答声称

a[1:] -= a[:-1]

被口译员翻译成类似

的东西

a.__setitem__(slice(1,None), a.__getitem__(slice(1,None)).
   __isub__(a.__getitem__(slice(None,-1))))

In [171]: a=np.arange(10)
In [172]: a[1:] -= a[:-1]
In [173]: a
Out[173]: array([0, 1, 1, 2, 2, 3, 3, 4, 4, 5])

计算结果与：

相同

In [175]: for i in range(9):
     ...:     a[i+1] = a[i+1]-a[i]

我可以看出它是如何从嵌套的 __setitem__ 等表达式中得到它的。我尝试用 np.nditer.

复制它

你提到的反面是

In [178]: for i in range(8,-1,-1):
     ...:     a[i+1] = a[i+1]-a[i]

我不认为 numpy 可以推断出需要这样的反向迭代。 __setitem__ 的第二个参数在前向迭代中的计算结果很好。缓冲该术语是唯一简单的解决方案。

引入 .at ufunc 方法是为了解决 a[idx] += b 等表达式中的缓冲问题。特别是当 idx 有重复时。对 a 的影响应该是累积的还是应该只适用最后一个实例。

在您的示例中，if 的行为就像 a[1:] - a[:-1]:

In [165]: a=np.arange(10)
In [166]: idx=np.arange(1,10)
In [167]: np.add.at(a, idx, -a[:-1])
In [168]: a
Out[168]: array([0, 1, 1, 1, 1, 1, 1, 1, 1, 1])

那是因为 add.at 的第 3 个参数在使用前已被完全评估。这是一个临时副本。我从其他测试中知道 add.at 比普通的 a[idx] += 慢。 [我对 'buffering' add.at 绕过了什么感到有点困惑；这与这里出现问题的明显缺乏缓冲有何不同？]

为什么要使用 += 表示法？只是为了让代码更紧凑？或者希望让它更快？但是，如果速度是目标，我们是否希望 numpy 加入一些额外的缓冲，以使其更安全？

nditer 相当于 a[1:] -= a[:-1]

In [190]: a=np.arange(10)
In [191]: it = np.nditer([a[1:],a[1:],a[:-1]], op_flags=['readwrite'])
In [192]: for i,j,k in it:
     ...:     print(i,j,k)
     ...:     i[...] = j-k
     ...:     print(i)

1 1 0
1
2 2 1
1
3 3 1
2
4 4 2
2
5 5 2
3
6 6 3
...

可以简化迭代

In [197]: it = np.nditer([a[1:],a[:-1]], op_flags=['readwrite'])
In [198]: for i,j in it:
     ...:     i[...] -= j

因为它是一个视图，所以来自 a[:-1] 的迭代值将反映在前一个循环中所做的更改。

我不确定数组 += 表达式中是否使用了 c 版本的 nditer，但 nditer 的目的是巩固迭代编码成一个统一的框架。

另一个有趣的观察是如果我定义

idx = array([1, 2, 3, 4, 5, 6, 7, 8, 9])

然后

a[idx] -= a[:-1]
a[1:] -= a[idx-1]
a[idx] -= a[idx-1]

全部给出想要的array([0, 1, 1, 1, 1, 1, 1, 1, 1, 1])。换句话说，-= 的两边都必须是 'views/slices'。那一定是 add.at 绕过的缓冲。 a[idx-1] 是一个副本是显而易见的。 a[idx]-= 抛出一个缓冲区，而 a[1:]-= 没有，不是很明显。

Answer 2

此行为以前未定义，但是 since NumPy 1.13.0，具有重叠输入和输出的操作现在表现得好像先复制输入一样。引用发行说明：

Operations where ufunc input and output operands have memory overlap produced undefined results in previous NumPy versions, due to data dependency issues. In NumPy 1.13.0, results from such operations are now defined to be the same as for equivalent operations where there is no memory overlap.

Operations affected now make temporary copies, as needed to eliminate data dependency. As detecting these cases is computationally expensive, a heuristic is used, which may in rare cases result to needless temporary copies. For operations where the data dependency is simple enough for the heuristic to analyze, temporary copies will not be made even if the arrays overlap, if it can be deduced copies are not necessary. As an example,np.add(a, b, out=a) will not involve copies.

Here is the relevant issue on GitHub.

具有重叠切片的 NumPy 就地操作

NumPy in-place operations with overlapping slices

python

numpy

in-place

更新：这个行为似乎改变了；现在他们给出了相同的结果。