为什么我的 Numpy 测试代码比 Matlab 慢 2 倍

Question

我一直在用 Python 开发基于菲涅耳系数的反射率求解器，但遇到了一些障碍，因为 Python + Numpy 的性能比 Matlab 慢 2 倍。我已将问题代码提炼成一个简单的示例，以显示在每种情况下执行的操作：

测试用例的 Python 代码：

import numpy as np
import time

def compare_fn(i):
    a = np.random.rand(400)
    vec = np.random.rand(400)
    t = time.time()
    for j in xrange(i):
        a = (2.3 + a * np.exp(2j*vec))/(1 + (2.3 * a * np.exp(2j*vec)))
    print (time.time()-t)
    return a

a = compare_fn(200000)

输出：10.7989997864

等效的 Matlab 代码：

function a = compare_fn(i)

    a = rand(1, 400);
    vec = rand(1, 400);
    tic
    for m = 1:i
        a = (2.3 + a .* exp(2j*vec))./(1 + (2.3 * a .* exp(2j*vec)));
    end
    toc

a = compare_fn(200000); 经过的时间是 5.644673 秒。

我被这个难住了。我已经安装了 MKL（Anaconda 学术许可证）。如果有任何帮助，我将非常感谢帮助我确定我的示例中的问题，以及如何使用 Numpy 获得同等的（如果不是更好的）性能。

一般来说，我不能并行化循环，因为求解多层的菲涅耳系数涉及递归计算，可以用上述循环的形式表示。

Answer 1

以下内容类似于 unutbu 的已删除答案，对于您的示例输入，在我的系统上运行速度提高了 3 倍。如果您在 Matlab 中像这样实现它，它可能也会运行得更快，但那是另一回事。为了能够使用 ipython 的 %timeit 功能，我将您的原始函数改写为：

def fn(a, vec, i):
    for j in xrange(i):
        a = (2.3 + a * np.exp(2j*vec))/(1 + (2.3 * a * np.exp(2j*vec)))
    return a

我已经通过从循环中删除指数计算来优化它：

def fn_bis(a, vec, n):
    exp_vec = np.exp(2j*vec)
    for j in xrange(n):
        a = (2.3 + a * exp_vec) / (1 + 2.3 * a * exp_vec)
    return a

采用两种方法进行试驾：

In [2]: a = np.random.rand(400)

In [3]: vec = np.random.rand(400)

In [9]: np.allclose(fn(a, vec, 100), fn_bis(a, vec, 100))
Out[9]: True

In [10]: %timeit fn(a, vec, 100)
100 loops, best of 3: 8.43 ms per loop

In [11]: %timeit fn_bis(a, vec, 100)
100 loops, best of 3: 2.57 ms per loop

In [12]: %timeit fn(a, vec, 200000)
1 loops, best of 3: 16.9 s per loop

In [13]: %timeit fn_bis(a, vec, 200000)
1 loops, best of 3: 5.25 s per loop

Answer 2

对于我最初的问题中的示例，我一直在做很多实验来尝试确定 Matlab 和 Python/Numpy 之间速度差异的来源。一些主要发现是：

1. Matlab 现在有一个 JIT 编译器，它在涉及循环的情况下提供了显着的好处。关闭它会使性能降低 2 倍，使其速度与原生 Python + Numpy 代码相似。

   功能加速关闭

   a = compare_fn(200000);

   经过的时间是 9.098062 秒。

2. 然后我开始探索使用 Numba 和 Cython 优化示例函数的选项，看看我能做得更好。对我来说，一个重要的发现是显式循环计算上的 Numba JIT 优化比 Numpy 数组上的原生矢量化数学运算更快。我不太明白为什么会这样，但我在下面包含了我的示例代码和测试时间。我也玩过 Cython（我不是专家），虽然它也更快，但 Numba 仍然比 Cython 快 2 倍，所以我最终还是坚持使用 Numba 进行测试。

