大系数数组的Python集成

Question

我正在尝试优化python中的简单集成

from scipy import integrate
import numpy as np
from scipy.special import kv
import time
#Example function 
def integrand(x, a, b, c):
    return a * (x ** (-b)) * (np.sqrt(x ** (c) + 1) - 1)
#Real Function that I want to calculate
def Bes(xx):
    return integrate.quad(lambda x: kv(5./3.,x), xx,np.inf)
def F(x,a,b,c,d,e,f):
    zx = 1/((x**2.+1)*a)
    feq = e*x**(f)
    if (x>c):
       feq *= c/x * np.exp(-(x/d)**2.)
    return b*Bes(zx)*feq*x**2.
start = time.time()
array_length = 10
a = np.random.rand(array_length)+3.
b = np.random.rand(array_length)+1.
c = np.random.rand(array_length)
d = (np.random.rand(array_length)+1)*100.
e = np.random.rand(array_length)*100.
f = np.random.rand(array_length)
inte = np.array([])
for i in range(array_length):
    result = integrate.quad(lambda x: F(x, a[i], b[i], c[i],d[i],e[i],f[i]),0.01,100000.)
    inte = np.append(inte,result[0])
print("For array length = %i" % array_length)
print("Time = %.2f [sec]" %(time.time()-start))

但我面临的问题是

• a、 b，c是长度大于10^7（相同长度）的数组
• x的积分范围从0.01开始，扩展到无穷大
• 小x处的积分（如[0.01，1]）非常重要，需要很小的步骤。你知道吗

我想对每个系数值积分这个函数，并高效地返回整个积分数组作为结果（长度~10^7）。你知道吗

我应该使用什么样的工具？你知道吗

（+）我刚刚将代码从简单的示例更改为需要解决的实际集成形式。抱歉搞混了。你知道吗

Answer 1

我怀疑这个积分会收敛于b和c的某些值，所以我试着用Sympy来计算它：

import sympy
sympy.init_printing()

a, b, c = sympy.symbols('a, b, c', positive=True)

x = sympy.Symbol('x', positive=True)

sympy.integrate(a*(x**(-b))*(sympy.sqrt(x**c+1)-1), (x, 0, sympy.oo))

这意味着只要您的系数通过了check函数，您就应该能够使用此代码获得正确的结果。你知道吗

from numpy import sqrt, pi

from scipy.special import gamma

def check(a, b, c):
    assert (-(-b + 1)/c < 1)
    assert (1/2 - (-b + 1)/c > 1)
    assert (1 - (-b + 1)/c > 1)

def result(a, b, c):
    return a*gamma(-b/c + 1 + 1/c)*gamma(b/c - 1/2 - 1/c)/(2*sqrt(pi)*(b - 1))