Why isn't this function type stable?

Let’s consider the function:

using Sobol

function fun_pdf(img1::Matrix{T}, n_skip) where T <: Real
    (nrows, ncols) = size(img1)
    
    indx = SobolSeq([0, 0], [ncols - 1, nrows - 1])
    indx = skip(indx, n_skip)

    n_samples = Int(floor(nrows * nrows * 0.1))

    for i in 1:n_samples
        ix = next!(indx)
        #xq, yq = ceil(Int, ix[1]), ceil(Int, ix[2])
        xq, yq = ceil(Int64, ix[1]), ceil(Int64, ix[2])
    end
end

Even with

xq, yq = ceil(Int64, ix[1]), ceil(Int64, ix[2])

instead of

xq, yq = ceil(Int, ix[1]), ceil(Int, ix[2])

the compiler can’t infer the type of xq and yq:

PDF = zeros(Float64,nf1,nf2)
@code_warntype  fun_pdf(PDF, 100)

Variables
  #self#::Core.Const(fun_pdf)
  img1::Matrix{Float64}
  n_skip::Int64
  @_4::Union{Nothing, Tuple{Int64, Int64}}
  @_5::Int64
  n_samples::Int64
  indx::ScaledSobolSeq{_A, Int64} where _A
  ncols::Int64
  nrows::Int64
  i::Int64
  yq::Any
  xq::Any
  ix::AbstractVector{var"#s7"} where var"#s7"<:AbstractFloat

Body::Nothing
1 ─ %1  = Main.size(img1)::Tuple{Int64, Int64}
β”‚   %2  = Base.indexed_iterate(%1, 1)::Core.PartialStruct(Tuple{Int64, Int64}, Any[Int64, Core.Const(2)])
β”‚         (nrows = Core.getfield(%2, 1))
β”‚         (@_5 = Core.getfield(%2, 2))
β”‚   %5  = Base.indexed_iterate(%1, 2, @_5::Core.Const(2))::Core.PartialStruct(Tuple{Int64, Int64}, Any[Int64, Core.Const(3)])
β”‚         (ncols = Core.getfield(%5, 1))
β”‚   %7  = Base.vect(0, 0)::Vector{Int64}
β”‚   %8  = (ncols - 1)::Int64
β”‚   %9  = (nrows - 1)::Int64
β”‚   %10 = Base.vect(%8, %9)::Vector{Int64}
β”‚         (indx = Main.SobolSeq(%7, %10))
β”‚         (indx = Main.skip(indx, n_skip))
β”‚   %13 = (nrows * nrows * 0.1)::Float64
β”‚   %14 = Main.floor(%13)::Float64
β”‚         (n_samples = Main.Int(%14))
β”‚   %16 = (1:n_samples)::Core.PartialStruct(UnitRange{Int64}, Any[Core.Const(1), Int64])
β”‚         (@_4 = Base.iterate(%16))
β”‚   %18 = (@_4 === nothing)::Bool
β”‚   %19 = Base.not_int(%18)::Bool
└──       goto #4 if not %19
2 β”„ %21 = @_4::Tuple{Int64, Int64}::Tuple{Int64, Int64}
β”‚         (i = Core.getfield(%21, 1))
β”‚   %23 = Core.getfield(%21, 2)::Int64
β”‚         (ix = Main.next!(indx))
β”‚   %25 = Base.getindex(ix, 1)::Any
β”‚   %26 = Main.ceil(Main.Int64, %25)::Any
β”‚   %27 = Base.getindex(ix, 2)::Any
β”‚   %28 = Main.ceil(Main.Int64, %27)::Any
β”‚         (xq = %26)
β”‚         (yq = %28)
β”‚         (@_4 = Base.iterate(%16, %23))
β”‚   %32 = (@_4 === nothing)::Bool
β”‚   %33 = Base.not_int(%32)::Bool
└──       goto #4 if not %33
3 ─       goto #2
4 β”„       return nothing

Any ideas?

The ScaledSobolSeq object you’re constructing needs to know the dimension (as part of it’s type information).

You can either:

  • construct indx outside the function, and use a function barrier, or
  • provide the type information directly yourself. I.e, write
    indx = ScaledSobolSeq{2, Float64}(Float64[0,0], Float64[ncols -1, nrows-1])
    

of course, there may be another option that I’m missing!

3 Likes