NaN in digamma caused by _cotpi

[Youjack]

The definition of _digamma(z::ComplexOrReal{Float64}) contains the formula

SpecialFunctions.jl/src/gamma.jl

Line 29 in cf35c58

ψ = -π * _cotpi(z)

in which _cotpi is defined as

SpecialFunctions.jl/src/gamma.jl

Line 63 in cf35c58

_cotpi(x) = @static VERSION >= v"1.10.0-DEV.525" ? inv(tanpi(x)) : cospi(x) / sinpi(x)

For large complex arguments, _cotpi produces NaN since cospi and sinpi produce Inf. For example,

julia> cospi(500 * exp(im * 3pi/4))
Inf - Inf*im

julia> sinpi(500 * exp(im * 3pi/4))
Inf + Inf*im

julia> cospi(500 * exp(im * 3pi/4)) / sinpi(500 * exp(im * 3pi/4))
NaN + NaN*im

Nevertheless, the ordinary cot gives the correct answer

julia> cot(pi * 500 * exp(im * 3pi/4))
0.0 - 1.0im

A way to avoid NaN is to define cotpi as in mpmath
https://github.com/mpmath/mpmath/blob/b600dbcabf4b7406a61e82b9e607f754a9f12ff9/mpmath/libfp.py#L153-L161

[stevengj]

Wouldn't it be better to fix tanpi in base?

[Youjack]

I have added a comment under a tanpi issue in the Julia repo.