for most functions if f(::Int64)::Float64 then f(::BigInt)::BigFloat

[oxinabox]

As I learnt when doing #17447 (comment)
the expectation is that, for a function f (eg gamma, exp2)

if f(::Int64)::Float64 then f(::BigInt)::BigFloat

And for the most part this is true (indeed I added in the PR test cases to ensure it).

I have found 3 exceptions to it, all functions from a part of math the I've not really dealt with much. zeta, digamma, and trigamma

digamma

julia> digamma(big"1")
-0.5772156649015315

julia> digamma(big"1.0")
-5.772156649015328606065120900824024310421593359399235988057672348848677267776685e-01

This is the simple case.
The fallback definition for BigInt, is falling down to digamma(z::Number) at special/gamma.jl:368
and so never ends up calling digamma(x::BigFloat) at mpfr.jl:450, which is the natural call back

zeta

julia> zeta(0)
-0.5

julia> zeta(big"0")
-0.5

julia> zeta(big"0.0")
-5.000000000000000000000000000000000000000000000000000000000000000000000000000000e-01

julia> zeta(big"1.0")
inf #this has type BigFloat

julia> zeta(big"1")
NaN #this has type Float64

julia> julia> zeta(1.0)
NaN

julia> zeta(1)
NaN

julia> zeta(complex(big"1",big"0"))
nan + 0.000000000000000000000000000000000000000000000000000000000000000000000000000000im

So it is working fine for complex numbers, which is what most people care about (From the little I know of the Riemann zeta function).

For BigInt it is hitting zeta(x::Integer) at special/gamma.jl:462
Which is similar to the issue for digamma

however we also see that zeta(big"1.0") != zeta(big("1")) (ignoring type).
NaN !=inf. which is more concerning.
Divining what the correct answer for zeta(1) should be, I leave to better informed minds than my own.
WolframAlpha says it is "Complex Infinity"

Trigamma

julia> trigamma(big"1")
1.6449340668482262

julia> trigamma(big"1.0")
1.6449340668482262

trigamma is not defined for BigFloat, or BigInt, so perhaps it doesn't belong here, but just pointing it out since i spotted it.

[oxinabox]

I made a script to search for these cases: https://gist.github.com/oxinabox/36e939544c2d93c1baa4d34fca776622

Ran it on 0.4.6 (my 0.5 is not accessibe right now, but I think the script should work on 0.5)
I'm not honestly sure that all of these are worth fixing.
But given I have now found them, I figured I would share my notes.

Returns a Float64

• cond -- I'm not entirely sure what this does (https://en.wikipedia.org/wiki/Condition_number). But for both BigFloat, BigInt, Int it always returns a 1.0 as a Float64
• digamma -- see above
• zeta -- see above
• trigamma -- see above
• eta -- hits eta(::Integer), and thus never calls eta(::BigFloat)
• invdigamma -- like trigamma, both BigInt and BigFloat inputs result in Float64 outputs. cos it just hits invdigamma.

Throw errors

And the following functions which throw errors on BigInts, but not on Int64s of the same value
Which from my understanding of them is not reasonable.
I haven't investigated why, but these are all mathematical functions that say they have a method that works on BigInt (or at least its super type),
And its not like they throw an error for values in the region I checked, since they works fine with a Int64

The error in all 6 of these cases is ERROR: MethodError: **f** has no method matching **f**(::BigFloat).
So it tries to fallback to BigFloat and fails.

• dawson
• erfcinv
• erfcx
• erfi
• erfinv

Reasonable/required

• ctime the result being float64 is because of how the filesystem stores time
• mtime ditto
• float64 deprecated casting syntax, obviously should do cast
• peakflops this is specified to return peak flop rate in double precision, and the int argument is only for how many BLAS threads

[simonbyrne]

Thanks for chasing this up. The main problem here is that we don't have BigFloat versions of most of those functions. However we should probably be consistent as to whether we give a Float64 or throw an error.

[simonbyrne]

It seems we do have methods for digamma, zeta and eta, so they could be changed.

We should probably change the behaviour of trigamma and invdigamma to throw an error on BigFloat arguments.

cond is a bit funny, as it is really intended for matrices, so it's treating the number as a 1x1 matrix. In that case it isn't really correct either (e.g. for Float32).

[oxinabox]

We should probably change the behaviour of trigamma and invdigamma to throw an error on BigFloat arguments

Ideally it would not be throwing errors, but instead not saying it has a method defined on BigInt in the first place.
The ones I mentioned above all through MethodErrors, but through them as errors saying the method was not defined for BigFloat.
It doesn't make too much a difference, but it does indicate mis-specification of the accepted types in the type signatures.

[stevengj]

If foo(x::Integer) = foo(float(x)), then it will automatically convert BigInt to BigFloat and will hence automatically throw a MethodError if there is no BigFloat method of foo.

Many of the f64(x) calls in gamma.jl seem problematic to me. Instead of foo(x::Number) = foo(f64(f)) methods, we should probably have something like:

foo{T<:Union{Integer,Rational}}(x::Union{T,Complex{T}}) = foo(float(x))
foo{T<:Union{Float16,Float32}}(x::Union{T,Complex{T}}) = oftype(x, foo(f64(x)))

[simonbyrne]

Instead of

foo{T<:Union{Integer,Rational}}(x::Union{T,Complex{T}}) = foo(float(x))

what's wrong with

foo(x::Number) = foo(float(x))

(which also works for Irrationals as well)

[simonbyrne]

Oh right, the mess in #14979...

We would also need

foo(x::FloatingPoint) = throw(MethodError(foo,(x,)))

[stevengj]

Probably to be absolutely sure you don't get circular method calls you would need

function foo(x::Number)
   y = float(x)
   typeof(y) === typeof(x) && throw(MethodError(foo, (x,)))
   return foo(y)
end

(Think about what float does for Complex{T}, Quaternion{T}, etc.)

[oxinabox]

Manually throwing MethodError seems generally problematic to me.

It is a bit weird, and potentially confusing if when doing some protyping in the REPL,
and I need a method foo, for a BigInt, so I check if it exists:

julia>  methods(foo, (BigInt,))
1-element Array{Any,1}:
 foo(x::number) at math.jl:115

then when I run

julia> foo(big"10000000000000000000000000000000000000")
ERROR: MethodError: `foo` has no method matching foo(::BigInt)

Why would the methods lie to me like that?
It said there was a method, but I got a MethodError.
And the docstring for MethodError says: "A method with the required type signature does not exist in the given generic function."

The reverse is also true, when ever I get a method error, my first step is to open the REPL and call methods -- often to find that I got the parameters the wrong way round.
If on calling methods I find a compatible method declared, then I will be left confused, not knowing what to do next. All I can do at that point to workout what is going on is dig through the julia source code.

In short, any time MethodError is thrown manually, it is making the methods function a lier, and thus hurting the interactive/REPL help system.

[stevengj]

@oxinabox, it is a common pattern in Julia to dispatch from a generic signature like foo(::Number) to more concrete signatures, which may or may not exist. Because of that, method_exists is not a reliable way to tell whether a method works for a given type. This is life with multiple dispatch.

[oxinabox]

I closed this in code that was eventually merged in:
JuliaMath/SpecialFunctions.jl#18