为什么 numpy.linalg.norm 在处理小数据时多次调用时会很慢？

Question

import numpy as np
from datetime import datetime
import math

def norm(l):
    s = 0
    for i in l:
        s += i**2
    return math.sqrt(s)

def foo(a, b, f):
    l = range(a)
    s = datetime.now()
    for i in range(b):
        f(l)
    e = datetime.now()
    return e-s

foo(10**4, 10**5, norm)
foo(10**4, 10**5, np.linalg.norm)
foo(10**2, 10**7, norm)
foo(10**2, 10**7, np.linalg.norm)

我得到以下输出：

0:00:43.156278
0:00:23.923239
0:00:44.184835
0:01:00.343875

似乎当针对小尺寸数据多次调用

np.linalg.norm

时，它的运行速度比我的标准函数慢。

这是什么原因？

Answer 1

首先：

datetime.now()

不适合衡量性能，它包括挂机时间，当高优先级进程运行或Python GC启动时，您可能会选择一个不好的时间（对于您的计算机），..

Python 中有专用的计时函数/模块：IPython/Jupyter 中的内置

timeit

模块或

%timeit

以及其他几个外部模块（如

perf

，...）

让我们看看如果我在您的数据上使用这些会发生什么：

import numpy as np
import math

def norm(l):
    s = 0
    for i in l:
        s += i**2
    return math.sqrt(s)

r1 = range(10**4)
r2 = range(10**2)

%timeit norm(r1)
3.34 ms ± 150 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)
%timeit np.linalg.norm(r1)
1.05 ms ± 3.92 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

%timeit norm(r2)
30.8 µs ± 1.53 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)
%timeit np.linalg.norm(r2)
14.2 µs ± 313 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)

对于短迭代来说，它并不慢，它仍然更快。但请注意，如果您已经拥有 NumPy 数组，那么 NumPy 函数的真正优势就会显现出来：

a1 = np.arange(10**4)
a2 = np.arange(10**2)

%timeit np.linalg.norm(a1)
18.7 µs ± 539 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)
%timeit np.linalg.norm(a2)
4.03 µs ± 157 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)

是的，现在速度快多了。 18.7us 与 1ms - 对于 10000 个元素几乎快了 100 倍。这意味着示例中

np.linalg.norm

的大部分时间都花在将

range

转换为

np.array

上。

Answer 2

你走在正确的道路上

np.linalg.norm

在小数组上有相当高的开销。在大型数组上，jit 编译函数和

np.linalg.norm

都在内存瓶颈中运行，这对于大多数时间执行简单乘法的函数来说是预期的。

如果从另一个 jitted 函数调用 jitted 函数，它可能会被内联，这可能比 numpy-norm 函数有更大的优势。

示例

import numba as nb
import numpy as np

@nb.njit(fastmath=True)
def norm(l):
    s = 0.
    for i in range(l.shape[0]):
        s += l[i]**2
    return np.sqrt(s)

性能

r1 = np.array(np.arange(10**2),dtype=np.int32)
Numba:0.42µs
linalg:4.46µs

r1 = np.array(np.arange(10**2),dtype=np.int32)
Numba:8.9µs
linalg:13.4µs

r1 = np.array(np.arange(10**2),dtype=np.float64)
Numba:0.35µs
linalg:3.71µs

r2 = np.array(np.arange(10**4), dtype=np.float64)
Numba:1.4µs
linalg:5.6µs

测量性能

• 在测量之前调用一次 jit-compiled 函数（第一次调用时有静态编译开销）
• 明确测量是否有效（由于小型数组保留在处理器缓存中，因此在实际示例中可能会出现超过 RAM 吞吐量的乐观结果。示例）

Answer 3

如果平方求和时的速度和潜在溢出是一个问题，但使用像

numba

这样的低级方法不是一个选择，那么

scipy.linalg.norm

可以用作优化低级例程的高级接口。 它依赖于 BLAS 例程

nrm2

（参见源代码），该例程使用缩放来避免溢出。

使用

perfplot

，我发现这个函数（禁用有限检查）比

numpy.linalg.norm

更快，当大小超过 10,000 时，这显然会改变计算/内存处理。对于大小高达 10,000 的 “small” 向量，

numpy

似乎主要受一些恒定开销的支配，并且仅在大小为 100,000 及以上时才达到与 BLAS 例程相同的水平。
根据应用程序的不同，100,000 不一定再小，但已经是一个相当大的向量了。

重现代码

### Imports ### 
import numpy as np
from scipy.linalg import norm as scipy_norm

import perfplot

### Performance comparison ###
perfplot.show(
    setup=lambda n: np.random.rand(n),
    kernels=[
        lambda a: np.linalg.norm(a),
        lambda a: scipy_norm(a, check_finite=False),
    ],
    labels=[
        "numpy",
        "scipy",
    ],
    n_range=[2**k for k in range(20)],
    xlabel="len(a)",
)