加快 scipy 中集成的函数评估

Question

我正在尝试将代码从 Matlab 移植到 SciPy。这是我到目前为止编写的代码的简化版本：https://gist.github.com/atmo/01b6e007be9ef90e402c。然而，Python 版本比 Matlab 慢得多。我在要点中包含了分析结果，它们表明几乎 90% 的时间 python 花在评估函数 f 上。除了用 C 或 Cython 重写之外，有什么方法可以加快其评估速度吗？

Answer 1

您的 numpy 版本在速度上可能与旧版 MATLAB 运行相当。但是新的 MATLAB 版本进行了各种形式的即时编译，大大加快了重复计算的速度。

我的猜测是，您可以蚕食 lambda 和 f 代码，并可能将它们的评估时间缩短一半。但真正的杀手是你打电话 f 这么多次。

首先，我会尝试预先计算 f 中的内容。例如定义 K1=K[1] 并在计算中使用 K1。这将减少索引调用的次数。指数是否重复？也许用常规 def 替换 lambda 定义，或者将其与 f.

组合

Answer 2

正如我在评论中提到的，如果您考虑到 ，您可以摆脱对 quad 的大约一半调用（以及因此复杂的函数 f）矩阵是对称的.

进一步的速度提升，仍然是纯粹的 python，将通过重写那个复杂的函数来获得。我在 sympy 中完成了大部分工作。

最后，我尝试使用 np.vectorize 向量化对 quad 的调用。

from scipy.integrate import quad
from scipy.special import jn as besselj
from scipy import exp, zeros, linspace
from scipy.linalg import norm
import numpy as np

def complicated_func(lmbd, a, n, k):
    u,v,w = 5, 3, 2
    x = a*lmbd
    fac = exp(2*x)
    comm = (2*w + x)
    part1 = ((v**2 + 4*w*(w + 2*x) + 2*x*(x - 1))*fac**5
        + 2*u*fac**4
        + (-v**2 - 4*(w*(3*w + 4*x + 1) + x*(x-2)) + 1)*fac**3
        + (-8*(w + x) + 2)*fac**2
        + (2*comm*(comm + 1) - 1)*fac)
    return part1/lmbd *besselj(n+1, lmbd) * besselj(k+1, lmbd)

def perform_quad(n, k, a):
    return quad(complicated_func, 0, np.inf, args=(a,n,k))[0]

def improved_main():
    sz = 20
    amatrix = np.zeros((sz,sz))
    ls = -np.linspace(1, 10, 20)/2
    inds = np.tril_indices(sz)
    myv3 = np.vectorize(perform_quad)
    res = myv3(inds[0], inds[1], ls.reshape(-1,1))
    results = np.empty(res.shape[0])
    for rowind, row in enumerate(res):
        amatrix[inds] = row
        symm_matrix = amatrix + amatrix.T - np.diag(amatrix.diagonal())
        results[rowind] = norm(symm_matrix)
    return results

计时结果显示我的速度提高了 5 倍（如果我只运行一次，你会原谅我，它需要足够长的时间）：

In [11]: %timeit -n1 -r1 improved_main()
1 loops, best of 1: 6.92 s per loop

In [12]: %timeit -n1 -r1 main()
1 loops, best of 1: 35.9 s per loop

如果您立即将 v 替换为它的平方，也会有微不足道的收获，因为这是它唯一一次用于那个复杂的函数：作为它的平方。

对 besselj 的调用中也存在大量重复，但我不知道如何避免这种情况，因为 quad 将确定 lmbd，因此您可以'轻松预先计算这些值，然后执行查找。

如果您分析 improved_main，您会发现对 complicated_func 的调用量几乎减少了 2 倍（仍然需要计算对角线）。所有其他速度提升都可以归因于 np.vectorize 和 complicated_func.

的改进

我的系统上没有Matlab，所以如果你在那里改进复杂的功能，我不能对它的速度增益做任何陈述。

加快 scipy 中集成的函数评估

Speed up function evaluation for integration in scipy

python

numpy

scipy