scipy 在数组上集成变量边界答案

【问题标题】：scipy integrate over array with variable boundsscipy 在数组上集成变量边界
【发布时间】：2013-12-05 15:56:23
【问题描述】：

我正在尝试在一个点列表上集成一个函数，并将整个数组传递给一个集成函数，以便对事物进行矢量化。对于初学者来说，调用 scipy.integrate.quad 太慢了，因为我有 10 000 000 个积分需要积分。使用 scipy.integrate.romberg 会更快，几乎是瞬间完成，而 quad 很慢，因为您必须循环它或对其进行矢量化。

我的函数相当复杂，但出于演示目的，假设我想将 x^2 从 a 积分到 b，但 x 是一个标量数组来评估 x。例如

将 numpy 导入为 np

from scipy.integrate import quad, romberg

def integrand(x, y):   

    return x**2 + y**2


quad(integrand, 0, 10, args=(10) # this fails since y is not a scalar

romberg(integrand, 0, 10)  # y works here, giving the integral over
                           # the entire range

但这仅适用于固定界限。有没有办法做类似的事情

z = np.arange(20,30)
romberg(integrand, 0, z)  # Fails since the function doesn't seem to
                          # support variable bounds

我看到的唯一方法是在 numpy 中重新实现算法本身并使用它来代替，这样我就可以拥有可变界限。任何支持这样的功能？还有 romb，你必须直接提供被积函数的值和一个 dx 区间，但这对于我的复杂函数来说太不精确了（marcum Q 函数，找不到任何实现，这可能是另一种方式来点它)。

【问题讨论】：

标签： python numpy scipy integrate

【解决方案1】：

尝试评估特殊函数的最佳方法是编写一个函数，该函数使用函数的属性在所有参数范围内快速准确地评估它。单一方法不太可能为所有参数范围提供准确（甚至稳定）的结果。在这种情况下，直接计算积分几乎肯定会在许多情况下失败。

话虽如此，在许多范围内评估积分的一般问题可以通过将积分转换为微分方程并求解来解决。粗略地说，步骤是

假设积分 I(t) 是函数 f(x) 从 0 到 t 的积分 [这可以推广到任意下限]，将其写为微分方程 dI/dt = f(x)。
在从 0 到 t 的某个时间范围内，使用 scipy.integrate.odeint() 在某些初始条件（此处为 I(0)）求解此微分方程。此范围应包含所有感兴趣的限制。采样的精细程度取决于函数以及需要评估的准确度。
结果将是我们输入的 t 集合从 0 到 t 的积分值。我们可以使用插值将其转换为“连续”函数。例如，使用样条曲线我们可以定义i = scipy.interpolate.InterpolatedUnivariateSpline(t,I)。
分别给定数组b 和a 中的一组上限和下限，然后我们可以将它们一次性评估为res=i(b)-i(a)。

这种方法是否适用于您的情况将需要您在参数范围内仔细研究它。另请注意，Marcum Q 函数涉及半无限积分。原则上这不是问题，只需将积分转换为有限范围内的积分即可。例如，考虑转换 x->1/x。不能保证这种方法对于您的问题在数值上是稳定的。

【讨论】：