【问题标题】:Incomplete gamma function does not accept complex value inputs不完整的 gamma 函数不接受复数值输入
【发布时间】:2018-07-31 08:38:11
【问题描述】:

我正在尝试计算upper incomplete gamma function 中定义的this post。如果我使用

from scipy.special import gamma,gammainc
from numpy import linspace

a = 0
z = (2+3j)*np.linspace(0,10)
gamma(a)*(1-gammainc(a,z))

z 是一个复数向量我得到一个错误

TypeError: ufunc 'gammainc' not supported for the input types, and the inputs could not be safely coerced to any supported types according to the casting rule ''safe''

是否有替代函数来进行计算?当我尝试使用 WolframAlpha 的 Gamma 函数执行此操作时,似乎没有问题。

【问题讨论】:

  • 不完全伽马函数是通过积分定义的,x(或z)是积分边界。你希望这个积分如何在复平面上进行?复平面中的积分通常取决于积分路径而不仅仅是边界。我不知道 WolframAlpha 是如何处理这种情况的,但这将是首先要弄清楚的事情。另请注意,对于您示例中的特定情况 (a = 0),积分有一个解析解,您可以使用 scipy.special.expi 函数来计算结果(对于实际值边界)。
  • @a_guest 至少对于具有正实部的参数,它是明确的:那里没有被积函数的极点,而因子 exp(-t) 在无穷大处提供了很好的衰减。在左半平面上,解析延拓使用沿负实轴切割的分支。

标签: python scipy gamma-function


【解决方案1】:

当 SciPy 不足以处理棘手的特殊功能时,mpmath 通常会提供帮助。比较

>>> mpmath.gammainc(0, 2+3j)
mpc(real='-0.024826207944199364', imag='0.020316674911044622')

同样来自Wolfram Alpha

使用 Python 编写,比 SciPy 慢;它也没有矢量化。所以用你的数据它会像

import mpmath
result = np.array([mpmath.gammainc(w) for w in z[1:]], dtype=np.complex)

请注意,我避免将 0 作为参数传递(它是一个极点)。 mpmath.gammainc 的返回类型是它自己的 mpc 对象类型,但可以如上转换回 NumPy。

【讨论】:

  • 谢谢,我还找到了一个 sympy 函数 sympy.functions.special.gamma_functions.uppergamma 来解决我的问题。在速度和准确性方面,您一般会推荐mpmathsympy 函数吗?我将在我的真实程序中输入一个包含 10^6 到 10^7 点的数组(不是帖子中的 MWE)。
  • SymPy 函数调用 mpmath 进行数值计算,它自己并不执行它们。没有理由通过 SymPy 访问 mpmath 函数,它增加了函数调用的开销。像这样的大型阵列需要时间。如果有帮助,请考虑降低 mpmath 的精度,mpmath.org/doc/1.0.0/basics.html#setting-the-precision
猜你喜欢
  • 2017-03-13
  • 1970-01-01
  • 2021-08-07
  • 2018-09-19
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 2013-05-12
  • 2016-03-10
相关资源
最近更新 更多