Cannot compute incomplete gamma function for some arguments

[araujoms]

gamma(Double64(z),x) works fine when x < 3, but if I call for example gamma(Double64(z),3) it chokes. The reason is that it goes on to a different branch of the code that uses logabsgamma(), which is not defined for Double64 arguments. I cannot do the extension like you did using Quadmath as it does not implement logabsgamma(), so I think the simplest solution would be defining

SpecialFunctions.logabsgamma(x::Double64) = Double64.(logabsgamma(BigFloat(x)))

A bit unrelated, but I'm confused by what you are doing in the file gamma_erf.jl. Why do you write

loggamma(x::Double64) = Double64Float128(lgamma, x)

instead of

SpecialFunctions.loggamma(x::Double64) = Double64Float128(lgamma, x)

In this way you're not extending the function, but shadowing it, so it conflicts with the version from SpecialFunctions.

[araujoms]

I can fix this myself but I need some guidance.

[JeffreySarnoff]

Glad you are into this! The shadowing is unintended. I will change that.
SpecialFunctions.logabsgamma(x::Double64) = Double64.(logabsgamma(BigFloat(x))) is simplest.
I would like to keep the current implmentation for x<3 and use your suggestion for x>=3.

[araujoms]

Great! Currently there's no implementation of logabsgamma, so adding it will only change the behaviour of gamma for x >= 3, where it currently fails. I'll do a PR later today then.

[JeffreySarnoff]

Thank you. When using BigFloats inside DoubleFloats, do this:

fn(x::DoubleFloat{T}) where {T<:IEEEFloat}
  keep_precision = precision(BigFloat)
  setprecision(BigFloat, 128)
  # do stuff, get result as a BigFloat
  setprecision(BigFloat, keep_precision)
  DoubleFloat{T}(result)
end