mrpy.base.special.polygamma

mrpy.base.special.polygamma(m, z)

Polygamma function.

Note

This is exactly as defined by mpmath, but modified to take array_like arguments, and return results with type float.

Notes

Gives the polygamma function of order m of z, psi^{(m)}(z). Special cases are known as the digamma function (psi^{(0)}(z)), the trigamma function (psi^{(1)}(z)), etc. The polygamma functions are defined as the logarithmic derivatives of the gamma function:

\[\psi^{(m)}(z) = \left(\frac{d}{dz}\right)^{m+1} \log \Gamma(z)\]

In particular, psi^{(0)}(z) = Gamma’(z)/Gamma(z). In the present implementation of psi(), the order m must be a nonnegative integer, while the argument z may be an arbitrary complex number (with exception for the polygamma function’s poles at z = 0, -1, -2, ldots).

Examples

For various rational arguments, the polygamma function reduces to a combination of standard mathematical constants:

>>> from mpmath import *
>>> mp.dps = 25; mp.pretty = True
>>> psi(0, 1), -euler
(-0.5772156649015328606065121, -0.5772156649015328606065121)
>>> psi(1, '1/4'), pi**2+8*catalan
(17.19732915450711073927132, 17.19732915450711073927132)
>>> psi(2, '1/2'), -14*apery
(-16.82879664423431999559633, -16.82879664423431999559633)

The polygamma functions are derivatives of each other:

>>> diff(lambda x: psi(3, x), pi), psi(4, pi)
(-0.1105749312578862734526952, -0.1105749312578862734526952)
>>> quad(lambda x: psi(4, x), [2, 3]), psi(3,3)-psi(3,2)
(-0.375, -0.375)

The digamma function diverges logarithmically as z to infty, while higher orders tend to zero:

>>> psi(0,inf), psi(1,inf), psi(2,inf)
(+inf, 0.0, 0.0)

Evaluation for a complex argument:

>>> psi(2, -1-2j)
(0.03902435405364952654838445 + 0.1574325240413029954685366j)

Evaluation is supported for large orders m and/or large arguments z:

>>> psi(3, 10**100)
2.0e-300
>>> psi(250, 10**30+10**20*j)
(-1.293142504363642687204865e-7010 + 3.232856260909107391513108e-7018j)

Application to infinite series

Any infinite series where the summand is a rational function of the index k can be evaluated in closed form in terms of polygamma functions of the roots and poles of the summand:

>>> a = sqrt(2)
>>> b = sqrt(3)
>>> nsum(lambda k: 1/((k+a)**2*(k+b)), [0, inf])
0.4049668927517857061917531
>>> (psi(0,a)-psi(0,b)-a*psi(1,a)+b*psi(1,a))/(a-b)**2
0.4049668927517857061917531

This follows from the series representation (m > 0)

\[\psi^{(m)}(z) = (-1)^{m+1} m! \sum_{k=0}^{\infty} \frac{1}{(z+k)^{m+1}}.\]

Since the roots of a polynomial may be complex, it is sometimes necessary to use the complex polygamma function to evaluate an entirely real-valued sum:

>>> nsum(lambda k: 1/(k**2-2*k+3), [0, inf])
1.694361433907061256154665
>>> nprint(polyroots([1,-2,3]))
[(1.0 - 1.41421j), (1.0 + 1.41421j)]
>>> r1 = 1-sqrt(2)*j
>>> r2 = r1.conjugate()
>>> (psi(0,-r2)-psi(0,-r1))/(r1-r2)
(1.694361433907061256154665 + 0.0j)