用于 SciPy 集成和插值的 Numba答案

【问题标题】：Numba for SciPy Integration and Interpolation用于 SciPy 集成和插值的 Numba
【发布时间】：2021-09-30 03:09:39
【问题描述】：

我使用 Numba 来加速我的代码。它工作得很好，并提供了 2-3 倍的改进。然而，我的代码中花费的主要时间（大约 90%）是在 scipy quad 积分和插值（线性和三次样条）中。我做了几百次这些集成，所以我认为这是 Numba 可以提升的东西。看起来 Numba 不支持这些？我听说 Numba-Scipy 可以让 Numba 识别 Scipy，但这似乎仍然不起作用。有没有办法让 Numba 优化我的集成/插值？

【问题讨论】：

您能提供到目前为止您尝试过的内容以及遇到的错误吗？
scipy 是一组未连接的包。最后，我查看了numba-scipy，只“知道”了“特殊功能”模块。 quad 基本上是一个迭代函数，多次调用你的函数（至少 20 次）。您可以通过使您的功能尽可能快来提高其速度，但除此之外您无能为力。一些集成/优化函数已经编译了组件，但它们仍然受到调用您自己的函数所需时间的限制。没有人编译你的代码或将整个过程变成一个编译包。
请提供一个实际输入的工作示例（可能是随机数）。通常要走的路是创建一个 LowLevel callable，它提供了显着的加速。但是开箱即用不支持线性或三次样条类。如果只是对样条曲线的评估，这应该可以在 Numba 中实现。替代方案是 C 或 Cython 中的 LowLevel Callables，但存在相同的问题。例如。使用 np.interpolate 和低级可调用对象的示例。 stackoverflow.com/a/58561573/4045774

标签： python scipy integration interpolation numba

【解决方案1】：

刚刚为cquadpack 写了一个包装器，名为NumbaQuadpack，它应该可以满足您的要求：https://github.com/Nicholaswogan/NumbaQuadpack。 cquadpack 是 Quadpack 的 C 版本，是 scipy.integrate.quad 使用的。

from NumbaQuadpack import quadpack_sig, dqags
import numpy as np
import numba as nb
import timeit

@nb.cfunc(quadpack_sig)
def f(x, data):
    return x**2 + 2 + np.log(x)
funcptr = f.address
a = 0
b = 1
sol, abserr, success = dqags(funcptr, a, b)
print(sol) # definite integral solution

# test speed
@nb.njit()
def timetest_nb():
    sol, abserr, success = dqags(funcptr, a, b)
timetest_nb()
n_time=10000
print(timeit.Timer(timetest_nb).timeit(number=n_time)/n_time)

在我的电脑上，这个小积分需要 4.2 µs，而当我用 scipy.integrate.quad 做同样的事情时，需要 68.1 µs。

对于插值，只需使用np.interp（一维插值）。它可以在 numba-jitted 函数中使用。

一般来说，任何 C/C++ 或 fortran 代码都可以用 ctypes 包装，并从 numba-jitted 函数中调用。

【讨论】：

您是否还找到了一种通过 double func(double x, void *user_data) 将多个输入（numpy.ndarrays、floats、int64）传递给您的实现的方法？比这样做更优雅的东西stackoverflow.com/a/58561573/4045774
@max9111 不，我没有。我只是硬编码在一个 float64 数组中。如果 numba 有一个结构，你可以得到一个指针，然后从一个 jitted 函数中取消引用，那就太好了……这可能吗？
哦，我知道这是您在其他 SO 帖子中设法做到的。凉爽的。这有点笨拙，但我猜它有效。
@nicholaswogan，请问您可以为 scipy 的 interpolate.RectBivariateSpline 函数写类似的东西吗？我正在训练这样一个RectBivariateSpline 对象，这一步的速度还可以，但是我查询训练有素的插值器 10^7 次，需要 1 小时左右......所以查询训练有素的 2D 插值器是我的问题。
@velenos14 ，interpolate.RectBivariateSpline 的更快包装器是可能的。它只是Fitpack fortran 库的包装器。你可以做到的。遵循我制作的其中一个包的结构。例如github.com/Nicholaswogan/NumbaMinpack