Possible regression using quadgk function and Unitful limits

[mikeingold]

I think I might have discovered a recent regression in QuadGK when using the quadgk function with ranges specified in Unitful types. I tried to dig into the source by following the stack trace, and it seems like it might be related to the recent changes to handle_infinities, but I'm not familiar enough with this package or Julia's underlying type mechanics to nail down the issue at the moment.

Steps to Reproduce

I was able to reproduce the issue using the following minimal code snippet. This snippet should integrate a constant 1 Joule value over the domain [0,10] seconds, which should result in a value of 10 Joule-seconds.

using QuadGK, Unitful
f(t) = 1.0u"J"
quadgk(f, 0.0u"s", 10.0u"s")

Test Configuration 1 (v2.6.0, Error)

System information for this configuration:

versioninfo()
# Julia Version 1.8.5
# Commit 17cfb8e65e (2023-01-08 06:45 UTC)
# Platform Info:
#   OS: Windows (x86_64-w64-mingw32)
#   CPU: 32 × AMD Ryzen 9 3950X 16-Core Processor
#   WORD_SIZE: 64
#   LIBM: libopenlibm
#   LLVM: libLLVM-13.0.1 (ORCJIT, znver2)
#   Threads: 1 on 32 virtual cores

import Pkg
Pkg.status(["Unitful", "QuadGK"])
# Status `C:\Users\mikei\.julia\environments\v1.8\Project.toml`
#   [1fd47b50] QuadGK v2.6.0
#   [1986cc42] Unitful v1.12.2

Running the code snippet above produces an error

using QuadGK, Unitful

f(t) = 1.0u"J"
# f (generic function with 1 method)

quadgk(f, 0.0u"s", 10.0u"s")
#=
ERROR: MethodError: no method matching Quantity{Float64, 𝐓 , Unitful.FreeUnits{(s,), 𝐓 , nothing}}(::Int64)
The applicable method may be too new: running in world age 32472, while current world is 32475.
Closest candidates are:
  Quantity{T, D, U}(::Number) where {T, D, U} at C:\Users\mikei\.julia\packages\Unitful\fbiNW\src\types.jl:151 (method too new to be called from this world context.)
  (::Type{T})(::T) where T<:Number at boot.jl:772
  (::Type{T})(::AbstractChar) where T<:Union{AbstractChar, Number} at char.jl:50
  ...
