QuadGK is type unstable for ULogFloat64

[vandenman]

I noticed a somewhat odd interaction with QuadGK.jl and a type instability for ULogFloat64.

Below is a MWE.

import LogarithmicNumbers, QuadGK

QuadGK.quadgk(x -> x^2, LogarithmicNumbers.LogFloat64(0),  LogarithmicNumbers.LogFloat64(1))
# (+exp(-1.0986122886681098), 0.0)
QuadGK.quadgk(x -> x^2, LogarithmicNumbers.ULogFloat64(0), LogarithmicNumbers.ULogFloat64(1))
# (+exp(-1.0986122886681098), 0.0)

both work, although I noticed that the error term has type Float64 instead of LogFloat64. However,

@code_warntype QuadGK.quadgk(x -> x^2, LogarithmicNumbers.ULogFloat64(0), LogarithmicNumbers.ULogFloat64(1))

shows type instabilities. Specifically,

import LinearAlgebra
f(x) = x^2
segs = (zero(LogarithmicNumbers.ULogFloat64), one(LogarithmicNumbers.ULogFloat64))
atol=nothing; rtol=nothing; maxevals=10^7; order=7; segbuf=nothing; eval_segbuf=nothing
@code_warntype QuadGK.handle_infinities(f, segs) do f, s, _
    QuadGK.do_quadgk(f, s, order, atol, rtol, maxevals, LinearAlgebra.norm, segbuf, eval_segbuf)
end

is type unstable:

MethodInstance for QuadGK.handle_infinities(::var"#19#20", ::typeof(f), ::Tuple{LogarithmicNumbers.ULogFloat64, LogarithmicNumbers.ULogFloat64})
  from handle_infinities(workfunc, f, s) @ QuadGK ~/.julia/packages/QuadGK/wnfvV/src/adapt.jl:163
Arguments
  #self#::Core.Const(QuadGK.handle_infinities)
  workfunc::Core.Const(var"#19#20"())
  f::Core.Const(f)
  s::Tuple{LogarithmicNumbers.ULogFloat64, LogarithmicNumbers.ULogFloat64}
Locals
  #14::QuadGK.var"#14#23"{LogarithmicNumbers.ULogFloat64}
  #13::QuadGK.var"#13#22"
  #12::QuadGK.var"#12#21"{typeof(f), LogarithmicNumbers.ULogFloat64}
  inf2::Bool
  inf1::Bool
  s2::LogarithmicNumbers.ULogFloat64
  s1::LogarithmicNumbers.ULogFloat64
  #20::QuadGK.var"#20#29"{LogarithmicNumbers.ULogFloat64, LogarithmicNumbers.ULogFloat64}
  #19::QuadGK.var"#19#28"{LogarithmicNumbers.ULogFloat64}
  #18::QuadGK.var"#18#27"{typeof(f), LogarithmicNumbers.ULogFloat64, LogarithmicNumbers.ULogFloat64}
  #17::QuadGK.var"#17#26"{LogarithmicNumbers.ULogFloat64, LogarithmicNumbers.ULogFloat64}
  #16::QuadGK.var"#16#25"{LogarithmicNumbers.ULogFloat64}
  #15::QuadGK.var"#15#24"{typeof(f), LogarithmicNumbers.ULogFloat64, LogarithmicNumbers.ULogFloat64}
  @_18::Int64
  si::LogarithmicNumbers.ULogFloat64
  s0::LogarithmicNumbers.ULogFloat64
  @_21::Tuple{LogarithmicNumbers.ULogFloat64, LogarithmicNumbers.ULogFloat64}
