Numpy:多种外部产品

Question

一般问题

假设我有一个ndarray v形状(nrow,ncols,3).我想计算outer_array形状的ndarray ,(nrow,ncols,3,3)其中包含(3)每个索引处的形状向量的所有外积(nrow,ncol).当然,这就是numpy.einsum存在的问题.现在,我尝试过的是:

outer_array = numpy.einsum("xyi,xyj->xyij",v,v.conjugate())

现在,我不确定这是否有效:尽管outer_array具有预期形状的事实,外部产品矩阵的元素与我期望的不一致.

我认为这是由于einsum表达式中标签的选择:产品应该被总结x并y因为索引被重复,但由于我在输出表达式中重用它们,因此总和的结果是以某种方式广播的.

另一方面,如果我写:

outer_array = numpy.einsum("xyi,uvj->...ij",v,v.conjugate())

numpy将为每对计算外部产品的所有可能组合,(x,y)并(u,v)产生一个形状阵列(ncols,nrow,ncols,nrow,3,3),其中对角线 (u,v) = (x,y)将包含所需的输出.

精确的问题

如何在einsum表示法中选择前两个索引以获得一个数组,在每个索引处x,y我得到vector的外积v而不必求助于第二个表达式？

显然编辑,这个表单似乎也有效:

outer_array = numpy.einsum("...i,...j->...ij",v,v.conjugate())

我只能欣赏numpy广播的有用性!

Answer 1

'xyi,xyj->xyij'工作的关键xy是在输出字符串中重复.

让我们使用一个更简单的数组:

x = np.arange(6).reshape(3,2)
np.einsum.einsum('ij->j',x)
# array([6, 9])  # sums on the 1st dimension of `x`

现在对于这个外部产品x:

In [20]: x[:,:,None]*x[:,None,:]  # shape (3,2,2)
Out[20]: 
array([[[ 0,  0],
        [ 0,  1]],

       [[ 4,  6],
        [ 6,  9]],

       [[16, 20],
        [20, 25]]])

这是numpy广播的一个例子(即添加维度并扩展它们)

In "...i,...j->...ij",...更多地作为现有但匿名维度的占位符.

相当于einsum:

np.einsum('ij,ik->ijk',x,x)

(我应该在最后2个维度中进行非对称的计算).

我已经找到了一个纯粹的Python工作einsum.重点是解析参数字符串,以及它如何为iter对象创建输入.它可以在github上找到:https://github.com/hpaulj/numpy-einsum/blob/master/einsum_py.py 欢迎您使用它.它有一个debug标志来显示中间步骤.

用我的einsum调试输出:

In [23]: einsum_py.myeinsum('ij,ik->ijk',x,x,debug=True)
# ... some parsing information
[('ij', [105, 106], 'NONE'), ('ik', [105, 107], 'NONE')]
('ijk', [105, 106, 107], 'NONE')
iter labels: [105, 106, 107],'ijk'

op_axes [[0, 1, -1], [0, -1, 1], [0, 1, 2]]

op_axes是用于创建a的关键参数iter,该对象迭代输入和输出数组的轴.它遍历所有数组的第1轴.对于第一个op和输出,第二个轴是1,对于第二个op,第二个轴是newaxis(-1).

随着ellipsis:

In [24]: einsum_py.myeinsum('...j,...k->...jk',x,x,debug=True)
...
iter labels: [0, 106, 107],'0jk'
op_axes [[0, 1, -1], [0, -1, 1], [0, 1, 2]]

这会生成op_axes相同的计算结果.