Internals of BigFloat getting changed in 1.12

[KristofferC]

See JuliaLang/julia#55906.

Now causes

Error During Test at /home/pkgeval/.julia/packages/SpecialFunctions/Zijv9/test/gamma_inc.jl:231
  Test threw exception
  Expression: gamma(a, x) ≈ gamma(big(a), big(x))
  not compatible with mpfr
  Stacktrace:
    [1] error(s::String)
      @ Base ./error.jl:44
    [2] unsafe_convert(::Type{Ref{BigFloat}}, x::Base.RefValue{BigFloat})
      @ Base.MPFR ./mpfr.jl:212
    [3] _gamma_big(a::BigInt, x::BigFloat)
      @ SpecialFunctions ~/.julia/packages/SpecialFunctions/Zijv9/src/gamma_inc.jl:1159
    [4] _gamma
      @ ~/.julia/packages/SpecialFunctions/Zijv9/src/gamma_inc.jl:1147 [inlined]
    [5] gamma(a::BigInt, x::BigFloat)
      @ SpecialFunctions ~/.julia/packages/SpecialFunctions/Zijv9/src/gamma_inc.jl:1120
    [6] macro expansion
      @ /opt/julia/share/julia/stdlib/v1.12/Test/src/Test.jl:676 [inlined]
    [7] (::var"#14#15")()
      @ Main ~/.julia/packages/SpecialFunctions/Zijv9/test/gamma_inc.jl:231
    [8] setprecision(f::var"#14#15", ::Type{BigFloat}, prec::Int64; base::Int64)
      @ Base.MPFR ./mpfr.jl:1180
    [9] setprecision(f::Function, ::Type{BigFloat}, prec::Int64)
      @ Base.MPFR ./mpfr.jl:1179
   [10] top-level scope
      @ ~/.julia/packages/SpecialFunctions/Zijv9/test/gamma_inc.jl:229
   [11] macro expansion
      @ /opt/julia/share/julia/stdlib/v1.12/Test/src/Test.jl:1771 [inlined]
   [12] macro expansion
      @ ~/.julia/packages/SpecialFunctions/Zijv9/test/gamma_inc.jl:229 [inlined]
   [13] include(mapexpr::Function, mod::Module, _path::String)
      @ Base ./Base.jl:308
   [14] macro expansion
      @ ~/.julia/packages/SpecialFunctions/Zijv9/test/runtests.jl:45 [inlined]
   [15] macro expansion
      @ /opt/julia/share/julia/stdlib/v1.12/Test/src/Test.jl:1771 [inlined]

[inkydragon]

Simple MWE:

julia> using SpecialFunctions
Precompiling SpecialFunctions finished.
  1 dependency successfully precompiled in 2 seconds. 17 already precompiled.

julia> gamma(big(1.), 0)
1.0

julia> gamma(big(1), 0)
ERROR: not compatible with mpfr
Stacktrace:

[stevengj]

The current implmentation is:

using Base.MPFR: MPFRRoundingMode, ROUNDING_MODE

function _gamma(x::BigFloat)
    isnan(x) && return x
    z = BigFloat()
    ccall((:mpfr_gamma, :libmpfr), Int32, (Ref{BigFloat}, Ref{BigFloat}, Int32), z, x, ROUNDING_MODE[])
    isnan(z) && throw(DomainError(x, "NaN result for non-NaN input."))
    return z
end

The current master still passes things as Ref{BigFloat}. It seems like the rounding mode should be passed as:

 ccall((:mpfr_gamma, :libmpfr), Int32, (Ref{BigFloat}, Ref{BigFloat}, MPFRRoundingMode), z, x, Base.Rounding.rounding_raw(BigFloat))

but this shouldn't change anything.

We should also change :libmpfr to Base.MPFR.libmpfr,

What else should change?

[giordano]

I think you're looking at the wrong method? The failing one is

SpecialFunctions.jl/src/gamma_inc.jl

Lines 1149 to 1161 in f8fd782

	function _gamma_big(a::Real,x::Real)
	if x < 0
	# MPFR returns NaN in this case
	if isinteger(a) && a > 0
	return invoke(_gamma, Tuple{Number,Number}, a, BigFloat(x))
	else
	throw(DomainError((a, x), "gamma will only return a complex result if called with a complex argument"))
	end
	end
	z = BigFloat()
	ccall((:mpfr_gamma_inc, :libmpfr), Int32, (Ref{BigFloat}, Ref{BigFloat}, Ref{BigFloat}, Int32), z, a, x, ROUNDING_MODE[])
	return z
	end

for the incomplete gamma function: #485 (comment)

[giordano]

I believe the issue is that x::Int is passed in the ccall as x::Ref{BigFloat} when x >= 0. Julia v1.11:

julia> ccall((:mpfr_gamma_inc, :libmpfr), Int32, (Ref{BigFloat}, Ref{BigFloat}, Ref{BigFloat}, Int32), BigFloat(0.0), BigFloat(1.0), 1.0, Base.MPFR.ROUNDING_MODE[])
-1

julia nightly:

julia> ccall((:mpfr_gamma_inc, :libmpfr), Int32, (Ref{BigFloat}, Ref{BigFloat}, Ref{BigFloat}, Int32), BigFloat(0.0), BigFloat(1.0), 1.0, Base.MPFR.ROUNDING_MODE[])
ERROR: not compatible with mpfr
Stacktrace:
 [1] error(s::String)
   @ Base ./error.jl:44
 [2] unsafe_convert(::Type{Ref{BigFloat}}, x::Base.RefValue{BigFloat})
   @ Base.MPFR ./mpfr.jl:215
 [3] top-level scope
   @ ./REPL[6]:1

[stevengj]

Shouldn't Base.cconvert take care of that? Base.cconvert(Ref{BigFloat}, 1.0) works on 1.11

[giordano]

Based on the stacktrace we're hitting this unsafe_convert method, introduced by the PR linked above. The Base.cconvert is the argument to Bsae.unsafe_convert

[giordano]

Basically, 1.11:

julia> Base.unsafe_convert(Ref{BigFloat}, Base.cconvert(Ref{BigFloat}, 1.0))
Ptr{BigFloat} @0x00007fd8cd79c460

master:

julia> Base.unsafe_convert(Ref{BigFloat}, Base.cconvert(Ref{BigFloat}, 1.0))
ERROR: not compatible with mpfr
Stacktrace:
 [1] error(s::String)
   @ Base ./error.jl:44
 [2] unsafe_convert(::Type{Ref{BigFloat}}, x::Base.RefValue{BigFloat})
   @ Base.MPFR ./mpfr.jl:215
 [3] top-level scope
   @ REPL[7]:1

[stevengj]

Why can't cconvert and/or unsafe_convert be updated to support this?

[giordano]

Since this was specifically introduced in the PR which reworked BigFloat this looks completely deliberate.

[stevengj]

In particular, couldn't we just have:

Base.cconvert(::Type{Ref{BigFloat}}, x::Number) = convert(BigFloat, x).d

?

This mirrors the behavior in earlier versions, where cconvert creates a BigFloat object.

[giordano]

This mirrors the behavior in earlier versions, where cconvert creates a BigFloat object.

It still does.

[stevengj]

The problem is that cconvert(Ref{BigFloat}, x) calls convert(Ref{BigFloat}, x) (via the default fallback method), and that is no longer the correct conversion in the context of a ccall.

Doing this right is the whole job of cconvert, so this seems like a bug in cconvert on master: JuliaLang/julia#57367 (which makes _gamma_big work again).

[KristofferC]

Fixed in JuliaLang/julia#57367