如何用c++编写numpy.tensordot

Question

我正试图在c++中复制numpy.tensordot。numpy文档中的示例显示了一个嵌套循环，我可以开始工作，但如果不是呢

c = np.tensordot(a,b, axes=([1,0],[0,1]))

我想做：

c = np.tensordot(a,b, axes=([1,2],[0,1]))

在python中，这个新的嵌套循环会是什么样子？在c++中，有没有一种更简单/更快的方法来进行此操作？现在我正在使用相同的嵌套"；对于"；用c++中的std:：vector循环。我看到了一些可能会有所帮助的库，但我正在尝试只使用c++标准库。

这是一个numpy示例，以及指向文档的链接：https://numpy.org/doc/stable/reference/generated/numpy.tensordot.html

Examples
A “traditional” example:
>>>
a = np.arange(60.).reshape(3,4,5)
b = np.arange(24.).reshape(4,3,2)
c = np.tensordot(a,b, axes=([1,0],[0,1]))
c.shape
(5, 2)
c
array([[4400., 4730.],
[4532., 4874.],
[4664., 5018.],
[4796., 5162.],
[4928., 5306.]])
# A slower but equivalent way of computing the same...
d = np.zeros((5,2))
for i in range(5):
for j in range(2):
for k in range(3):
for n in range(4):
d[i,j] += a[k,n,i] * b[n,k,j]
c == d
array([[ True,  True],
[ True,  True],
[ True,  True],
[ True,  True],
[ True,  True]])

谢谢

Answer 1

我发现首先重写到np.einsum是有帮助的，因为循环代码的结果在概念上看起来非常相似：

a = np.random.rand(16, 8, 2)
b = np.random.rand(8, 2, 1)

c =  np.tensordot(a, b, axes=([1,2],[0,1]))
# same thing written with einsum
c_ein = np.einsum("ijk,jko->io", a, b)
# same thing done with for loops, 
# notice how we can use the same letters and indexing as einsum
c_manual = np.zeros((16, 1))
for i in range(16):
for o in range(1):
# j and k are summed since they don't appear in output
total = 0
for j in range(8):
for k in range(2):
total += a[i, j, k] * b[j, k, o]
c_manual[i, o] = total
assert np.allclose(c, c_ein, c_manual)