python - 在 Cython 中更快地计算平方范数

Question

我想计算平方范数，可以写成

W 是一个矩阵，有 V 行。你是一个向量。我有两个 numpy 对象 W 和 u。

import numpy as np
import numpy.random as npr
V = 10
W = npr.normal(size=(V, 3))
u = npr.normal(size=(1,3))

如果我逐行计算，我可以这样做:

res = np.zeros(V)
for v in range(V):
    res[v] = (W[v] - u).dot((W[v] - u).transpose())

但是一旦 V 变大(比如 5000)，它可能会很慢，我需要一次又一次地重新计算它。所以我尝试了矩阵乘法，但是没有成功，因为它不是逐行乘法。

((W - u).transpose()).dot(W - u)

如何在 Numpy 中快速计算平方范数？

我打算使用 Cython，那么循环中的逐行乘法是否可以比 Numpy 更快？我了解到 Cython 有 parallelization ，但如果我在内部使用 Numpy 对象，Cython 似乎无法并行化 for 循环(添加于 11 月 20 日:可能我不能在 prange 内部使用 Python 对象，但我可以使用 Numpy 对象)。

Answer 1

方法#1

您可以利用基于快速 BLAS 的 np.dot直接使用 NumPy，无需任何循环，就像这样 -

res = (W**2).sum(1) + (u**2).sum(1) -2*W.dot(u.ravel())

引入np.einsum为 u -

获取 W 和 np.inner 的逐行求和

res = np.einsum('ij,ij->i',W,W) + np.inner(u,u).ravel() -2*W.dot(u.ravel())

方法说明

(W[v] - u).dot((W[v] - u) 在每次迭代中，我们都在做内点积。因为，我们对所有行都这样做W，利用广播转换为 ((W - u)**2).sum(1)。

现在，

(Xik-Yjk)**2 = Xik**2 + Yjk**2 - 2*Xik*Yjk

因此，

sum_k((Xik-Yjk)**2) = sum_k(Xik**2) + sum_k(Yjk**2) - 2*sum_k(Xik*Yjk)

RHS 上的最后一项基本上是矩阵乘法，我们正在利用 np.dot。

方法#2

或者，使用更多的np.einsum，它会像这样-

d = W -u
res = np.einsum('ij,ij->i',d,d)