如何使用 numpy 沿右轴做标量积并将过程矢量化
How to do a scalar product along the right axes with numpy and vectorize the process
我有维度为 (100, 100, 16, 16) 的 numpy 数组 'test'，它为 100x100 网格上的点提供了一个不同的 16x16 数组。
我还有一些特征值和向量，其中 vals 的维数为 (100, 100, 16)，vecs 的维数为 (100, 100, 16, 16)，其中 vecs[x, y, :, i] 是矩阵的第 i 个特征向量第i个特征值对应的点(x, y) vals[x, y, i].
现在我想在网格上的所有点上取数组的第一个特征向量，用测试矩阵做矩阵乘积，然后用数组的所有其他特征向量做结果向量的标量积网格上的所有点并将它们相加。
生成的数组应具有维度 (100, 100)。在此之后，我想取数组的第二个特征向量，将它与 test 矩阵相乘，然后将结果的标量乘以所有不是第二个的特征向量，依此类推，所以最后我有 16 (100 , 100) 或者更确切地说是 (100, 100, 16) 数组。到目前为止，我只成功了很多我想避免的 for 循环，但是使用 tensordot 给了我错误的维度，我不知道如何选择为 np.dot 函数矢量化的轴。
我听说 einsum 可能适合这项任务，但不依赖于 python 循环的一切对我来说都很好。
import numpy as np
from numpy import linalg as la
test = np.arange(16*16*100*100).reshape((100, 100, 16, 16))
vals, vecs = la.eig(test + 1)
np.tensordot(vecs, test, axes=[2, 3]).shape
>>> (10, 10, 16, 10, 10, 16)
编辑：好的，所以我使用 np.einsum 得到了正确的中间结果。
np.einsum('ijkl, ijkm -> ijlm', vecs, test)
但在下一步中，我只想对 vec 的所有其他条目执行标量积。我可以在这个 einsum 形式主义中实现一些反 Kronecker delta 吗？或者我现在应该切换回通常的 numpy 吗？
好的，我试了一下 np.einsum 我找到了一种方法来完成上述操作。 einsum 的一个很好的特性是，如果你在 'output' 中重复出现两次的索引（所以在 '->'-thing 的右边）你可以沿着一些轴进行元素乘法和沿着其他一些轴的收缩（一些东西你没有手写的张量代数符号）。
result = np.einsum('ijkl, ijlm -> ijkm', np.einsum('ijkl, ijkm -> ijlm', vecs, test), vecs)
这几乎可以解决问题。现在只需要去掉对角线项。我们可以通过像这样减去对角线项来做到这一点：
result = result - result * np.eye(np.shape(test)[-1])[None, None, ...]
我有维度为 (100, 100, 16, 16) 的 numpy 数组 'test'，它为 100x100 网格上的点提供了一个不同的 16x16 数组。 我还有一些特征值和向量，其中 vals 的维数为 (100, 100, 16)，vecs 的维数为 (100, 100, 16, 16)，其中 vecs[x, y, :, i] 是矩阵的第 i 个特征向量第i个特征值对应的点(x, y) vals[x, y, i].
现在我想在网格上的所有点上取数组的第一个特征向量，用测试矩阵做矩阵乘积，然后用数组的所有其他特征向量做结果向量的标量积网格上的所有点并将它们相加。 生成的数组应具有维度 (100, 100)。在此之后，我想取数组的第二个特征向量，将它与 test 矩阵相乘，然后将结果的标量乘以所有不是第二个的特征向量，依此类推，所以最后我有 16 (100 , 100) 或者更确切地说是 (100, 100, 16) 数组。到目前为止，我只成功了很多我想避免的 for 循环，但是使用 tensordot 给了我错误的维度，我不知道如何选择为 np.dot 函数矢量化的轴。 我听说 einsum 可能适合这项任务，但不依赖于 python 循环的一切对我来说都很好。
import numpy as np
from numpy import linalg as la
test = np.arange(16*16*100*100).reshape((100, 100, 16, 16))
vals, vecs = la.eig(test + 1)
np.tensordot(vecs, test, axes=[2, 3]).shape
>>> (10, 10, 16, 10, 10, 16)
编辑：好的，所以我使用 np.einsum 得到了正确的中间结果。
np.einsum('ijkl, ijkm -> ijlm', vecs, test)
但在下一步中，我只想对 vec 的所有其他条目执行标量积。我可以在这个 einsum 形式主义中实现一些反 Kronecker delta 吗？或者我现在应该切换回通常的 numpy 吗？
好的，我试了一下 np.einsum 我找到了一种方法来完成上述操作。 einsum 的一个很好的特性是，如果你在 'output' 中重复出现两次的索引（所以在 '->'-thing 的右边）你可以沿着一些轴进行元素乘法和沿着其他一些轴的收缩（一些东西你没有手写的张量代数符号）。
result = np.einsum('ijkl, ijlm -> ijkm', np.einsum('ijkl, ijkm -> ijlm', vecs, test), vecs)
这几乎可以解决问题。现在只需要去掉对角线项。我们可以通过像这样减去对角线项来做到这一点：
result = result - result * np.eye(np.shape(test)[-1])[None, None, ...]