这是 3 个等效函数的代码。第一个是 Numba 优化函数，具有显式循环来执行元素计算。第二个函数是一个 Python+Numpy 函数，它依赖于 Numpy 向量化来执行计算。第三个函数尝试使用 Numba 来优化向量化的 Numpy 代码（如您在结果中所见，未能改进）。最后，我已经包含了 Cython 代码，尽管我只测试了一个案例。

import numpy as np
import numba as nb

@nb.jit(nb.complex128[:](nb.int16, nb.int16))
def compare_fn_jit(i, j):
    a = np.asarray(np.random.rand(j), dtype=np.complex128)
    vec = np.random.rand(j)
    exp_term = np.exp(2j*vec)

    for k in xrange(i):
        for l in xrange(j):
            a[l] = (2.3 + a[l] * exp_term[l])/(1 + (2.3 * a[l] * exp_term[l]))
    return a

def compare_fn(i, j):
    a = np.asarray(np.random.rand(j), dtype=np.complex128)
    vec = np.random.rand(j)
    exp_term = np.exp(2j*vec)
    for k in xrange(i):
        a = (2.3 + a * exp_term)/(1 + (2.3 * a * exp_term))
    return a

compare_fn_jit2 = nb.jit(nb.complex128[:](nb.int16, nb.int16))(compare_fn)


import numpy as np
cimport numpy as np
cimport cython
@cython.boundscheck(False)
def compare_fn_cython(int i, int j):
    cdef int k, l
    cdef np.ndarray[np.complex128_t, ndim=1] a, vec, exp_term
    a = np.asarray(np.random.rand(j), dtype=np.complex128)
    vec = np.asarray(np.random.rand(j), dtype=np.complex128)
    exp_term = np.exp(2j*vec)

    for k in xrange(i):
        for l in xrange(j):
            a[l] = (2.3 + a[l] * exp_term[l])/(1 + (2.3 * a[l] * exp_term[l]))
    return a

计时结果：

我。单个外循环的计时 - 展示矢量化计算的效率

%timeit -n 1 -r 10 compare_fn_jit(1,1000000) 1 个循环，10 次中最好的：352 每个循环的毫秒数

%timeit -n 1 -r 10 compare_fn(1,1000000) 1 次循环，最好的 10 次：498 毫秒 每个循环

%timeit -n 1 -r 10 compare_fn_jit2(1,1000000) 1 个循环，最好的 10：497 每个循环的毫秒数

%timeit -n 1 -r 10 compare_fn_cython(1,1000000) 1 次循环，最好的 10： 每个循环 424 毫秒

二。对短数组进行计算的大循环极端情况下的计时（预计 Numpy+Python 性能不佳）

%timeit -n 1 -r 5 compare_fn_jit(1000000,40) 1 个循环，最好的 5 个：1.44 每个循环的秒数

%timeit -n 1 -r 5 compare_fn(1000000,40) 1 次循环，最好的 5 次：28.2 秒 每个循环

%timeit -n 1 -r 5 compare_fn_jit2(1000000,40) 1 个循环，最好的 5: 29 s 每个循环

三。测试上述两种情况的中间位置

%timeit -n 1 -r 5 compare_fn_jit(100000,400) 1 次循环，最好的 5 次：1.4 秒 每个循环

%timeit -n 1 -r 5 compare_fn(100000,400) 1 次循环，最好的 5 次：5.26 秒 每个循环

%timeit -n 1 -r 5 compare_fn_jit2(100000,400) 1 个循环，最好的 5 个：5.34 每个循环的秒数

如您所见，对于这种特殊情况，使用 Numba 可以将效率提高 1.5 倍到 30 倍不等。与 Cython 相比，它的效率以及使用和实施的容易程度给我留下了深刻的印象。

Answer 3

我不知道 numpypy 是否足以满足你正在做的事情，但你可以尝试一下。

http://buildbot.pypy.org/numpy-status/latest.html