Stacktrace:
  [1] convert(#unused#::Type{Quantity{Float64, 𝐓, Unitful.FreeUnits{(s,), 𝐓, nothing}}}  , x::Int64)
    @ Base .\number.jl:7
  [2] one(#unused#::Type{Quantity{Float64, 𝐓, Unitful.FreeUnits{(s,), 𝐓, nothing}}}  )
    @ Base .\number.jl:334
  [3] #s14#1
    @ C:\Users\mikei\.julia\packages\QuadGK\kf0xA\src\gausskronrod.jl:260 [inlined]
  [4] var"#s14#1"(T::Any, ::Any, #unused#::Any, n::Any)
    @ QuadGK .\none:0
  [5] (::Core.GeneratedFunctionStub)(::Any, ::Vararg{Any})
    @ Core .\boot.jl:582
  [6] do_quadgk(f::typeof(f), s::Tuple{Quantity{Float64, 𝐓, Unitful.FreeUnits{(s,), 𝐓, nothing}}, Quantity{Float64, 𝐓
, Unitful.FreeUnits{(s,), 𝐓, nothing}}} , n::Int64, atol::Nothing, rtol::Nothing, maxevals::Int64, nrm::typeof(LinearAlgebra.norm), segbuf::Nothing)
    @ QuadGK C:\Users\mikei\.julia\packages\QuadGK\kf0xA\src\adapt.jl:7
  [7] #28
    @ C:\Users\mikei\.julia\packages\QuadGK\kf0xA\src\adapt.jl:186 [inlined]
  [8] handle_infinities(workfunc::QuadGK.var"#28#29"{Nothing, Nothing, Int64, Int64, typeof(LinearAlgebra.norm), Nothing}, f::typeof(f), s::Tuple{Quantity{Float64, 𝐓, Unitful.FreeUnits{(s,), 𝐓, nothing}}, Quantity{Float64, 𝐓, Unitful.Fre
eUnits{(s,), 𝐓, nothing}}} )
    @ QuadGK C:\Users\mikei\.julia\packages\QuadGK\kf0xA\src\adapt.jl:115
  [9] quadgk(::Function, ::Quantity{Float64, 𝐓, Unitful.FreeUnits{(s,), 𝐓, nothing}}  , ::Vararg{Quantity{Float64, 𝐓, U
nitful.FreeUnits{(s,), 𝐓, nothing}}} ; atol::Nothing, rtol::Nothing, maxevals::Int64, order::Int64, norm::Function, segbuf::Nothing)
    @ QuadGK C:\Users\mikei\.julia\packages\QuadGK\kf0xA\src\adapt.jl:185
 [10] quadgk(::Function, ::Quantity{Float64, 𝐓, Unitful.FreeUnits{(s,), 𝐓, nothing}}  , ::Vararg{Quantity{Float64, 𝐓, U
nitful.FreeUnits{(s,), 𝐓, nothing}}} )
    @ QuadGK C:\Users\mikei\.julia\packages\QuadGK\kf0xA\src\adapt.jl:183
 [11] top-level scope
    @ REPL[3]:1
=#

Test Configuration 2 (v.2.5.0, Working)

For the next configuration, I simply pinned QuadGK to version 2.5.0 in Julia's Pkg REPL and re-ran the test. System information for this configuration:

versioninfo()
# Julia Version 1.8.5
# Commit 17cfb8e65e (2023-01-08 06:45 UTC)
# Platform Info:
#   OS: Windows (x86_64-w64-mingw32)
#   CPU: 32 × AMD Ryzen 9 3950X 16-Core Processor
#   WORD_SIZE: 64
#   LIBM: libopenlibm
#   LLVM: libLLVM-13.0.1 (ORCJIT, znver2)
#   Threads: 1 on 32 virtual cores

import Pkg
Pkg.status(["Unitful", "QuadGK"])
# ⌃ [1fd47b50] QuadGK v2.5.0 ⚲
#   [1986cc42] Unitful v1.12.2
# Info Packages marked with ⌃ have new versions available and may be upgradable.

In this configuration, the test runs as expected.

using QuadGK, Unitful

f(t) = 1.0u"J"
# f (generic function with 1 method)

quadgk(f, 0.0u"s", 10.0u"s")
# (10.0 J s, 0.0 J s)

[mikeingold]

I just found a related bug report on the Unitful.jl package repo.

[stevengj]

The funny thing is that it is failing on:

  [2] one(#unused#::Type{Quantity{Float64, 𝐓, Unitful.FreeUnits{(s,), 𝐓, nothing}}}  )
    @ Base .\number.jl:334

which is this fallback method of one(T). That seems like Unitful's fault — it should really be dispatching to a one method in Unitful.

[stevengj]

I just tried it with the latest versions (QuadGK 2.6.0 and Unitful v1.12.2) in Julia 1.8.2, and it works fine for me:

julia> using Unitful, QuadGK

julia> f(t) = 1.0u"J"
f (generic function with 1 method)

julia> quadgk(f, 0.0u"s", 10.0u"s")
(10.0 J s, 0.0 J s)

So something is weird on your system?

[stevengj]

Hmm, interesting. I can sometimes reproduce your error. If I enter the code above line-by-line, it works. But if I paste it in all at once, then I get the error. This seems like a bug in Julia, TBH.

[stevengj]

Aha, I just noticed that the difference between your code and mine: I was doing using Unitful, QuadGK whereas you were doing using QuadGK, Unitful.

It seems like this is a bug in Julia, since it dispatch shouldn't depend on the order of the imports: JuliaLang/julia#48295

[mikeingold]

Fascinating. I can reproduce the import-order issue that you detected in both a prior install of Julia v1.7.3 and a fresh install of Julia Version 1.9.0-beta2, Commit 7daffeecb8 (2022-12-29 07:45 UTC) with only those two packages added.

[stevengj]

Should be fixed now.

[mikeingold]

Awesome. Thanks!

[giordano]

Maybe could add a test?