Python vs. MATLAB 计算无穷大的积分，结果不同，替代方案（即将 Gauss-Legendre 正交扩展到 -x-> Infinity）？

Question

对于从 (-x 或 0) -> 无穷大的积分，我在 MATLAB 的 quadgk 和 Python 的 quad 例程之间得到不一致的结果。我相信 MATLAB 版本是正确的（基于将 flag 参数从 1 切换到 -1 的意义检查），而 Python 版本给出错误的结果，在本例中为 0。MATLAB 产生 0.1022。 integrands 是相同的，我已经列出了每一步，甚至将 MATLAB 的 quadgk 生成的 x 值插入到 Python 中（这导致 Python 版本生成与 MATLAB 相同的值，只是将它们传递给 integrand 函数）。在这一点上，我希望使用另一个例程而不是 SciPy，例如此处的 Gauss-Legendre 正交 https://sourceforge.net/projects/fastgausslegendrequadrature/，但我不确定如何将其从 a/b 范围扩展到 -a ->infinity（我见过这些方法只能达到有限的数量：

whereas b=np.Inf results in NaN. Also not sure how to setup the integration from the returned nodes and weights, although I've been reading the transformation but only for a and b finite ranges: https://pomax.github.io/bezierinfo/legendre-gauss.html 或者如果有人知道可以处理这个的 Python 库 - 我真的不喜欢 quad 不是矢量化的事实，并且可能会在 Cython 中编写代码，因为我必须快速集成 600,000 个函数（即 link 到上面的 C++ 库 link）。这里真正奇怪的是，我通过将 vol 输入上移任何位置 >= 0.39 来获得相同的结果，低于 Python 的结果在 0 处崩溃。非常混乱。感谢任何帮助，微积分已经有好几年了……这里是 Python 代码：

from scipy.stats import norm, lognorm
from scipy.integrate import quad
import numpy as np

def integrand(x, flag, F, K, vol, T2, T1):
    d1 = (np.log(x / (x+K)) + 0.5 * (vol**2) * (T2-T1)) / (vol * np.sqrt(T2 - T1))
    d2 = d1 - vol*np.sqrt(T2 - T1)
    mu = np.log(F) - 0.5 *vol **2 * T1
    sigma = vol * np.sqrt(T1)
    value = lognorm.pdf(x, scale=np.exp(mu), s=sigma) * (flag * x*norm.cdf(flag * d1) - flag * (x+K)*norm.cdf(flag * d2))
    return value

if __name__ == '__main__':
        flag = 1
        F = 54.31
        K = 1.1967
        vol = 0.1328
        T2 = 0.0411
        T1 = 0.0137
        quad(integrand, 0, np.Inf, args=(flag, F, K, vol, T2, T1), epsabs=1e-12)[0]

这是 MATLAB 代码（必须将 integrand 保存为 .M 然后可以在命令 window 中输入脚本）：

function value = integrand(x, flag, F,K,vol,T2,T1)
d1 = (log(x ./ (x+K)) + 0.5 .* (vol.^2) .* (T2-T1)) ./ (vol .* sqrt(T2 - T1));
d2 = d1 - vol.*sqrt(T2 - T1);
mu = log(F) - 0.5 .*vol .^2 .* T1;
sigma = vol .* sqrt(T1);
value = lognpdf(x, mu, sigma) .* (flag .* x.*normcdf(flag .* d1) - flag .* (x+K).*normcdf(flag .* d2));
end

% 脚本部分

flag = 1
F = 54.31
K = 1.1967
vol = 0.1328
T2 = 0.0411
T1 = 0.0137
quadgk(@(x) integrand(x,flag, F, K, vol, T2, T1), 0, Inf, 'AbsTol',1e-12)

我应该注意到 MATLAB 和 Python 在传递这些输入时使用 quad 生成相同的结果（转置以上变量）：

 current_opt = [  -1.0000    1.2075    0.1251    0.4300    0.0685    0.0411     
                1.0000    1.2075    0.0512    0.5600    0.0685    0.0411]

Answer 1

好的，这很有趣。除非 epsabs 变量设置得离谱地高，否则集成会分崩离析。我已经成功地使用 epsabs=-1e1000 在 MATLAB 和 Python 之间复制了结果。尽可能慢，但至少它有效。

Python vs. MATLAB 计算无穷大的积分，结果不同，替代方案（即将 Gauss-Legendre 正交扩展到 -x-> Infinity）？

Python vs. MATLAB computing an integral to infinity with different results, alternative (i.e. expand Gauss-Legendre quadrature to -x-> Infinity)?

python

matlab

cython

numerical-integration

quad