JianghuiDu
JianghuiDu

Reputation: 319

Julia function type instability

Still a bit confused about to to make functions type stable in Julia. Here is my function

function tran7(C::Array{Float64,1}, flux_up::Float64, D::Array{Float64,1},
    v::Array{Float64,1}, AFDW::Array{Float64,1}, VF::Array{Float64,1},
    VF_mid::Array{Float64,1})

    C_up = Float64[];
    flux = Array{Float64,1}(N+1); dC = Array{Float64, 1}(N)

    C_up = (flux_up + VF[1]*(D[1]/dx_aux[1] + (1-AFDW[1])*v[1])*C[1])/
    (VF[1]*(D[1]/dx_aux[1] + AFDW[1]*v[1]))

    flux = -VF.*D.*(vcat(C_up,C,C[N])[2:(N+2)] -
            vcat(C_up,C,C[N])[1:(N+1)])./dx_aux +
            VF.*v.*(AFDW.*vcat(C_up,C) + (1-AFDW).*vcat(C,C[N]))
    flux[1] = flux_up
    dC = -(flux[2:(N+1)] - flux[1:N])./(VF_mid.*dx)
    return dC
end

And here is the @code_warntype results

Variables:
  #self# <optimized out>
  C::Array{Float64,1}
  flux_up::Float64
  D::Array{Float64,1}
  v::Array{Float64,1}
  AFDW::Array{Float64,1}
  VF::Array{Float64,1}
  VF_mid::Array{Float64,1}
  #200::##200#203
  #201::##201#204
  #202::##202#205
  C_down <optimized out>
  C_up::ANY
  flux::ANY
  dC::ANY
  #temp#@_16::Core.MethodInstance
  #temp#@_17::ANY
  #temp#@_18::Core.MethodInstance
  #temp#@_19::ANY
  #temp#@_20 <optimized out>

