## `psi` -- the digamma/polygamma function

#### Introduction

`psi``(x)` represents the digamma function, i.e., the logarithmic derivative `diff(ln(gamma(x)),x)` of the `gamma` function.

`psi``(x,n)` represents the n-th polygamma function, i.e., the n-th derivative `diff(psi(x),x$n)`.

#### Call(s)

`psi(x)`
`psi(x, n)`

#### Parameters

`x` | - | an arithmetical expression
`n` | - | a nonnegative integer

#### Returns

an arithmetical expression.

### Overloadable:

`x`

#### Side Effects

When called with a floating point value `x`, the function is sensitive to the environment variable `DIGITS` which determines the numerical working precision.

#### Related Functions

`fact`, `gamma`, `zeta`

#### Details

- `psi(x, 0)` is equivalent to `psi(x)`.
- The digamma/polygamma function is defined for all complex arguments apart from the singular points 0,-1,-2,.. .
- If `x` is a floating point value, then a floating point value is returned.

  If `x` is a positive integer smaller than 1000, or if `x` is an odd integer multiple of 1/2 of absolute value smaller than 1000, then the relation

             psi(x+1,n)=psi(x,n)+(-1)^n*n!/x^(n+1)

  is applied. In conjunction with

  `psi(1) = -``EULER`,

  `psi(1,n) = (-1)^(n+1)*n!*``zeta``(n+1)`, `n>0`,

  `psi(1/2) = -2*ln(2) -` `EULER`,

  `psi(1/2,n) = (-1)^(n+1)*n!*(2^(n+1)-1)*``zeta``(n+1)`, `n>0`

  an explicit expression for the value of `psi` is returned.

  The special values `psi(infinity) = psi(infinity,0) = infinity` and `psi(infinity,n) = 0` for `n>0` are implemented.

  For all other arguments, a symbolic function call of `psi` is returned.

- The float attribute of the digamma function `psi(x)` is a kernel function, i.e., floating point evaluation is fast. The float attribute of the polygamma function `psi(x,n)` with `n>0` is a library function. Note that `psi``(float(x))` and `psi``(float(x),n)` rather than `float(psi(x))` and `float(psi(x,n))` should be used for float evaluation, because for integers and odd integer multiples of 1/2, the computation of the symbolic result `psi``(x)`, `psi``(x,n)` is costly and its float evaluation may be numerically unstable.
- The `expand` attribute uses

             psi(x+1,n)=psi(x,n)+(-1)^n*n!/x^(n+1)

  to rewrite `psi(x)`. For numerical `x`, this formula is used to shift the argument to the range `0<x<1`. Cf. examples 3 and 4.

#### Example 1

We demonstrate some calls with exact and symbolic input data:

    >> psi(-3/2), psi(4, 1), psi(3/2, 2)

                                         2
                                       PI
                8/3 - 2 ln(2) - EULER, --- - 49/36, 16 - 14 zeta(3)
                                        6

    >> psi(x + sqrt(2), 4), psi(infinity, 5)

                                         1/2
                                psi(x + 2   , 4), 0

Floating point values are computed for floating point arguments:

    >> psi(-5.2), psi(1.0, 3), psi(2.0 + 3.0*I, 10)

              6.065773152, 6.493939402, 0.7526409593 - 2.299472238 I

#### Example 2

`psi` is singular for nonpositive integers:

    >> psi(-2)

          Error: singularity [psi]

#### Example 3

For positive integers and odd integer multiples of 1/2, the result is expressed in terms of `EULER`, `PI`, `ln`, and `zeta`, respectively, if the absolute value of the argument is smaller than 1000:

    >> psi(-5/2), psi(-3/2, 1), psi(4, 3), psi(9/2, 2)

                                     2           4
                                   PI          PI
          46/15 - 2 ln(2) - EULER, --- + 40/9, --- - 1393/216,
                                    2          15

             19410176/1157625 - 14 zeta(3)

For larger arguments, the `expand` attribute can be used to obtain such expressions:

    >> psi(1000, 1)

                                   psi(1000, 1)

    >> expand(%)

            2
          PI
          --- -
           6

          835458876624295851523752364295.../50820720104325812617835292...

#### Example 4

The functions `diff`, `expand`, `float`, `limit`, and `series` handle expressions involving `psi`:

    >> diff(psi(x^2 + 1, 3), x), float(ln(3 + psi(sqrt(PI))))

                                   2
                          2 x psi(x  + 1, 4), 1.183103343

    >> expand(psi(15/4)), expand(psi(x + 3, 2))

                                             2       2          2
             psi(3/4) + 524/231, psi(x, 2) + -- + -------- + --------
                                              3          3          3
                                             x    (x + 1)    (x + 2)

    >> limit(x*psi(x), x = 0), limit(psi(x, 3), x = infinity)

                                       -1, 0

    >> series(psi(x), x = 0), series(psi(x, 3), x = infinity, 3)

                            2                 3   4
            1           x PI     2           x  PI       4
          - - - EULER + ----- - x  zeta(3) + ------ + O(x ),
            x             6                    90

             2    3    2     / 1  \
             -- + -- + -- + O| -- |
              3    4    5    |  6 |
             x    x    x     \ x  /

#### Changes

- Explicit expressions are now returned for all positive integers and odd integer multiples of 1/2 of size |x|<1000. The `expand` attribute now also rewrites symbolic calls with numerical arguments. The special values `psi(infinity)` = `psi(infinity,0)` =0 and `psi(infinity,n)` =0 for `n>0` were implemented.

Do you have questions or comments?
Copyright © SciFace Software GmbH & Co. KG 2000