4

我正在使用 python 库 mpmath,特别是评估不完整的 gamma 函数。这是寻根例程的一部分,但对于复值参数的某些组合,它的评估非常慢。

import mpmath as mp
from mpmath import gammainc
def Gamma(a,z0,z1):
    return gammainc(a,a=z0,b=z1,regularized=False)

这里函数的评估mpmath.gammainc卡住了:

>> Gamma(mp.mpc(12.5+17.5j), mp.mpf(0.0), mp.mpf(-12.5))

另一方面,Mathematica 几乎立即返回结果:

In[1]:= Gamma[12.5 + 17.5 I, 0, -12.5]
Out[1]:= 2.38012*10^-7 + 5.54827*10^-7 I

在其他情况下,对于不同的参数mpmathMathematica返回相同的输出:

数学

In[2]: Gamma[3.5 I, 0, 10]
Out[2]:= 0.0054741 + 0.000409846 I

Python mpmath

>> Gamma(3.5j,0,10)
mpc(real='0.0054741038497352953', imag='0.00040984640431357779')

您对这种行为的原因有一些了解吗?这可以被认为mpmath是一个正交问题还是一个数学问题?不幸的是scipy,它没有提供gamma复杂参数的函数实现,所以它不是一个选项。

4

1 回答 1

3

gammainc显然,在某些情况下,mpmath 中的一个错误会导致它在评估时进入无限循环。值得报告mpmath 跟踪器。但至少对于您提到的情况,一种解决方法是将三参数不完全伽马函数编写为两个不完全伽马函数的差(参考)。上不完全伽马函数是通过传递两个参数来计算的gammainc(即,z1隐含地认为是正无穷大)。

def Gamma(a, z0, z1):
    return gammainc(a, z0) - gammainc(a, z1)

print(Gamma(12.5+17.5j, 0.0, -12.5)) 

打印(2.3801203496987e-7 + 5.54827238374855e-7j)与 WolframAlpha 一致。

于 2018-08-26T00:22:15.377 回答