2

我已经通读了einsum手册和ajcr的基本介绍

我在非编码环境中对爱因斯坦求和的经验为零，尽管我试图通过一些互联网研究来解决这个问题（会提供链接，但还没有超过两个的声誉）。我还尝试使用 einsum 在 python 中进行实验，看看是否可以更好地处理事情。

然而，我仍然不清楚这样做是否既可行又有效：

在长度 (3) 和高度 (n) 相等的数组 (a 和 b) 的两个数组上，逐行产生 (row i: a on b) 的外积加上 (row i: b on a )，然后将所有外部乘积矩阵相加以输出一个最终矩阵。

我知道 'i,j->ij' 会产生一个向量对另一个向量的外积——这是让我失去兴趣的下一步。（'ijk,jik->ij' 绝对不是）

我的另一个可用选项是遍历数组并从我用 cython 编写的函数调用基本函数（双外积和矩阵加法）（使用内置的 numpy 外部和 sum 函数不是一个选项，它太慢了）。我很可能最终也会将循环本身移动到 cython。

所以：

  1. 我怎样才能概括地表达我上面描述的过程？

  2. 它会比在 cython 中做所有事情提供真正的收益吗？还是有其他我不知道的选择？（包括我使用 numpy 的效率低于我可能的可能性......）

谢谢。

编辑示例：

A=np.zeros((3,3))
arrays_1=np.array([[1,0,0],[1,2,3],[0,1,0],[3,2,1]])
arrays_2=np.array([[1,2,3],[0,1,0],[1,0,0],[3,2,1]])
for i in range(len(arrays_1)):
  A=A+(np.outer(arrays_1[i], arrays_2[i])+np.outer(arrays_2[i],arrays_1[i]))

（但是请注意，实际上我们正在处理长度更大的数组（即每个内部成员的长度仍然为 3，但最多有几千个这样的成员），并且这部分代码（不可避免地）被调用了很多次）

如果它有帮助，这里是两个外部产品相加的cython：

def outer_product_sum(np.ndarray[DTYPE_t, ndim=1] a_in, np.ndarray[DTYPE_t, ndim=1] b_in):
    cdef double *a = <double *>a_in.data
    cdef double *b = <double *>b_in.data
    return np.array([
[a[0]*b[0]+a[0]*b[0], a[0]*b[1]+a[1]*b[0], a[0] * b[2]+a[2] * b[0]],
[a[1]*b[0]+a[0]*b[1], a[1]*b[1]+a[1]*b[1], a[1] * b[2]+a[2] * b[1]],
[a[2]*b[0]+a[0]*b[2], a[2]*b[1]+a[1]*b[2], a[2] * b[2]+a[2] * b[2]]])

现在，我从“i in range(len(array))”循环中调用，如上所示。

4

2 回答 2

3

爱因斯坦求和只能用于问题的乘法部分（即外积）。幸运的是，求和不必逐元素执行，但您可以在归约矩阵上执行此操作。使用示例中的数组：

arrays_1 = np.array([[1,0,0],[1,2,3],[0,1,0],[3,2,1]])
arrays_2 = np.array([[1,2,3],[0,1,0],[1,0,0],[3,2,1]])
A = np.einsum('ki,kj->ij', arrays_1, arrays_2) + np.einsum('ki,kj->ij', arrays_2, arrays_1)

输入数组的形状为 (4,3)，求和发生在第一个索引（名为'k'）上。如果要在第二个索引上进行求和，请将下标字符串更改为'ik,jk->ij'.

于 2016-02-22T10:53:45.450 回答
1

无论你能做什么np.einsum，你通常可以更快地使用np.dot. 在这种情况下，A是两个点积的和：

arrays_1 = np.array([[1,0,0],[1,2,3],[0,1,0],[3,2,1]])
arrays_2 = np.array([[1,2,3],[0,1,0],[1,0,0],[3,2,1]])

A1 = (np.einsum('ki,kj->ij', arrays_1, arrays_2) +
      np.einsum('ki,kj->ij', arrays_2, arrays_1))

A2 = arrays_1.T.dot(arrays_2) + arrays_2.T.dot(arrays_1)

print(np.allclose(A1, A2))
# True

%timeit (np.einsum('ki,kj->ij', arrays_1, arrays_2) +
         np.einsum('ki,kj->ij', arrays_2, arrays_1))
# 100000 loops, best of 3: 7.51 µs per loop

%timeit arrays_1.T.dot(arrays_2) + arrays_2.T.dot(arrays_1)
# 100000 loops, best of 3: 4.51 µs per loop
于 2016-02-22T12:01:12.907 回答