Body:
  begin 
      C_up::ANY = $(Expr(:foreigncall, :(:jl_alloc_array_1d), 
Array{Float64,1}, svec(Any, Int64), Array{Float64,1}, 0, 0, 0)) # line 5:
      $(Expr(:foreigncall, :(:jl_alloc_array_1d), Array{Float64,1}, svec(Any, Int64), Array{Float64,1}, 0, 0, 0)) # line 6:
  flux::ANY = (Array{Float64,1})((Main.N + 1)::ANY)::Array{Float64,1} # line 6:
  dC::ANY = (Array{Float64,1})(Main.N)::Array{Float64,1} # line 8:
  SSAValue(21) = (Base.arrayref)(VF::Array{Float64,1}, 1)::Float64
  SSAValue(20) = (((Base.arrayref)(D::Array{Float64,1}, 1)::Float64 / (Main.getindex)(Main.dx_aux, 1)::ANY)::ANY + (Base.mul_float)((Base.sub_float)((Base.sitofp)(Float64, 1)::Float64, (Base.arrayref)(AFDW::Array{Float64,1}, 1)::Float64)::Float64, (Base.arrayref)(v::Array{Float64,1}, 1)::Float64)::Float64)::ANY
  SSAValue(19) = (Base.arrayref)(C::Array{Float64,1}, 1)::Float64
  $(Expr(:inbounds, false))
  # meta: location operators.jl * 424
  SSAValue(22) = ((SSAValue(21) * SSAValue(20))::ANY * SSAValue(19))::ANY
  # meta: pop location
  $(Expr(:inbounds, :pop))
  C_up::ANY = ((flux_up::Float64 + SSAValue(22))::ANY / ((Base.arrayref)(VF::Array{Float64,1}, 1)::Float64 * (((Base.arrayref)(D::Array{Float64,1}, 1)::Float64 / (Main.getindex)(Main.dx_aux, 1)::ANY)::ANY + (Base.mul_float)((Base.arrayref)(AFDW::Array{Float64,1}, 1)::Float64, (Base.arrayref)(v::Array{Float64,1}, 1)::Float64)::Float64)::ANY)::ANY)::ANY # line 11:
  #200::##200#203 = $(Expr(:new, :(Main.##200#203)))
  SSAValue(1) = $(Expr(:invoke, MethodInstance for -(::Array{Float64,1}), :(Main.-), :(VF)))
  SSAValue(3) = (Main.vcat)(C_up::ANY, C::Array{Float64,1}, (Main.getindex)(C::Array{Float64,1}, Main.N)::ANY)::ANY
  SSAValue(4) = (Main.getindex)(SSAValue(3), (Main.colon)(2, (Main.N + 2)::ANY)::ANY)::ANY
  SSAValue(5) = (Main.vcat)(C_up::ANY, C::Array{Float64,1}, (Main.getindex)(C::Array{Float64,1}, Main.N)::ANY)::ANY
  SSAValue(6) = (Main.getindex)(SSAValue(5), (Main.colon)(1, (Main.N + 1)::ANY)::ANY)::ANY
  SSAValue(7) = (SSAValue(4) - SSAValue(6))::ANY
  SSAValue(33) = D::Array{Float64,1}
  SSAValue(34) = Main.dx_aux
  $(Expr(:inbounds, false))
  # meta: location broadcast.jl broadcast 434
  SSAValue(26) = SSAValue(7)
  SSAValue(27) = SSAValue(34)
  # meta: location broadcast.jl containertype 34
  SSAValue(25) = (Base.Broadcast.promote_containertype)(Array, (Base.Broadcast.promote_containertype)((Base.Broadcast._containertype)((Base.Broadcast.typeof)(SSAValue(26))::DataType)::ANY, (Base.Broadcast._containertype)((Base.Broadcast.typeof)(SSAValue(27))::DataType)::ANY)::ANY)::UNION{TYPE{ARRAY}, TYPE{BASE.SPARSEARRAYS.HIGHERORDERFNS.PROMOTETOSPARSE}}
  unless (SSAValue(25) isa Type{Array})::Bool goto 37
  #temp#@_16::Core.MethodInstance = MethodInstance for promote_containertype(::Type{Array}, ::Type{Array})
  goto 46
  37: 
  unless (SSAValue(25) isa Type{Base.SparseArrays.HigherOrderFns.PromoteToSparse})::Bool goto 41
  #temp#@_16::Core.MethodInstance = MethodInstance for promote_containertype(::Type{Array}, ::Type{Base.SparseArrays.HigherOrderFns.PromoteToSparse})
  goto 46
  41: 
  goto 43
  43: 
  #temp#@_17::ANY = (Base.Broadcast.promote_containertype)(Array, SSAValue(25))::UNION{TYPE{ARRAY}, TYPE{BASE.SPARSEARRAYS.HIGHERORDERFNS.PROMOTETOSPARSE}}
  goto 48
  46: 
  #temp#@_17::ANY = $(Expr(:invoke, :(#temp#@_16), :(Base.Broadcast.promote_containertype), Array, SSAValue(25)))
  48: 
  # meta: pop location
  # meta: pop location
  $(Expr(:inbounds, :pop))
  SSAValue(8) = (Base.Broadcast.broadcast_c)(#200::##200#203, #temp#@_17::ANY, SSAValue(1), SSAValue(33), SSAValue(7), SSAValue(34))::ANY
  #201::##201#204 = $(Expr(:new, :(Main.##201#204)))
  SSAValue(12) = ((Base.broadcast)(Main.*, AFDW::Array{Float64,1}, (Main.vcat)(C_up::ANY, C::Array{Float64,1})::ANY)::ANY + (Base.broadcast)(Main.*, $(Expr(:invoke, MethodInstance for -(::Int64, ::Array{Float64,1}), :(Main.-), 1, :(AFDW))), (Main.vcat)(C::Array{Float64,1}, (Main.getindex)(C::Array{Float64,1}, Main.N)::ANY)::ANY)::ANY)::ANY
  $(Expr(:inbounds, false))
  # meta: location broadcast.jl broadcast 434
  # meta: location broadcast.jl containertype 34
  SSAValue(28) = (Base.Broadcast.promote_containertype)(Array, (Base.Broadcast._containertype)((Base.Broadcast.typeof)(SSAValue(12))::DataType)::ANY)::UNION{TYPE{ARRAY}, TYPE{BASE.SPARSEARRAYS.HIGHERORDERFNS.PROMOTETOSPARSE}}
  unless (SSAValue(28) isa Type{Array})::Bool goto 62
  #temp#@_18::Core.MethodInstance = MethodInstance for promote_containertype(::Type{Array}, ::Type{Array})
  goto 71
  62: 
  unless (SSAValue(28) isa Type{Base.SparseArrays.HigherOrderFns.PromoteToSparse})::Bool goto 66
  #temp#@_18::Core.MethodInstance = MethodInstance for promote_containertype(::Type{Array}, ::Type{Base.SparseArrays.HigherOrderFns.PromoteToSparse})
  goto 71
  66: 
  goto 68
  68: 
  #temp#@_19::ANY = (Base.Broadcast.promote_containertype)(Array, SSAValue(28))::UNION{TYPE{ARRAY}, TYPE{BASE.SPARSEARRAYS.HIGHERORDERFNS.PROMOTETOSPARSE}}
  goto 73
  71: 
  #temp#@_19::ANY = $(Expr(:invoke, :(#temp#@_18), :(Base.Broadcast.promote_containertype), Array, SSAValue(28)))
  73: 
  # meta: pop location
  # meta: pop location
  $(Expr(:inbounds, :pop))
  SSAValue(13) = (Base.Broadcast.broadcast_c)(#201::##201#204, #temp#@_19::ANY, VF::Array{Float64,1}, v::Array{Float64,1}, SSAValue(12))::ANY
  flux::ANY = (SSAValue(8) + SSAValue(13))::ANY # line 14:
  (Main.setindex!)(flux::ANY, flux_up::Float64, 1)::ANY # line 15:
  #202::##202#205 = $(Expr(:new, :(Main.##202#205)))
  SSAValue(15) = -(((Main.getindex)(flux::ANY, (Main.colon)(2, (Main.N + 1)::ANY)::ANY)::ANY - (Main.getindex)(flux::ANY, (Main.colon)(1, Main.N)::ANY)::ANY)::ANY)::ANY
  SSAValue(35) = VF_mid::Array{Float64,1}
  SSAValue(36) = Main.dx
  $(Expr(:inbounds, false))
  # meta: location broadcast.jl broadcast 434
  SSAValue(32) = (Base.Broadcast.promote_containertype)((Base.Broadcast._containertype)((Base.Broadcast.typeof)(SSAValue(15))::DataType)::ANY, (Base.Broadcast.promote_containertype)(Array, (Base.Broadcast._containertype)((Base.Broadcast.typeof)(SSAValue(36))::DataType)::ANY)::UNION{TYPE{ARRAY}, TYPE{BASE.SPARSEARRAYS.HIGHERORDERFNS.PROMOTETOSPARSE}})::ANY
  # meta: pop location
  $(Expr(:inbounds, :pop))
  dC::ANY = (Base.Broadcast.broadcast_c)(#202::##202#205, SSAValue(32), SSAValue(15), SSAValue(35), SSAValue(36))::ANY # line 16:
  return dC::ANY
  end::ANY

I have a few questions. (1) why do the variables that were preallocated (dC, C_up, flux) still have type any? (2) if I have a parameter (constant) defined outside of the function, such as dx_aux in here, is there an advantage of including it as an input? (3) what's the best way to deal with complex intermediary computational results?

Upvotes: 1

Views: 490

Answers (1)

Chris Rackauckas
Chris Rackauckas

Reputation: 19132

C_up, flux, and dC are not defined in the function header, so I assume they are globals. Globals can always change types because

x = 3
x = 2.0

works in the REPL, so the functions do not specialize on the type of global variables out of fear of getting the type wrong. However, if it's a constant,

const x = 3

then the compiler can make the assumption it won't change, and thus specialize on the type and fix this issue. Note that the value of arrays is their reference, so

const x = [3,4,5]
x[2] = 6

is fine because the value of x is just the reference, so it doesn't change even when internal values of the array change.

(3) what's the best way to deal with complex intermediary computational results?

Cache arrays like you have here, or callable types with a cache, or ..., this is far too broad to answer well.

I'll leave with two links:

http://www.stochasticlifestyle.com/7-julia-gotchas-handle/

That describes that global variable problem and why it happens. But also:

https://github.com/JuliaLang/julia/issues/8870

suggests a new syntax:

x::Int = 3

so that non-constant globals can have a declared type that is inferred in further functions. This seems like it won't make 1.0, but it will be a 1.x.

Upvotes: 3

Related Questions