Body::Tuple
1 ──        Core.NewvarNode(:(#14))
│           Core.NewvarNode(:(#13))
│           Core.NewvarNode(:(#12))
│           Core.NewvarNode(:(inf2))
│           Core.NewvarNode(:(inf1))
│    %6   = Base.getindex(s, 1)::LogarithmicNumbers.ULogFloat64
│    %7   = Base.lastindex(s)::Core.Const(2)
│    %8   = Base.getindex(s, %7)::LogarithmicNumbers.ULogFloat64
│           (s1 = %6)
│           (s2 = %8)
│    %11  = QuadGK.realone(s1)::Core.Const(true)
└───        goto #16 if not %11
2 ── %13  = QuadGK.realone(s2)::Core.Const(true)
└───        goto #16 if not %13
3 ── %15  = QuadGK.isinf(s1)::Bool
│    %16  = QuadGK.isinf(s2)::Bool
│           (inf1 = %15)
│           (inf2 = %16)
└───        goto #5 if not inf1
4 ──        goto #6
5 ──        goto #16 if not inf2
6 ┄─        goto #9 if not inf1
7 ──        goto #9 if not inf2
8 ── %24  = QuadGK.:(var"#12#21")::Core.Const(QuadGK.var"#12#21")
│    %25  = Core.typeof(f)::Core.Const(typeof(f))
│    %26  = Core.typeof(s1)::Core.Const(LogarithmicNumbers.ULogFloat64)
│    %27  = Core.apply_type(%24, %25, %26)::Core.Const(QuadGK.var"#12#21"{typeof(f), LogarithmicNumbers.ULogFloat64})
│           (#12 = %new(%27, f, s1))
│    %29  = #12::QuadGK.var"#12#21"{typeof(f), LogarithmicNumbers.ULogFloat64}
│           (#13 = %new(QuadGK.:(var"#13#22")))
│    %31  = #13::Core.Const(QuadGK.var"#13#22"())
│    %32  = QuadGK.map(%31, s)::Tuple{Union{LogarithmicNumbers.ULogFloat64, LogarithmicNumbers.LogFloat64}, Union{LogarithmicNumbers.ULogFloat64, LogarithmicNumbers.LogFloat64}}
│    %33  = QuadGK.:(var"#14#23")::Core.Const(QuadGK.var"#14#23")
│    %34  = Core.typeof(s1)::Core.Const(LogarithmicNumbers.ULogFloat64)
│    %35  = Core.apply_type(%33, %34)::Core.Const(QuadGK.var"#14#23"{LogarithmicNumbers.ULogFloat64})
│           (#14 = %new(%35, s1))
│    %37  = #14::QuadGK.var"#14#23"{LogarithmicNumbers.ULogFloat64}
│    %38  = (workfunc)(%29, %32, %37)::Tuple
└───        return %38
9 ┄─        goto #11 if not inf1
10 ─        (@_21 = Core.tuple(s2, s1))
└───        goto #12
11 ─        (@_21 = Core.tuple(s1, s2))
12 ┄ %44  = @_21::Tuple{LogarithmicNumbers.ULogFloat64, LogarithmicNumbers.ULogFloat64}
│           Core.NewvarNode(:(#20))
│           Core.NewvarNode(:(#19))
│           Core.NewvarNode(:(#18))
│           Core.NewvarNode(:(#17))
│           Core.NewvarNode(:(#16))
│           Core.NewvarNode(:(#15))
│    %51  = Base.indexed_iterate(%44, 1)::Core.PartialStruct(Tuple{LogarithmicNumbers.ULogFloat64, Int64}, Any[LogarithmicNumbers.ULogFloat64, Core.Const(2)])
│           (s0 = Core.getfield(%51, 1))
│           (@_18 = Core.getfield(%51, 2))
│    %54  = Base.indexed_iterate(%44, 2, @_18::Core.Const(2))::Core.PartialStruct(Tuple{LogarithmicNumbers.ULogFloat64, Int64}, Any[LogarithmicNumbers.ULogFloat64, Core.Const(3)])
│           (si = Core.getfield(%54, 1))
│    %56  = si::LogarithmicNumbers.ULogFloat64
│    %57  = QuadGK.zero(si)::Core.Const(exp(-Inf))
│    %58  = (%56 < %57)::Bool
└───        goto #15 if not %58
13 ─ %60  = QuadGK.:(var"#15#24")::Core.Const(QuadGK.var"#15#24")
│    %61  = Core.typeof(f)::Core.Const(typeof(f))
│    %62  = Core.typeof(s1)::Core.Const(LogarithmicNumbers.ULogFloat64)
│    %63  = Core.typeof(s0)::Core.Const(LogarithmicNumbers.ULogFloat64)
│    %64  = Core.apply_type(%60, %61, %62, %63)::Core.Const(QuadGK.var"#15#24"{typeof(f), LogarithmicNumbers.ULogFloat64, LogarithmicNumbers.ULogFloat64})
│    %65  = s1::LogarithmicNumbers.ULogFloat64
│           (#15 = %new(%64, f, %65, s0))
│    %67  = #15::QuadGK.var"#15#24"{typeof(f), LogarithmicNumbers.ULogFloat64, LogarithmicNumbers.ULogFloat64}
│    %68  = QuadGK.:(var"#16#25")::Core.Const(QuadGK.var"#16#25")
│    %69  = Core.typeof(s0)::Core.Const(LogarithmicNumbers.ULogFloat64)
│    %70  = Core.apply_type(%68, %69)::Core.Const(QuadGK.var"#16#25"{LogarithmicNumbers.ULogFloat64})
│           (#16 = %new(%70, s0))
│    %72  = #16::QuadGK.var"#16#25"{LogarithmicNumbers.ULogFloat64}
│    %73  = QuadGK.map(%72, s)::Tuple{LogarithmicNumbers.LogFloat64, LogarithmicNumbers.LogFloat64}
│    %74  = QuadGK.reverse(%73)::Tuple{LogarithmicNumbers.LogFloat64, LogarithmicNumbers.LogFloat64}
│    %75  = QuadGK.:(var"#17#26")::Core.Const(QuadGK.var"#17#26")
│    %76  = Core.typeof(s1)::Core.Const(LogarithmicNumbers.ULogFloat64)
│    %77  = Core.typeof(s0)::Core.Const(LogarithmicNumbers.ULogFloat64)
│    %78  = Core.apply_type(%75, %76, %77)::Core.Const(QuadGK.var"#17#26"{LogarithmicNumbers.ULogFloat64, LogarithmicNumbers.ULogFloat64})
│    %79  = s1::LogarithmicNumbers.ULogFloat64
│           (#17 = %new(%78, %79, s0))
│    %81  = #17::QuadGK.var"#17#26"{LogarithmicNumbers.ULogFloat64, LogarithmicNumbers.ULogFloat64}
│    %82  = (workfunc)(%67, %74, %81)::Tuple
└───        return %82
14 ─        Core.Const(:(goto %108))
15 ┄ %85  = QuadGK.:(var"#18#27")::Core.Const(QuadGK.var"#18#27")
│    %86  = Core.typeof(f)::Core.Const(typeof(f))
│    %87  = Core.typeof(s1)::Core.Const(LogarithmicNumbers.ULogFloat64)
│    %88  = Core.typeof(s0)::Core.Const(LogarithmicNumbers.ULogFloat64)
│    %89  = Core.apply_type(%85, %86, %87, %88)::Core.Const(QuadGK.var"#18#27"{typeof(f), LogarithmicNumbers.ULogFloat64, LogarithmicNumbers.ULogFloat64})
│    %90  = s1::LogarithmicNumbers.ULogFloat64
│           (#18 = %new(%89, f, %90, s0))
│    %92  = #18::QuadGK.var"#18#27"{typeof(f), LogarithmicNumbers.ULogFloat64, LogarithmicNumbers.ULogFloat64}
│    %93  = QuadGK.:(var"#19#28")::Core.Const(QuadGK.var"#19#28")
│    %94  = Core.typeof(s0)::Core.Const(LogarithmicNumbers.ULogFloat64)
│    %95  = Core.apply_type(%93, %94)::Core.Const(QuadGK.var"#19#28"{LogarithmicNumbers.ULogFloat64})
│           (#19 = %new(%95, s0))
│    %97  = #19::QuadGK.var"#19#28"{LogarithmicNumbers.ULogFloat64}
│    %98  = QuadGK.map(%97, s)::Tuple{LogarithmicNumbers.LogFloat64, LogarithmicNumbers.LogFloat64}
│    %99  = QuadGK.:(var"#20#29")::Core.Const(QuadGK.var"#20#29")
│    %100 = Core.typeof(s1)::Core.Const(LogarithmicNumbers.ULogFloat64)
│    %101 = Core.typeof(s0)::Core.Const(LogarithmicNumbers.ULogFloat64)
│    %102 = Core.apply_type(%99, %100, %101)::Core.Const(QuadGK.var"#20#29"{LogarithmicNumbers.ULogFloat64, LogarithmicNumbers.ULogFloat64})
│    %103 = s1::LogarithmicNumbers.ULogFloat64
│           (#20 = %new(%102, %103, s0))
│    %105 = #20::QuadGK.var"#20#29"{LogarithmicNumbers.ULogFloat64, LogarithmicNumbers.ULogFloat64}
│    %106 = (workfunc)(%92, %98, %105)::Tuple
└───        return %106
16 ┄ %108 = (workfunc)(f, s, QuadGK.identity)::Tuple
└───        return %108

the highlighting on GitHub doesn't quite show it but all output from workfunc returns an untyped Tuple.

I played around with some ways to resolve this but it's not entirely clear to me what the best way to fix this would be. In addition, I had a few questions.

1. Why is it that Logarithmic <: Real, instead of Logarithmic <: AbstractFloat?
2. Why isn't float defined as float(x::Logarithmic) = x?

As a learning project, I ported the R package Brobdingag to Julia (which does essentially the same as this package). There, I did make those decisions and then QuadGK worked without issues, although I'm sure there are reasons for doing otherwise.