# How to investigate output of type instability checks?

I am trying to eliminate the type instabilities that might cause my code (shown below) of an iterative solver to some non linear system of equations to become slow:

``````using Printf,Plots, TimerOutputs
include("ccd-helper.jl")
function ccd_by_hand(maxitr)
initialize_cc_variables()
nv::Int64 = deserialize("nv.jlbin")
nocc::Int64 = deserialize("nocc.jlbin")
erhf::Float64 = deserialize("erhf.jlbin")
R2 = TensorMap(zeros(Float64,nv,nv,nocc,nocc),ℝ^nv ⊗ ℝ^nv , ℝ^nocc ⊗ ℝ^nocc)
println("R2 is of type \$(typeof(R2))")
R2u = TensorMap(zeros(Float64,nv,nv,nocc,nocc),ℝ^nv ⊗ ℝ^nv , ℝ^nocc ⊗ ℝ^nocc)
R_iter = copy(R2)
Scaled_R2 = copy(R2)
T2_old::TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}= initialize_t2_only()
println("T2 is of type \$(typeof(T2_old))")
normtol=1.0e-8
e_new::Float64 = calculate_ECCD(T2_old)
e_old::Float64 = copy(0.0)
etol = 1.0e-10
p_max =  6# Maximum number of previous iterations to store for DIIS
p_min = 2 # DIIS will start after this many iterations
R_iter_storage = Array{typeof(R2)}(undef, 0)
T2_storage::Vector{typeof(T2_old)} = Array{typeof(T2_old)}(undef, 0)
R2_storage = Array{typeof(R2)}(undef, 0)
println("R2_storage is of type \$(typeof(R2_storage))")
println("Starting CCD Iterations with Convergence Criteria as ||R₂||<\$normtol , ΔE<\$etol and max iterations=\$maxitr")
println("-----------------------------------------------------------------------------------------------")
for i in 1:maxitr
@printf("\n-------------------------\nStarting Iteration: %d  Current Total Energy: %.8f\n", i, e_new+erhf)
println("Energy before updating = \$(e_old+erhf)")
# R2 = calculate_residual(T2_old,R2u,R2)
R2 = calculate_residual_memeff(T2_old)
println("R2 is of type \$(typeof(R2))")
T2_old = update_T2(T2_old,R2,Scaled_R2)
println("T2 is of type \$(typeof(T2_old))")
conv::Bool,r2norm::Float64 = check_convergence(R2,normtol,e_old,e_new,etol)
if conv
@printf("CCD Converged in %d iterations, no further updates required,current ||R2||=%.10f and ΔE = %.12f\n\n-------------------------\nE_Total = %.8f\n", i,r2norm,abs(e_old-e_new),e_new+erhf)
break
end
R_iter = calculate_R_iter(R_iter,R2)
@printf("CCD Not Converged in %d iterations, current ||R2|| = %.10f and ΔE = %.12f \n", i,r2norm,abs(e_old-e_new))
if i >= p_min   #DIIS  starts
println("DIIS is being implemented in iteration \$i")
# if i >= p_min + p_max - 1
if length(R_iter_storage)>=p_max  # Adjusted condition to start popping correctly
popfirst!(R_iter_storage)
popfirst!(T2_storage)
popfirst!(R2_storage)
end
push!(R_iter_storage, copy(R_iter))
push!(T2_storage, copy(T2_old))
push!(R2_storage, copy(R2))
#update T2 via DIIS
p::Int64 = length(R_iter_storage)
diis_result = PerformDiis(R_iter_storage,p,T2_storage,R2_storage)
if diis_result==false
# T2_old = update_T2(T2_old,R2,fvv,foo)
nothing
else
T2_old = diis_result
# println("DIIS has updated T2")
# display(T2_old)
end
elseif i == p_min - 1 #DIIS starts next iteration so add current stuff to memory
println("DIIS will start in next iteration")
push!(R_iter_storage, copy(R_iter))
push!(T2_storage, copy(T2_old))
push!(R2_storage, copy(R2))
just_show_Bmatrix(R_iter_storage,length(R_iter_storage),T2_storage,R2_storage)
else
println("DIIS start is still far away")
# No change to memory
end
e_old = e_new
e_new = calculate_ECCD(T2_old)
println("Energy after updating = \$(e_new+erhf)\n-------------------------\n")
end

end
maxitr = 100;
@code_warntype ccd_by_hand(maxitr);
ccd_by_hand(maxitr)
``````

On running the `@code_warntype` macro on the function `ccd_by_hand(..)` , I get the following output :

``````MethodInstance for ccd_by_hand(::Int64)
from ccd_by_hand(maxitr) @ Main ~/Desktop/Andreas Köhn/phase 1.5/code/tensorkit/ccd-memeff.jl:3
Arguments
#self#::Core.Const(ccd_by_hand)
maxitr::Int64
Locals
@_3::Union{Nothing, Tuple{Int64, Int64}}
R2_storage::Vector{TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}}
T2_storage::Vector{TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}}
R_iter_storage::Vector{TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}}
p_min::Int64
p_max::Int64
etol::Float64
e_old::Float64
e_new::Float64
normtol::Float64
T2_old::TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}
Scaled_R2::TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}
R_iter::TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}
R2u::TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}
R2::TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}
erhf::Float64
nocc::Int64
nv::Int64
@_21::Int64
i::Int64
diis_result::TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}
p::Int64
r2norm::Float64
conv::Bool
@_27::Any
@_28::Any
@_29::Any
@_30::TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}
@_31::Any
@_32::Float64
@_33::Vector{TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}}
@_34::TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}
@_35::Bool
@_36::Any
@_37::Int64
@_38::TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}
@_39::Float64
@_40::Any
Body::Nothing
1 ──        Core.NewvarNode(:(@_3))
│           Core.NewvarNode(:(R2_storage))
│           Core.NewvarNode(:(T2_storage))
│           Core.NewvarNode(:(R_iter_storage))
│           Core.NewvarNode(:(p_min))
│           Core.NewvarNode(:(p_max))
│           Core.NewvarNode(:(etol))
│           Core.NewvarNode(:(e_old))
│           Core.NewvarNode(:(e_new))
│           Core.NewvarNode(:(normtol))
│           Core.NewvarNode(:(T2_old))
│           Core.NewvarNode(:(Scaled_R2))
│           Core.NewvarNode(:(R_iter))
│           Core.NewvarNode(:(R2u))
│           Core.NewvarNode(:(R2))
│           Core.NewvarNode(:(erhf))
│           Core.NewvarNode(:(nocc))
│           Core.NewvarNode(:(nv))
│           Main.initialize_cc_variables()
│    %20  = Main.deserialize("nv.jlbin")::Any
│           (@_27 = %20)
│    %22  = (@_27 isa Main.Int64)::Bool
└───        goto #3 if not %22
2 ──        goto #4
3 ── %25  = Base.convert(Main.Int64, @_27)::Any
└───        (@_27 = Core.typeassert(%25, Main.Int64))
4 ┄─        (nv = @_27::Int64)
│    %28  = Main.deserialize("nocc.jlbin")::Any
│           (@_28 = %28)
│    %30  = (@_28 isa Main.Int64)::Bool
└───        goto #6 if not %30
5 ──        goto #7
6 ── %33  = Base.convert(Main.Int64, @_28)::Any
└───        (@_28 = Core.typeassert(%33, Main.Int64))
7 ┄─        (nocc = @_28::Int64)
│    %36  = Main.deserialize("erhf.jlbin")::Any
│           (@_29 = %36)
│    %38  = (@_29 isa Main.Float64)::Bool
└───        goto #9 if not %38
8 ──        goto #10
9 ── %41  = Base.convert(Main.Float64, @_29)::Any
└───        (@_29 = Core.typeassert(%41, Main.Float64))
10 ┄        (erhf = @_29::Float64)
│    %44  = Main.zeros(Main.Float64, nv, nv, nocc, nocc)::Array{Float64, 4}
│    %45  = (Main.ℝ ^ nv)::CartesianSpace
│    %46  = (Main.ℝ ^ nv)::CartesianSpace
│    %47  = (%45 ⊗ %46)::ProductSpace{CartesianSpace, 2}
│    %48  = (Main.ℝ ^ nocc)::CartesianSpace
│    %49  = (Main.ℝ ^ nocc)::CartesianSpace
│    %50  = (%48 ⊗ %49)::ProductSpace{CartesianSpace, 2}
│           (R2 = Main.TensorMap(%44, %47, %50))
│    %52  = Main.typeof(R2)::Core.Const(TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}})
│    %53  = Base.string("R2 is of type ", %52)::String
│           Main.println(%53)
│    %55  = Main.zeros(Main.Float64, nv, nv, nocc, nocc)::Array{Float64, 4}
│    %56  = (Main.ℝ ^ nv)::CartesianSpace
│    %57  = (Main.ℝ ^ nv)::CartesianSpace
│    %58  = (%56 ⊗ %57)::ProductSpace{CartesianSpace, 2}
│    %59  = (Main.ℝ ^ nocc)::CartesianSpace
│    %60  = (Main.ℝ ^ nocc)::CartesianSpace
│    %61  = (%59 ⊗ %60)::ProductSpace{CartesianSpace, 2}
│           (R2u = Main.TensorMap(%55, %58, %61))
│           (R_iter = Main.copy(R2))
│           (Scaled_R2 = Main.copy(R2))
│    %65  = Main.initialize_t2_only()::TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}
│    %66  = Main.TrivialTensorMap::Core.Const(TrivialTensorMap)
│    %67  = Main.CartesianSpace::Core.Const(CartesianSpace)
│    %68  = Core.apply_type(Main.Matrix, Main.Float64)::Core.Const(Matrix{Float64})
│    %69  = Core.apply_type(%66, %67, 2, 2, %68)::Core.Const(TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}})
│           (@_30 = %65)
│    %71  = (@_30 isa %69)::Core.Const(true)
└───        goto #12 if not %71
11 ─        goto #13
12 ─        Core.Const(:(Base.convert(%69, @_30)))
└───        Core.Const(:(@_30 = Core.typeassert(%74, %69)))
13 ┄        (T2_old = @_30)
│    %77  = Main.typeof(T2_old)::Core.Const(TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}})
│    %78  = Base.string("T2 is of type ", %77)::String
│           Main.println(%78)
│           (normtol = 1.0e-8)
│    %81  = Main.calculate_ECCD(T2_old)::Any
│           (@_31 = %81)
│    %83  = (@_31 isa Main.Float64)::Bool
└───        goto #15 if not %83
14 ─        goto #16
15 ─ %86  = Base.convert(Main.Float64, @_31)::Any
└───        (@_31 = Core.typeassert(%86, Main.Float64))
16 ┄        (e_new = @_31::Float64)
│    %89  = Main.copy(0.0)::Core.Const(0.0)
│           (@_32 = %89)
│    %91  = (@_32::Core.Const(0.0) isa Main.Float64)::Core.Const(true)
└───        goto #18 if not %91
17 ─        goto #19
18 ─        Core.Const(:(Base.convert(Main.Float64, @_32)))
└───        Core.Const(:(@_32 = Core.typeassert(%94, Main.Float64)))
19 ┄        (e_old = @_32::Core.Const(0.0))
│           (etol = 1.0e-10)
│           (p_max = 6)
│           (p_min = 2)
│    %100 = Main.Array::Core.Const(Array)
│    %101 = Main.typeof(R2)::Core.Const(TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}})
│    %102 = Core.apply_type(%100, %101)::Core.Const(Array{TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}})
│    %103 = Main.undef::Core.Const(UndefInitializer())
│           (R_iter_storage = (%102)(%103, 0))
│    %105 = Main.Array::Core.Const(Array)
│    %106 = Main.typeof(T2_old)::Core.Const(TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}})
│    %107 = Core.apply_type(%105, %106)::Core.Const(Array{TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}})
│    %108 = Main.undef::Core.Const(UndefInitializer())
│    %109 = (%107)(%108, 0)::Vector{TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}}
│    %110 = Main.Vector::Core.Const(Vector)
│    %111 = Main.typeof(T2_old)::Core.Const(TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}})
│    %112 = Core.apply_type(%110, %111)::Core.Const(Vector{TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}})
│           (@_33 = %109)
│    %114 = (@_33 isa %112)::Core.Const(true)
└───        goto #21 if not %114
20 ─        goto #22
21 ─        Core.Const(:(Base.convert(%112, @_33)))
└───        Core.Const(:(@_33 = Core.typeassert(%117, %112)))
22 ┄        (T2_storage = @_33)
│    %120 = Main.Array::Core.Const(Array)
│    %121 = Main.typeof(R2)::Core.Const(TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}})
│    %122 = Core.apply_type(%120, %121)::Core.Const(Array{TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}})
│    %123 = Main.undef::Core.Const(UndefInitializer())
│           (R2_storage = (%122)(%123, 0))
│    %125 = Main.typeof(R2_storage)::Core.Const(Vector{TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}})
│    %126 = Base.string("R2_storage is of type ", %125)::String
│           Main.println(%126)
│    %128 = Base.string("Starting CCD Iterations with Convergence Criteria as ||R₂||<", normtol::Core.Const(1.0e-8), " , ΔE<", etol::Core.Const(1.0e-10), " and max iterations=", maxitr)::String
│           Main.println(%128)
│           Main.println("-----------------------------------------------------------------------------------------------")
│    %131 = (1:maxitr)::Core.PartialStruct(UnitRange{Int64}, Any[Core.Const(1), Int64])
│           (@_3 = Base.iterate(%131))
│    %133 = (@_3 === nothing)::Bool
│    %134 = Base.not_int(%133)::Bool
└───        goto #57 if not %134
23 ┄        Core.NewvarNode(:(@_21))
│           Core.NewvarNode(:(diis_result))
│           Core.NewvarNode(:(p))
│           Core.NewvarNode(:(r2norm))
│           Core.NewvarNode(:(conv))
│    %141 = @_3::Tuple{Int64, Int64}
│           (i = Core.getfield(%141, 1))
│    %143 = Core.getfield(%141, 2)::Int64
│    %144 = Base.getproperty(Printf, :format)::Core.Const(Printf.format)
│    %145 = Main.stdout::IO
│    %146 = i::Int64
│    %147 = (e_new + erhf)::Float64
│           (%144)(%145, Printf.Format{Base.CodeUnits{UInt8, String}, Tuple{Printf.Spec{Val{'d'}}, Printf.Spec{Val{'f'}}}}(UInt8[0x0a, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x0a, 0x53, 0x74, 0x61, 0x72, 0x74, 0x69, 0x6e, 0x67, 0x20, 0x49, 0x74, 0x65, 0x72, 0x61, 0x74, 0x69, 0x6f, 0x6e, 0x3a, 0x20, 0x25, 0x64, 0x20, 0x20, 0x43, 0x75, 0x72, 0x72, 0x65, 0x6e, 0x74, 0x20, 0x54, 0x6f, 0x74, 0x61, 0x6c, 0x20, 0x45, 0x6e, 0x65, 0x72, 0x67, 0x79, 0x3a, 0x20, 0x25, 0x2e, 0x38, 0x66, 0x0a], UnitRange{Int64}[1:47, 50:73, 78:78], (%.0d, %.8f), 2), %146, %147)
│    %149 = (e_old + erhf)::Float64
│    %150 = Base.string("Energy before updating = ", %149)::String
│           Main.println(%150)
│           (R2 = Main.calculate_residual_memeff(T2_old))
│    %153 = Main.typeof(R2)::Core.Const(TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}})
│    %154 = Base.string("R2 is of type ", %153)::String
│           Main.println(%154)
│    %156 = Main.update_T2(T2_old, R2, Scaled_R2)::TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}
│    %157 = Main.TrivialTensorMap::Core.Const(TrivialTensorMap)
│    %158 = Main.CartesianSpace::Core.Const(CartesianSpace)
│    %159 = Core.apply_type(Main.Matrix, Main.Float64)::Core.Const(Matrix{Float64})
│    %160 = Core.apply_type(%157, %158, 2, 2, %159)::Core.Const(TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}})
│           (@_34 = %156)
│    %162 = (@_34 isa %160)::Core.Const(true)
└───        goto #25 if not %162
24 ─        goto #26
25 ─        Core.Const(:(Base.convert(%160, @_34)))
└───        Core.Const(:(@_34 = Core.typeassert(%165, %160)))
26 ┄        (T2_old = @_34)
│    %168 = Main.typeof(T2_old)::Core.Const(TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}})
│    %169 = Base.string("T2 is of type ", %168)::String
│           Main.println(%169)
│    %171 = Main.check_convergence(R2, normtol::Core.Const(1.0e-8), e_old, e_new, etol::Core.Const(1.0e-10))::Tuple{Bool, Any}
│    %172 = Base.indexed_iterate(%171, 1)::Core.PartialStruct(Tuple{Bool, Int64}, Any[Bool, Core.Const(2)])
│    %173 = Core.getfield(%172, 1)::Bool
│           (@_21 = Core.getfield(%172, 2))
│    %175 = Base.indexed_iterate(%171, 2, @_21::Core.Const(2))::Core.PartialStruct(Tuple{Any, Int64}, Any[Any, Core.Const(3)])
│    %176 = Core.getfield(%175, 1)::Any
│           (@_35 = %173)
│    %178 = (@_35 isa Main.Bool)::Core.Const(true)
└───        goto #28 if not %178
27 ─        goto #29
28 ─        Core.Const(:(Base.convert(Main.Bool, @_35)))
└───        Core.Const(:(@_35 = Core.typeassert(%181, Main.Bool)))
29 ┄        (conv = @_35)
│           (@_36 = %176)
│    %185 = (@_36 isa Main.Float64)::Bool
└───        goto #31 if not %185
30 ─        goto #32
31 ─ %188 = Base.convert(Main.Float64, @_36)::Any
└───        (@_36 = Core.typeassert(%188, Main.Float64))
32 ┄        (r2norm = @_36::Float64)
└───        goto #34 if not conv
33 ─ %192 = Base.getproperty(Printf, :format)::Core.Const(Printf.format)
│    %193 = Main.stdout::IO
│    %194 = i::Int64
│    %195 = r2norm::Float64
│    %196 = (e_old - e_new)::Float64
│    %197 = Main.abs(%196)::Float64
│    %198 = (e_new + erhf)::Float64
│           (%192)(%193, Printf.Format{Base.CodeUnits{UInt8, String}, Tuple{Printf.Spec{Val{'d'}}, Printf.Spec{Val{'f'}}, Printf.Spec{Val{'f'}}, Printf.Spec{Val{'f'}}}}(UInt8[0x43, 0x43, 0x44, 0x20, 0x43, 0x6f, 0x6e, 0x76, 0x65, 0x72, 0x67, 0x65, 0x64, 0x20, 0x69, 0x6e, 0x20, 0x25, 0x64, 0x20, 0x69, 0x74, 0x65, 0x72, 0x61, 0x74, 0x69, 0x6f, 0x6e, 0x73, 0x2c, 0x20, 0x6e, 0x6f, 0x20, 0x66, 0x75, 0x72, 0x74, 0x68, 0x65, 0x72, 0x20, 0x75, 0x70, 0x64, 0x61, 0x74, 0x65, 0x73, 0x20, 0x72, 0x65, 0x71, 0x75, 0x69, 0x72, 0x65, 0x64, 0x2c, 0x63, 0x75, 0x72, 0x72, 0x65, 0x6e, 0x74, 0x20, 0x7c, 0x7c, 0x52, 0x32, 0x7c, 0x7c, 0x3d, 0x25, 0x2e, 0x31, 0x30, 0x66, 0x20, 0x61, 0x6e, 0x64, 0x20, 0xce, 0x94, 0x45, 0x20, 0x3d, 0x20, 0x25, 0x2e, 0x31, 0x32, 0x66, 0x0a, 0x0a, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x0a, 0x45, 0x5f, 0x54, 0x6f, 0x74, 0x61, 0x6c, 0x20, 0x3d, 0x20, 0x25, 0x2e, 0x38, 0x66, 0x0a], UnitRange{Int64}[1:17, 20:75, 81:91, 97:134, 139:139], (%.0d, %.10f, %.12f, %.8f), 4), %194, %195, %197, %198)
└───        goto #57
34 ─        (R_iter = Main.calculate_R_iter(R_iter, R2))
│    %202 = Base.getproperty(Printf, :format)::Core.Const(Printf.format)
│    %203 = Main.stdout::IO
│    %204 = i::Int64
│    %205 = r2norm::Float64
│    %206 = (e_old - e_new)::Float64
│    %207 = Main.abs(%206)::Float64
│           (%202)(%203, Printf.Format{Base.CodeUnits{UInt8, String}, Tuple{Printf.Spec{Val{'d'}}, Printf.Spec{Val{'f'}}, Printf.Spec{Val{'f'}}}}(UInt8[0x43, 0x43, 0x44, 0x20, 0x4e, 0x6f, 0x74, 0x20, 0x43, 0x6f, 0x6e, 0x76, 0x65, 0x72, 0x67, 0x65, 0x64, 0x20, 0x69, 0x6e, 0x20, 0x25, 0x64, 0x20, 0x69, 0x74, 0x65, 0x72, 0x61, 0x74, 0x69, 0x6f, 0x6e, 0x73, 0x2c, 0x20, 0x63, 0x75, 0x72, 0x72, 0x65, 0x6e, 0x74, 0x20, 0x7c, 0x7c, 0x52, 0x32, 0x7c, 0x7c, 0x20, 0x3d, 0x20, 0x25, 0x2e, 0x31, 0x30, 0x66, 0x20, 0x61, 0x6e, 0x64, 0x20, 0xce, 0x94, 0x45, 0x20, 0x3d, 0x20, 0x25, 0x2e, 0x31, 0x32, 0x66, 0x20, 0x0a], UnitRange{Int64}[1:21, 24:53, 59:69, 75:76], (%.0d, %.10f, %.12f), 3), %204, %205, %207)
│    %209 = (i >= p_min::Core.Const(2))::Bool
└───        goto #46 if not %209
35 ─ %211 = Base.string("DIIS is being implemented in iteration ", i)::String
│           Main.println(%211)
│    %213 = Main.length(R_iter_storage)::Int64
│    %214 = (%213 >= p_max::Core.Const(6))::Bool
└───        goto #37 if not %214
36 ─        Main.popfirst!(R_iter_storage)
│           Main.popfirst!(T2_storage)
└───        Main.popfirst!(R2_storage)
37 ┄ %219 = R_iter_storage::Vector{TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}}
│    %220 = Main.copy(R_iter)::TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}
│           Main.push!(%219, %220)
│    %222 = T2_storage::Vector{TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}}
│    %223 = Main.copy(T2_old)::TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}
│           Main.push!(%222, %223)
│    %225 = R2_storage::Vector{TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}}
│    %226 = Main.copy(R2)::TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}
│           Main.push!(%225, %226)
│    %228 = Main.length(R_iter_storage)::Int64
│           (@_37 = %228)
│    %230 = (@_37 isa Main.Int64)::Core.Const(true)
└───        goto #39 if not %230
38 ─        goto #40
39 ─        Core.Const(:(Base.convert(Main.Int64, @_37)))
└───        Core.Const(:(@_37 = Core.typeassert(%233, Main.Int64)))
40 ┄        (p = @_37)
│           (diis_result = Main.PerformDiis(R_iter_storage, p, T2_storage, R2_storage))
│    %237 = (diis_result == false)::Core.Const(false)
└───        goto #42 if not %237
41 ─        Core.Const(:(Main.nothing))
└───        Core.Const(:(goto %253))
42 ┄ %241 = diis_result::TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}
│    %242 = Main.TrivialTensorMap::Core.Const(TrivialTensorMap)
│    %243 = Main.CartesianSpace::Core.Const(CartesianSpace)
│    %244 = Core.apply_type(Main.Matrix, Main.Float64)::Core.Const(Matrix{Float64})
│    %245 = Core.apply_type(%242, %243, 2, 2, %244)::Core.Const(TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}})
│           (@_38 = %241)
│    %247 = (@_38 isa %245)::Core.Const(true)
└───        goto #44 if not %247
43 ─        goto #45
44 ─        Core.Const(:(Base.convert(%245, @_38)))
└───        Core.Const(:(@_38 = Core.typeassert(%250, %245)))
45 ┄        (T2_old = @_38)
└───        goto #49
46 ─ %254 = i::Int64
│    %255 = (p_min::Core.Const(2) - 1)::Core.Const(1)
│    %256 = (%254 == %255)::Bool
└───        goto #48 if not %256
47 ─        Main.println("DIIS will start in next iteration")
│    %259 = R_iter_storage::Vector{TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}}
│    %260 = Main.copy(R_iter)::TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}
│           Main.push!(%259, %260)
│    %262 = T2_storage::Vector{TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}}
│    %263 = Main.copy(T2_old)::TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}
│           Main.push!(%262, %263)
│    %265 = R2_storage::Vector{TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}}
│    %266 = Main.copy(R2)::TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}
│           Main.push!(%265, %266)
│    %268 = R_iter_storage::Vector{TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}}
│    %269 = Main.length(R_iter_storage)::Int64
│    %270 = T2_storage::Vector{TrivialTensorMap{CartesianSpace, 2, 2, Matrix{Float64}}}
│           Main.just_show_Bmatrix(%268, %269, %270, R2_storage)
└───        goto #49
48 ─        Main.println("DIIS start is still far away")
49 ┄ %274 = e_new::Float64
│           (@_39 = %274)
│    %276 = (@_39 isa Main.Float64)::Core.Const(true)
└───        goto #51 if not %276
50 ─        goto #52
51 ─        Core.Const(:(Base.convert(Main.Float64, @_39)))
└───        Core.Const(:(@_39 = Core.typeassert(%279, Main.Float64)))
52 ┄        (e_old = @_39)
│    %282 = Main.calculate_ECCD(T2_old)::Any
│           (@_40 = %282)
│    %284 = (@_40 isa Main.Float64)::Bool
└───        goto #54 if not %284
53 ─        goto #55
54 ─ %287 = Base.convert(Main.Float64, @_40)::Any
└───        (@_40 = Core.typeassert(%287, Main.Float64))
55 ┄        (e_new = @_40::Float64)
│    %290 = (e_new + erhf)::Float64
│    %291 = Base.string("Energy after updating = ", %290, "\n-------------------------\n")::String
│           Main.println(%291)
│           (@_3 = Base.iterate(%131, %143))
│    %294 = (@_3 === nothing)::Bool
│    %295 = Base.not_int(%294)::Bool
└───        goto #57 if not %295
56 ─        goto #23
57 ┄        return nothing
``````

Now I know that whenever any variable has the type `::Any` , in the output, that indicates type instability and must be eliminated. However, there seem to be variables name `@_27` and `@_28` which are type unstable, however these are not any variables that I have declared or used ? Are they some kind of intermediates ? How can I figure out where exactly these instabilitites are happening and fix them?

I would suggest, first, split your code into smaller functions. Separate the data input/generation, from the actual calculations. Data input and generation are many times type-unstable, and passing the final data to functions that actually require performance is the way to not propagate these instabilities.

Also, by splitting into smaller functions, it will be easier to interpret the output of `@code_warntype` for each function.

Also: note that `initialize_ccvariables()` does not have any argument and its result is not being assigned to any variable, thus suggesting that it is modifying some global variables. This may be the cause of instabilities (not sure, because we cannot see that code).

4 Likes

I tried to do the follows by first analyzing the type instabilities of one of the smaller functions related to input generation given below

``````function initialize_cc_variables()
main("water_dump.fcidump")
K = deserialize("K.jlbin")
f = deserialize("f.jlbin")
g = deserialize("g.jlbin")
nocc = deserialize("nocc.jlbin")
nv = deserialize("nv.jlbin")
g_vvoo::Array{Float64,4} = zeros(nv,nv,nocc,nocc)
g_voov::Array{Float64,4} = zeros(nv,nocc,nocc,nv)
g_vovo::Array{Float64,4} = zeros(nv,nocc,nv,nocc)
g_oovv::Array{Float64,4} = zeros(nocc,nocc,nv,nv)
g_oooo::Array{Float64,4} = zeros(nocc,nocc,nocc,nocc)
g_vvvv::Array{Float64,4} = zeros(nv,nv,nv,nv)
fvv = zeros(Float64,nv,nv)
foo = zeros(Float64,nocc,nocc)
for a in 1:nv, b in 1:nv
fvv[a,b] = f[nocc+a,nocc+b]

for i in 1:nocc, j in 1:nocc
g_vvoo[a,b,i,j] = g[nocc+a,i,nocc+b,j]
g_voov[a,i,j,b] = g[nocc+a,j,i,nocc+b]
g_vovo[a,i,b,j] = g[nocc+a,nocc+b,i,j]
g_oovv[i,j,a,b] = g[i,nocc+a,j,nocc+b]
foo[i,j] = f[i,j]
end
end
for i in 1:nocc, j in 1:nocc, k in 1:nocc, l in 1:nocc
g_oooo[i,j,k,l] = g[i,k,j,l]
end
for a in 1:nv, b in 1:nv, c in 1:nv, d in 1:nv
g_vvvv[a,b,c,d] = g[nocc+a,nocc+c,nocc+b,nocc+d]
end
serialize("g_vvoo.jlbin", g_vvoo)
serialize("g_voov.jlbin", g_voov)
serialize("g_vovo.jlbin", g_vovo)
serialize("g_oovv.jlbin", g_oovv)
serialize("g_oooo.jlbin", g_oooo)
serialize("g_vvvv.jlbin", g_vvvv)
serialize("K.jlbin", K)
serialize("fvv.jlbin", fvv)
serialize("foo.jlbin", foo)
end
``````

Here as you can see, I have explicitly declared the type of `g_vvoo` as `Array{Float64}` , and yet when I analyze the `@code_warntype` I am presented with

``````MethodInstance for initialize_cc_variables()
from initialize_cc_variables() @ Main ~/Desktop/Andreas Köhn/phase 1.5/code/tensorkit/ccd-helper.jl:180
Arguments
#self#::Core.Const(initialize_cc_variables)
Locals
@_2::Union{}
@_3::Union{}
@_4::Union{}
foo::Union{}
fvv::Union{}
g_vvvv::Union{}
g_oooo::Union{}
g_oovv::Union{}
g_vovo::Union{}
g_voov::Union{}
g_vvoo::Union{}
nv::Union{}
nocc::Union{}
g::Union{}
f::Union{}
K::Union{}
@_18::Union{}
a@_19::Union{}
@_20::Union{}
b@_21::Union{}
a@_22::Union{}
@_23::Union{}
i@_24::Union{}
j@_25::Union{}
i@_26::Union{}
@_27::Union{}
i@_28::Union{}
@_29::Union{}
j@_30::Union{}
@_31::Union{}
k@_32::Union{}
i@_33::Union{}
j@_34::Union{}
l::Union{}
k@_36::Union{}
@_37::Union{}
a@_38::Union{}
@_39::Union{}
b@_40::Union{}
@_41::Union{}
c@_42::Union{}
a@_43::Union{}
b@_44::Union{}
d::Union{}
c@_46::Union{}
@_47::Union{}
@_48::Union{}
@_49::Union{}
@_50::Union{}
@_51::Union{}
@_52::Union{}
Body::Union{}
1 ─     Core.NewvarNode(:(@_2))
│       Core.NewvarNode(:(@_3))
│       Core.NewvarNode(:(@_4))
│       Core.NewvarNode(:(foo))
│       Core.NewvarNode(:(fvv))
│       Core.NewvarNode(:(g_vvvv))
│       Core.NewvarNode(:(g_oooo))
│       Core.NewvarNode(:(g_oovv))
│       Core.NewvarNode(:(g_vovo))
│       Core.NewvarNode(:(g_voov))
│       Core.NewvarNode(:(g_vvoo))
│       Main.main("water_dump.fcidump")
│       Core.Const(:(K = Main.deserialize("K.jlbin")))
│       Core.Const(:(f = Main.deserialize("f.jlbin")))
│       Core.Const(:(g = Main.deserialize("g.jlbin")))
│       Core.Const(:(nocc = Main.deserialize("nocc.jlbin")))
│       Core.Const(:(nv = Main.deserialize("nv.jlbin")))
│       Core.Const(:(Main.zeros(nv, nv, nocc, nocc)))
│       Core.Const(:(Core.apply_type(Main.Array, Main.Float64, 4)))
│       Core.Const(:(@_47 = %18))
│       Core.Const(:(@_47 isa %19))
│       Core.Const(:(goto %24 if not %21))
│       Core.Const(:(goto %26))
│       Core.Const(:(Base.convert(%19, @_47)))
│       Core.Const(:(@_47 = Core.typeassert(%24, %19)))
│       Core.Const(:(g_vvoo = @_47))
│       Core.Const(:(Main.zeros(nv, nocc, nocc, nv)))
│       Core.Const(:(Core.apply_type(Main.Array, Main.Float64, 4)))
│       Core.Const(:(@_48 = %27))
│       Core.Const(:(@_48 isa %28))
│       Core.Const(:(goto %33 if not %30))
│       Core.Const(:(goto %35))
│       Core.Const(:(Base.convert(%28, @_48)))
│       Core.Const(:(@_48 = Core.typeassert(%33, %28)))
│       Core.Const(:(g_voov = @_48))
│       Core.Const(:(Main.zeros(nv, nocc, nv, nocc)))
│       Core.Const(:(Core.apply_type(Main.Array, Main.Float64, 4)))
│       Core.Const(:(@_49 = %36))
│       Core.Const(:(@_49 isa %37))
│       Core.Const(:(goto %42 if not %39))
│       Core.Const(:(goto %44))
│       Core.Const(:(Base.convert(%37, @_49)))
│       Core.Const(:(@_49 = Core.typeassert(%42, %37)))
│       Core.Const(:(g_vovo = @_49))
│       Core.Const(:(Main.zeros(nocc, nocc, nv, nv)))
│       Core.Const(:(Core.apply_type(Main.Array, Main.Float64, 4)))
│       Core.Const(:(@_50 = %45))
│       Core.Const(:(@_50 isa %46))
│       Core.Const(:(goto %51 if not %48))
│       Core.Const(:(goto %53))
│       Core.Const(:(Base.convert(%46, @_50)))
│       Core.Const(:(@_50 = Core.typeassert(%51, %46)))
│       Core.Const(:(g_oovv = @_50))
│       Core.Const(:(Main.zeros(nocc, nocc, nocc, nocc)))
│       Core.Const(:(Core.apply_type(Main.Array, Main.Float64, 4)))
│       Core.Const(:(@_51 = %54))
│       Core.Const(:(@_51 isa %55))
│       Core.Const(:(goto %60 if not %57))
│       Core.Const(:(goto %62))
│       Core.Const(:(Base.convert(%55, @_51)))
│       Core.Const(:(@_51 = Core.typeassert(%60, %55)))
│       Core.Const(:(g_oooo = @_51))
│       Core.Const(:(Main.zeros(nv, nv, nv, nv)))
│       Core.Const(:(Core.apply_type(Main.Array, Main.Float64, 4)))
│       Core.Const(:(@_52 = %63))
│       Core.Const(:(@_52 isa %64))
│       Core.Const(:(goto %69 if not %66))
│       Core.Const(:(goto %71))
│       Core.Const(:(Base.convert(%64, @_52)))
│       Core.Const(:(@_52 = Core.typeassert(%69, %64)))
│       Core.Const(:(g_vvvv = @_52))
│       Core.Const(:(fvv = Main.zeros(Main.Float64, nv, nv)))
│       Core.Const(:(foo = Main.zeros(Main.Float64, nocc, nocc)))
│       Core.Const(:(1:nv))
│       Core.Const(:(@_4 = Base.iterate(%74)))
│       Core.Const(:(@_4 === nothing))
│       Core.Const(:(Base.not_int(%76)))
│       Core.Const(:(goto %163 if not %77))
│       Core.Const(:(@_4))
│       Core.Const(:(a@_19 = Core.getfield(%79, 1)))
│       Core.Const(:(Core.getfield(%79, 2)))
│       Core.Const(:(1:nv))
│       Core.Const(:(@_18 = Base.iterate(%82)))
│       Core.Const(:(@_18 === nothing))
│       Core.Const(:(Base.not_int(%84)))
│       Core.Const(:(goto %158 if not %85))
│       Core.Const(:(a@_19))
│       Core.Const(:(a@_22 = %87))
│       Core.Const(:(@_18))
│       Core.Const(:(b@_21 = Core.getfield(%89, 1)))
│       Core.Const(:(Core.getfield(%89, 2)))
│       Core.Const(:(f))
│       Core.Const(:(nocc + a@_22))
│       Core.Const(:(nocc + b@_21))
│       Core.Const(:(Base.getindex(%92, %93, %94)))
│       Core.Const(:(Base.setindex!(fvv, %95, a@_22, b@_21)))
│       Core.Const(:(1:nocc))
│       Core.Const(:(@_20 = Base.iterate(%97)))
│       Core.Const(:(@_20 === nothing))
│       Core.Const(:(Base.not_int(%99)))
│       Core.Const(:(goto %153 if not %100))
│       Core.Const(:(@_20))
│       Core.Const(:(i@_24 = Core.getfield(%102, 1)))
│       Core.Const(:(Core.getfield(%102, 2)))
│       Core.Const(:(1:nocc))
│       Core.Const(:(@_23 = Base.iterate(%105)))
│       Core.Const(:(@_23 === nothing))
│       Core.Const(:(Base.not_int(%107)))
│       Core.Const(:(goto %148 if not %108))
│       Core.Const(:(i@_24))
│       Core.Const(:(i@_26 = %110))
│       Core.Const(:(@_23))
│       Core.Const(:(j@_25 = Core.getfield(%112, 1)))
│       Core.Const(:(Core.getfield(%112, 2)))
│       Core.Const(:(g))
│       Core.Const(:(nocc + a@_22))
│       Core.Const(:(i@_26))
│       Core.Const(:(nocc + b@_21))
│       Core.Const(:(Base.getindex(%115, %116, %117, %118, j@_25)))
│       Core.Const(:(Base.setindex!(g_vvoo, %119, a@_22, b@_21, i@_26, j@_25)))
│       Core.Const(:(g))
│       Core.Const(:(nocc + a@_22))
│       Core.Const(:(j@_25))
│       Core.Const(:(i@_26))
│       Core.Const(:(nocc + b@_21))
│       Core.Const(:(Base.getindex(%121, %122, %123, %124, %125)))
│       Core.Const(:(Base.setindex!(g_voov, %126, a@_22, i@_26, j@_25, b@_21)))
│       Core.Const(:(g))
│       Core.Const(:(nocc + a@_22))
│       Core.Const(:(nocc + b@_21))
│       Core.Const(:(i@_26))
│       Core.Const(:(Base.getindex(%128, %129, %130, %131, j@_25)))
│       Core.Const(:(Base.setindex!(g_vovo, %132, a@_22, i@_26, b@_21, j@_25)))
│       Core.Const(:(g))
│       Core.Const(:(i@_26))
│       Core.Const(:(nocc + a@_22))
│       Core.Const(:(j@_25))
│       Core.Const(:(nocc + b@_21))
│       Core.Const(:(Base.getindex(%134, %135, %136, %137, %138)))
│       Core.Const(:(Base.setindex!(g_oovv, %139, i@_26, j@_25, a@_22, b@_21)))
│       Core.Const(:(Base.getindex(f, i@_26, j@_25)))
│       Core.Const(:(Base.setindex!(foo, %141, i@_26, j@_25)))
│       Core.Const(:(@_23 = Base.iterate(%105, %114)))
│       Core.Const(:(@_23 === nothing))
│       Core.Const(:(Base.not_int(%144)))
│       Core.Const(:(goto %148 if not %145))
│       Core.Const(:(goto %110))
│       Core.Const(:(@_20 = Base.iterate(%97, %104)))
│       Core.Const(:(@_20 === nothing))
│       Core.Const(:(Base.not_int(%149)))
│       Core.Const(:(goto %153 if not %150))
│       Core.Const(:(goto %102))
│       Core.Const(:(@_18 = Base.iterate(%82, %91)))
│       Core.Const(:(@_18 === nothing))
│       Core.Const(:(Base.not_int(%154)))
│       Core.Const(:(goto %158 if not %155))
│       Core.Const(:(goto %87))
│       Core.Const(:(@_4 = Base.iterate(%74, %81)))
│       Core.Const(:(@_4 === nothing))
│       Core.Const(:(Base.not_int(%159)))
│       Core.Const(:(goto %163 if not %160))
│       Core.Const(:(goto %79))
│       Core.Const(:(1:nocc))
│       Core.Const(:(@_3 = Base.iterate(%163)))
│       Core.Const(:(@_3 === nothing))
│       Core.Const(:(Base.not_int(%165)))
│       Core.Const(:(goto %223 if not %166))
│       Core.Const(:(@_3))
│       Core.Const(:(i@_28 = Core.getfield(%168, 1)))
│       Core.Const(:(Core.getfield(%168, 2)))
│       Core.Const(:(1:nocc))
│       Core.Const(:(@_27 = Base.iterate(%171)))
│       Core.Const(:(@_27 === nothing))
│       Core.Const(:(Base.not_int(%173)))
│       Core.Const(:(goto %218 if not %174))
│       Core.Const(:(@_27))
│       Core.Const(:(j@_30 = Core.getfield(%176, 1)))
│       Core.Const(:(Core.getfield(%176, 2)))
│       Core.Const(:(1:nocc))
│       Core.Const(:(@_29 = Base.iterate(%179)))
│       Core.Const(:(@_29 === nothing))
│       Core.Const(:(Base.not_int(%181)))
│       Core.Const(:(goto %213 if not %182))
│       Core.Const(:(@_29))
│       Core.Const(:(k@_32 = Core.getfield(%184, 1)))
│       Core.Const(:(Core.getfield(%184, 2)))
│       Core.Const(:(1:nocc))
│       Core.Const(:(@_31 = Base.iterate(%187)))
│       Core.Const(:(@_31 === nothing))
│       Core.Const(:(Base.not_int(%189)))
│       Core.Const(:(goto %208 if not %190))
│       Core.Const(:(i@_28))
│       Core.Const(:(i@_33 = %192))
│       Core.Const(:(j@_30))
│       Core.Const(:(j@_34 = %194))
│       Core.Const(:(k@_32))
│       Core.Const(:(k@_36 = %196))
│       Core.Const(:(@_31))
│       Core.Const(:(l = Core.getfield(%198, 1)))
│       Core.Const(:(Core.getfield(%198, 2)))
│       Core.Const(:(Base.getindex(g, i@_33, k@_36, j@_34, l)))
│       Core.Const(:(Base.setindex!(g_oooo, %201, i@_33, j@_34, k@_36, l)))
│       Core.Const(:(@_31 = Base.iterate(%187, %200)))
│       Core.Const(:(@_31 === nothing))
│       Core.Const(:(Base.not_int(%204)))
│       Core.Const(:(goto %208 if not %205))
│       Core.Const(:(goto %192))
│       Core.Const(:(@_29 = Base.iterate(%179, %186)))
│       Core.Const(:(@_29 === nothing))
│       Core.Const(:(Base.not_int(%209)))
│       Core.Const(:(goto %213 if not %210))
│       Core.Const(:(goto %184))
│       Core.Const(:(@_27 = Base.iterate(%171, %178)))
│       Core.Const(:(@_27 === nothing))
│       Core.Const(:(Base.not_int(%214)))
│       Core.Const(:(goto %218 if not %215))
│       Core.Const(:(goto %176))
│       Core.Const(:(@_3 = Base.iterate(%163, %170)))
│       Core.Const(:(@_3 === nothing))
│       Core.Const(:(Base.not_int(%219)))
│       Core.Const(:(goto %223 if not %220))
│       Core.Const(:(goto %168))
│       Core.Const(:(1:nv))
│       Core.Const(:(@_2 = Base.iterate(%223)))
│       Core.Const(:(@_2 === nothing))
│       Core.Const(:(Base.not_int(%225)))
│       Core.Const(:(goto %288 if not %226))
│       Core.Const(:(@_2))
│       Core.Const(:(a@_38 = Core.getfield(%228, 1)))
│       Core.Const(:(Core.getfield(%228, 2)))
│       Core.Const(:(1:nv))
│       Core.Const(:(@_37 = Base.iterate(%231)))
│       Core.Const(:(@_37 === nothing))
│       Core.Const(:(Base.not_int(%233)))
│       Core.Const(:(goto %283 if not %234))
│       Core.Const(:(@_37))
│       Core.Const(:(b@_40 = Core.getfield(%236, 1)))
│       Core.Const(:(Core.getfield(%236, 2)))
│       Core.Const(:(1:nv))
│       Core.Const(:(@_39 = Base.iterate(%239)))
│       Core.Const(:(@_39 === nothing))
│       Core.Const(:(Base.not_int(%241)))
│       Core.Const(:(goto %278 if not %242))
│       Core.Const(:(@_39))
│       Core.Const(:(c@_42 = Core.getfield(%244, 1)))
│       Core.Const(:(Core.getfield(%244, 2)))
│       Core.Const(:(1:nv))
│       Core.Const(:(@_41 = Base.iterate(%247)))
│       Core.Const(:(@_41 === nothing))
│       Core.Const(:(Base.not_int(%249)))
│       Core.Const(:(goto %273 if not %250))
│       Core.Const(:(a@_38))
│       Core.Const(:(a@_43 = %252))
│       Core.Const(:(b@_40))
│       Core.Const(:(b@_44 = %254))
│       Core.Const(:(c@_42))
│       Core.Const(:(c@_46 = %256))
│       Core.Const(:(@_41))
│       Core.Const(:(d = Core.getfield(%258, 1)))
│       Core.Const(:(Core.getfield(%258, 2)))
│       Core.Const(:(g))
│       Core.Const(:(nocc + a@_43))
│       Core.Const(:(nocc + c@_46))
│       Core.Const(:(nocc + b@_44))
│       Core.Const(:(nocc + d))
│       Core.Const(:(Base.getindex(%261, %262, %263, %264, %265)))
│       Core.Const(:(Base.setindex!(g_vvvv, %266, a@_43, b@_44, c@_46, d)))
│       Core.Const(:(@_41 = Base.iterate(%247, %260)))
│       Core.Const(:(@_41 === nothing))
│       Core.Const(:(Base.not_int(%269)))
│       Core.Const(:(goto %273 if not %270))
│       Core.Const(:(goto %252))
│       Core.Const(:(@_39 = Base.iterate(%239, %246)))
│       Core.Const(:(@_39 === nothing))
│       Core.Const(:(Base.not_int(%274)))
│       Core.Const(:(goto %278 if not %275))
│       Core.Const(:(goto %244))
│       Core.Const(:(@_37 = Base.iterate(%231, %238)))
│       Core.Const(:(@_37 === nothing))
│       Core.Const(:(Base.not_int(%279)))
│       Core.Const(:(goto %283 if not %280))
│       Core.Const(:(goto %236))
│       Core.Const(:(@_2 = Base.iterate(%223, %230)))
│       Core.Const(:(@_2 === nothing))
│       Core.Const(:(Base.not_int(%284)))
│       Core.Const(:(goto %288 if not %285))
│       Core.Const(:(goto %228))
│       Core.Const(:(Main.serialize("g_vvoo.jlbin", g_vvoo)))
│       Core.Const(:(Main.serialize("g_voov.jlbin", g_voov)))
│       Core.Const(:(Main.serialize("g_vovo.jlbin", g_vovo)))
│       Core.Const(:(Main.serialize("g_oovv.jlbin", g_oovv)))
│       Core.Const(:(Main.serialize("g_oooo.jlbin", g_oooo)))
│       Core.Const(:(Main.serialize("g_vvvv.jlbin", g_vvvv)))
│       Core.Const(:(Main.serialize("K.jlbin", K)))
│       Core.Const(:(Main.serialize("fvv.jlbin", fvv)))
│       Core.Const(:(Main.serialize("foo.jlbin", foo)))
└──     Core.Const(:(return %296))
``````

Where the same variable is given as `Union{}` data type. Is this normal ? Or am I missing something ?

What is this function call? What it does?

I don’t know that `serialize` returns as a type, but it is likely type-unstable. But this is not a big issue if you can contain these instabilities in the function that initializes the data.

1 Like

The purpose of this function (given below) :

Function Definition
``````function main(pathtofcidump)
linenum=0
norb=0
nelec=0
begin # Finds out from where the main integral data starts
file = open(pathtofcidump, "r")
for line in eachline(file)
linenum += 1
if strip(line) == "/"
# println("The line number where integral data starts is: ", linenum)
break
end
end
end

begin
# Open the file
file = open(pathtofcidump, "r")

# Loop over each line in the file
for line in eachline(file)
# Check if the line starts with "&FCI NORB="
if startswith(line, " &FCI NORB=")
line = split(line, ",")
words1 = split(line[1], " ")
words2 = split(line[2], " ")
norb = parse(Int, words1[end])
nelec = parse(Int, words2[end])
break
end
end

# Close the file
close(file)
end
data = readdlm(pathtofcidump, skipstart=linenum)
hnuc::Float64 = data[end,1]

data = copy(data[1:end-1, :])
# h = Array{Union{Missing, Float64}}(missing, norb, norb);
# g = Array{Union{Missing, Float64}}(missing, norb, norb, norb, norb);
h::Array{Float64,2} = fill(0.0, norb, norb)
g::Array{Float64,4} = fill(0.0, norb, norb, norb, norb)
l::Int64 = length(data[:,1])
non_redundant_indices = []
for i in 1:l
if(data[i,4]== 0 && data[i,5]==0)
I = round(Int,data[i,2])
J = round(Int,data[i,3])
h[I,J] = data[i,1]
else
I = round(Int,data[i,2])
J = round(Int,data[i,3])
K = round(Int,data[i,4])
L = round(Int,data[i,5])
push!(non_redundant_indices, [I, J, K, L])
g[I,J,K,L] = data[i,1]
end
end
for (I,J,K,L) in non_redundant_indices
open("non_redundant_indices.txt", "a") do f
println(f, I,J,K,L,"   ",g[I,J,K,L])
end
g[K,L,I,J]=g[I,J,K,L]
g[J,I,L,K]=g[I,J,K,L]
g[L,K,J,I]=g[I,J,K,L]
g[J,I,K,L]=g[I,J,K,L]
g[L,K,I,J]=g[I,J,K,L]
g[I,J,L,K]=g[I,J,K,L]
g[K,L,J,I]=g[I,J,K,L]

end
for I in 1:norb
for J in 1:I
h[J,I] = h[I,J]
end
end
nocc = round(Int,nelec/2)
nv = norb - nocc
h = convert(Array{Float64}, h)
g = convert(Array{Float64}, g)
K = zeros(nv,nv,nocc,nocc)
for a in 1:nv, b in 1:nv, i in 1:nocc, j in 1:nocc
K[a,b,i,j] = g[nocc+a,i,nocc+b,j]
end
H1::Float64 = 0.0;
H2::Float64 = 0.0;
for i in 1:nocc
H1 = H1+  h[i,i]
for j in 1:nocc
H2 = H2 +  2 * g[i,i,j,j] - g[i,j,i,j]
end
end
erhf = 2*H1 + H2 + hnuc
f = zeros(size(h)) # fock matrix Initialization
for p in 1:norb , q in 1:norb
s_pq=0
for i in 1:nocc
s_pq += 2*g[p,q,i,i]-g[p,i,i,q]
end
f[p,q] = h[p,q] + s_pq
end
serialize("h.jlbin", h)
serialize("g.jlbin", g)
serialize("hnuc.jlbin", hnuc)
serialize("norbs.jlbin", norb)
serialize("nelec.jlbin", nelec)
serialize("nocc.jlbin", nocc)
serialize("nv.jlbin", nv)
serialize("K.jlbin", K)
serialize("erhf.jlbin", erhf)
serialize("f.jlbin", f)
end
``````

Is to simply parse a text file and initialize some large tensors using the data read from the file. It then writes the data to disc using `serialize(...)`. Also, could you please explain what is meant by “containing the instabilities within the function”?

https://docs.julialang.org/en/v1/manual/performance-tips/#kernel-functions

2 Likes

Yes these seem to be intermediaries from the IO.

That’s why separating initialization from computation into different functions is usually a good idea. However you wrote explicit type annotation, so I don’t think you suffer performance loss due to type instability. At least not in the code you have shown here. Note that you would not need these type annotations if you put the initialization in a different function.

Also this back and forth serialization looks a bit odd to me. Why not keep the values in memory? E.g. create a struct that holds the inputs to the algorithm (and serialize that if you need), then a function that is concerned with initializing all the values and creating the struct and then a third function where the iterative computation happens. Then it would be also much easier to benchmark the computation part in order to diagnose bottlenecks and optimise.

That is precisely what I am trying to avoid. These tensors `g_xxxx` will become very large and we want to read these files from the disk only when the part where they would be needed comes up. I know the code attached may not be very clear on this, but my other post talks about dealing with this issue.

But these tensors don’t appear in the loop of the code you shared? Only in this initialization function. Or are some function within this iterative process reading and writing files again? Then I can see why things are slow…

Since this whole post is essentially an “please help optimize the runtime of code”-type of post and it appears that the code you shared is not the whole story, could you share a bit more about your problem (like typical matrix sizes etc.), the algorithm and your current approach. To optimize for performance it is important to know the boundary conditions as optimization typically comes with trade-offs.

I don’t think that the code you shared suffers from type-instabilities. Not because there are none, but because you already contained them by explicit type annotations. But I somewhat doubt that even without the annotations this code would much slower, since I think the true bottleneck is elsewhere. Please profile your code and try to identify which part is slowest and then ideally try to reduce it to some minimal code that showcases the problem and is ideally executable for us, so we can play with it

1 Like

Sure, attached
ccd.jl (3.7 KB)
ccd-helper.jl (23.8 KB)

is the current code after making some minor modifications as suggested in this and the other thread. The first file is the main executable with the second having all the important function definitions. Please note that many of the functions are not needed and have been left there for legacy purposes.

And below is the:

@code_warntype output
``````MethodInstance for ccd_by_hand(::Int64)
from ccd_by_hand(maxitr) @ Main ~/Desktop/Andreas Köhn/phase 1.5/code/tensoroperations/ccd.jl:3
Arguments
#self#::Core.Const(ccd_by_hand)
maxitr::Int64
Locals
@_3::Union{Nothing, Tuple{Int64, Int64}}
R2_storage::Vector{Array{Float64, 4}}
T2_storage::Vector{Array{Float64, 4}}
R_iter_storage::Vector{Array{Float64, 4}}
p_min::Int64
p_max::Int64
etol::Float64
earr::Vector{Any}
e_old::Float64
e_new::Float64
normtol::Float64
T2_old::Array{Float64, 4}
Scaled_R2::Array{Float64, 4}
R_iter::Array{Float64, 4}
R2::Array{Float64, 4}
R2u::Array{Float64, 4}
erhf::Float64
nocc::Int64
nv::Int64
@_22::Int64
i::Int64
diis_result::Array{Float64, 4}
p::Int64
r2norm::Float64
conv::Bool
@_28::Any
@_29::Any
@_30::Any
Body::Nothing
1 ──        Core.NewvarNode(:(@_3))
│           Core.NewvarNode(:(R2_storage))
│           Core.NewvarNode(:(T2_storage))
│           Core.NewvarNode(:(R_iter_storage))
│           Core.NewvarNode(:(p_min))
│           Core.NewvarNode(:(p_max))
│           Core.NewvarNode(:(etol))
│           Core.NewvarNode(:(earr))
│           Core.NewvarNode(:(e_old))
│           Core.NewvarNode(:(e_new))
│           Core.NewvarNode(:(normtol))
│           Core.NewvarNode(:(T2_old))
│           Core.NewvarNode(:(Scaled_R2))
│           Core.NewvarNode(:(R_iter))
│           Core.NewvarNode(:(R2))
│           Core.NewvarNode(:(R2u))
│           Core.NewvarNode(:(erhf))
│           Core.NewvarNode(:(nocc))
│           Core.NewvarNode(:(nv))
│           Main.initialize_cc_variables()
│    %21  = Main.deserialize("nv.jlbin")::Any
│           (@_28 = %21)
│    %23  = (@_28 isa Main.Int64)::Bool
└───        goto #3 if not %23
2 ──        goto #4
3 ── %26  = Base.convert(Main.Int64, @_28)::Any
└───        (@_28 = Core.typeassert(%26, Main.Int64))
4 ┄─        (nv = @_28::Int64)
│    %29  = Main.deserialize("nocc.jlbin")::Any
│           (@_29 = %29)
│    %31  = (@_29 isa Main.Int64)::Bool
└───        goto #6 if not %31
5 ──        goto #7
6 ── %34  = Base.convert(Main.Int64, @_29)::Any
└───        (@_29 = Core.typeassert(%34, Main.Int64))
7 ┄─        (nocc = @_29::Int64)
│    %37  = Main.deserialize("erhf.jlbin")::Any
│           (@_30 = %37)
│    %39  = (@_30 isa Main.Float64)::Bool
└───        goto #9 if not %39
8 ──        goto #10
9 ── %42  = Base.convert(Main.Float64, @_30)::Any
└───        (@_30 = Core.typeassert(%42, Main.Float64))
10 ┄        (erhf = @_30::Float64)
│    %45  = Main.zeros(Main.Float64, nv, nv, nocc, nocc)::Array{Float64, 4}
│    %46  = Main.zeros(Main.Float64, nv, nv, nocc, nocc)::Array{Float64, 4}
│           (R2u = %45)
│           (R2 = %46)
│    %49  = Main.Float64::Core.Const(Float64)
│    %50  = Main.size(R2)::NTuple{4, Int64}
│           (R_iter = Main.zeros(%49, %50))
│    %52  = Main.Float64::Core.Const(Float64)
│    %53  = Main.size(R2)::NTuple{4, Int64}
│           (Scaled_R2 = Main.zeros(%52, %53))
│           (T2_old = Main.initialize_t2_only())
│           (normtol = 1.0e-8)
│           (e_new = Main.calculate_ECCD(T2_old))
│           (e_old = Main.copy(0.0))
│           (earr = Base.vect())
│           (etol = 1.0e-10)
│           (p_max = 6)
│           (p_min = 2)
│    %63  = Main.Array::Core.Const(Array)
│    %64  = Main.typeof(R2)::Core.Const(Array{Float64, 4})
│    %65  = Core.apply_type(%63, %64)::Core.Const(Array{Array{Float64, 4}})
│    %66  = Main.undef::Core.Const(UndefInitializer())
│           (R_iter_storage = (%65)(%66, 0))
│    %68  = Main.Array::Core.Const(Array)
│    %69  = Main.typeof(T2_old)::Core.Const(Array{Float64, 4})
│    %70  = Core.apply_type(%68, %69)::Core.Const(Array{Array{Float64, 4}})
│    %71  = Main.undef::Core.Const(UndefInitializer())
│           (T2_storage = (%70)(%71, 0))
│    %73  = Main.Array::Core.Const(Array)
│    %74  = Main.typeof(R2)::Core.Const(Array{Float64, 4})
│    %75  = Core.apply_type(%73, %74)::Core.Const(Array{Array{Float64, 4}})
│    %76  = Main.undef::Core.Const(UndefInitializer())
│           (R2_storage = (%75)(%76, 0))
│    %78  = earr::Vector{Any}
│    %79  = Base.getindex(e_new, 1)::Float64
│    %80  = (%79 + erhf)::Float64
│           Main.push!(%78, %80)
│    %82  = Base.string("Starting CCD Iterations with Convergence Criteria as ||R₂||<", normtol::Core.Const(1.0e-8), " , ΔE<", etol::Core.Const(1.0e-10), " and max iterations=", maxitr)::String
│           Main.println(%82)
│           Main.println("-----------------------------------------------------------------------------------------------")
│    %85  = (1:maxitr)::Core.PartialStruct(UnitRange{Int64}, Any[Core.Const(1), Int64])
│           (@_3 = Base.iterate(%85))
│    %87  = (@_3 === nothing)::Bool
│    %88  = Base.not_int(%87)::Bool
└───        goto #24 if not %88
11 ┄        Core.NewvarNode(:(diis_result))
│           Core.NewvarNode(:(p))
│    %92  = @_3::Tuple{Int64, Int64}
│           (i = Core.getfield(%92, 1))
│    %94  = Core.getfield(%92, 2)::Int64
│    %95  = Base.getproperty(Printf, :format)::Core.Const(Printf.format)
│    %96  = Main.stdout::IO
│    %97  = i::Int64
│    %98  = Base.getindex(e_new, 1)::Float64
│    %99  = (%98 + erhf)::Float64
│           (%95)(%96, Printf.Format{Base.CodeUnits{UInt8, String}, Tuple{Printf.Spec{Val{'d'}}, Printf.Spec{Val{'f'}}}}(UInt8[0x0a, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x0a, 0x53, 0x74, 0x61, 0x72, 0x74, 0x69, 0x6e, 0x67, 0x20, 0x49, 0x74, 0x65, 0x72, 0x61, 0x74, 0x69, 0x6f, 0x6e, 0x3a, 0x20, 0x25, 0x64, 0x20, 0x20, 0x43, 0x75, 0x72, 0x72, 0x65, 0x6e, 0x74, 0x20, 0x54, 0x6f, 0x74, 0x61, 0x6c, 0x20, 0x45, 0x6e, 0x65, 0x72, 0x67, 0x79, 0x3a, 0x20, 0x25, 0x2e, 0x38, 0x66, 0x0a], UnitRange{Int64}[1:47, 50:73, 78:78], (%.0d, %.8f), 2), %97, %99)
│    %101 = Base.getindex(e_old, 1)::Float64
│    %102 = (%101 + erhf)::Float64
│    %103 = Base.string("Energy before updating = ", %102)::String
│           Main.println(%103)
│           (R2 = Main.calculate_residual(T2_old, R2u, R2))
│           (T2_old = Main.update_T2(T2_old, R2, Scaled_R2))
│    %107 = Main.check_convergence(R2, normtol::Core.Const(1.0e-8), e_old, e_new, etol::Core.Const(1.0e-10))::Tuple{Bool, Float64}
│    %108 = Base.indexed_iterate(%107, 1)::Core.PartialStruct(Tuple{Bool, Int64}, Any[Bool, Core.Const(2)])
│           (conv = Core.getfield(%108, 1))
│           (@_22 = Core.getfield(%108, 2))
│    %111 = Base.indexed_iterate(%107, 2, @_22::Core.Const(2))::Core.PartialStruct(Tuple{Float64, Int64}, Any[Float64, Core.Const(3)])
│           (r2norm = Core.getfield(%111, 1))
└───        goto #13 if not conv
12 ─ %114 = Base.getproperty(Printf, :format)::Core.Const(Printf.format)
│    %115 = Main.stdout::IO
│    %116 = i::Int64
│    %117 = r2norm::Float64
│    %118 = Base.getindex(e_old, 1)::Float64
│    %119 = Base.getindex(e_new, 1)::Float64
│    %120 = (%118 - %119)::Float64
│    %121 = Main.abs(%120)::Float64
│    %122 = Base.getindex(e_new, 1)::Float64
│    %123 = (%122 + erhf)::Float64
│           (%114)(%115, Printf.Format{Base.CodeUnits{UInt8, String}, Tuple{Printf.Spec{Val{'d'}}, Printf.Spec{Val{'f'}}, Printf.Spec{Val{'f'}}, Printf.Spec{Val{'f'}}}}(UInt8[0x43, 0x43, 0x44, 0x20, 0x43, 0x6f, 0x6e, 0x76, 0x65, 0x72, 0x67, 0x65, 0x64, 0x20, 0x69, 0x6e, 0x20, 0x25, 0x64, 0x20, 0x69, 0x74, 0x65, 0x72, 0x61, 0x74, 0x69, 0x6f, 0x6e, 0x73, 0x2c, 0x20, 0x6e, 0x6f, 0x20, 0x66, 0x75, 0x72, 0x74, 0x68, 0x65, 0x72, 0x20, 0x75, 0x70, 0x64, 0x61, 0x74, 0x65, 0x73, 0x20, 0x72, 0x65, 0x71, 0x75, 0x69, 0x72, 0x65, 0x64, 0x2c, 0x63, 0x75, 0x72, 0x72, 0x65, 0x6e, 0x74, 0x20, 0x7c, 0x7c, 0x52, 0x32, 0x7c, 0x7c, 0x3d, 0x25, 0x2e, 0x31, 0x30, 0x66, 0x20, 0x61, 0x6e, 0x64, 0x20, 0xce, 0x94, 0x45, 0x20, 0x3d, 0x20, 0x25, 0x2e, 0x31, 0x32, 0x66, 0x0a, 0x0a, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x2d, 0x0a, 0x45, 0x5f, 0x54, 0x6f, 0x74, 0x61, 0x6c, 0x20, 0x3d, 0x20, 0x25, 0x2e, 0x38, 0x66, 0x0a], UnitRange{Int64}[1:17, 20:75, 81:91, 97:134, 139:139], (%.0d, %.10f, %.12f, %.8f), 4), %116, %117, %121, %123)
└───        goto #24
13 ─        (R_iter = Main.calculate_R_iter(R_iter, R2))
│    %127 = Base.getproperty(Printf, :format)::Core.Const(Printf.format)
│    %128 = Main.stdout::IO
│    %129 = i::Int64
│    %130 = r2norm::Float64
│    %131 = Base.getindex(e_old, 1)::Float64
│    %132 = Base.getindex(e_new, 1)::Float64
│    %133 = (%131 - %132)::Float64
│    %134 = Main.abs(%133)::Float64
│           (%127)(%128, Printf.Format{Base.CodeUnits{UInt8, String}, Tuple{Printf.Spec{Val{'d'}}, Printf.Spec{Val{'f'}}, Printf.Spec{Val{'f'}}}}(UInt8[0x43, 0x43, 0x44, 0x20, 0x4e, 0x6f, 0x74, 0x20, 0x43, 0x6f, 0x6e, 0x76, 0x65, 0x72, 0x67, 0x65, 0x64, 0x20, 0x69, 0x6e, 0x20, 0x25, 0x64, 0x20, 0x69, 0x74, 0x65, 0x72, 0x61, 0x74, 0x69, 0x6f, 0x6e, 0x73, 0x2c, 0x20, 0x63, 0x75, 0x72, 0x72, 0x65, 0x6e, 0x74, 0x20, 0x7c, 0x7c, 0x52, 0x32, 0x7c, 0x7c, 0x20, 0x3d, 0x20, 0x25, 0x2e, 0x31, 0x30, 0x66, 0x20, 0x61, 0x6e, 0x64, 0x20, 0xce, 0x94, 0x45, 0x20, 0x3d, 0x20, 0x25, 0x2e, 0x31, 0x32, 0x66, 0x20, 0x0a], UnitRange{Int64}[1:21, 24:53, 59:69, 75:76], (%.0d, %.10f, %.12f), 3), %129, %130, %134)
│    %136 = (i >= p_min::Core.Const(2))::Bool
└───        goto #19 if not %136
14 ─ %138 = Base.string("DIIS is being implemented in iteration ", i)::String
│           Main.println(%138)
│    %140 = Main.length(R_iter_storage)::Int64
│    %141 = (%140 >= p_max::Core.Const(6))::Bool
└───        goto #16 if not %141
15 ─        Main.popfirst!(R_iter_storage)
│           Main.popfirst!(T2_storage)
└───        Main.popfirst!(R2_storage)
16 ┄ %146 = R_iter_storage::Vector{Array{Float64, 4}}
│    %147 = Main.copy(R_iter)::Array{Float64, 4}
│           Main.push!(%146, %147)
│    %149 = T2_storage::Vector{Array{Float64, 4}}
│    %150 = Main.copy(T2_old)::Array{Float64, 4}
│           Main.push!(%149, %150)
│    %152 = R2_storage::Vector{Array{Float64, 4}}
│    %153 = Main.copy(R2)::Array{Float64, 4}
│           Main.push!(%152, %153)
│           (p = Main.length(R_iter_storage))
│           (diis_result = Main.PerformDiis(R_iter_storage, p, T2_storage, R2_storage))
│    %157 = (diis_result == false)::Core.Const(false)
└───        goto #18 if not %157
17 ─        Core.Const(:(Main.nothing))
└───        Core.Const(:(goto %162))
18 ┄        (T2_old = diis_result)
└───        goto #22
19 ─ %163 = i::Int64
│    %164 = (p_min::Core.Const(2) - 1)::Core.Const(1)
│    %165 = (%163 == %164)::Bool
└───        goto #21 if not %165
20 ─        Main.println("DIIS will start in next iteration")
│    %168 = R_iter_storage::Vector{Array{Float64, 4}}
│    %169 = Main.copy(R_iter)::Array{Float64, 4}
│           Main.push!(%168, %169)
│    %171 = T2_storage::Vector{Array{Float64, 4}}
│    %172 = Main.copy(T2_old)::Array{Float64, 4}
│           Main.push!(%171, %172)
│    %174 = R2_storage::Vector{Array{Float64, 4}}
│    %175 = Main.copy(R2)::Array{Float64, 4}
│           Main.push!(%174, %175)
│    %177 = R_iter_storage::Vector{Array{Float64, 4}}
│    %178 = Main.length(R_iter_storage)::Int64
│    %179 = T2_storage::Vector{Array{Float64, 4}}
│           Main.just_show_Bmatrix(%177, %178, %179, R2_storage)
└───        goto #22
21 ─        Main.println("DIIS start is still far away")
22 ┄        (e_old = e_new)
│           (e_new = Main.calculate_ECCD(T2_old))
│    %185 = Base.getindex(e_new, 1)::Float64
│    %186 = (%185 + erhf)::Float64
│    %187 = Base.string("Energy after updating = ", %186, "\n-------------------------\n")::String
│           Main.println(%187)
│    %189 = earr::Vector{Any}
│    %190 = Base.getindex(e_new, 1)::Float64
│    %191 = (%190 + erhf)::Float64
│           Main.push!(%189, %191)
│           (@_3 = Base.iterate(%85, %94))
│    %194 = (@_3 === nothing)::Bool
│    %195 = Base.not_int(%194)::Bool
└───        goto #24 if not %195
23 ─        goto #11
24 ┄        return nothing
``````

Lowered code will insert local variables to represent intermediate results in expanded code. The type annotations on left hand variables expand to additional `convert` and `typeassert` steps during any assignment, that’ll explain the existence of `@_` 27-29, 31, 36, 40, the intermediates annotated with `::Any`.

As pointed out already, 27-29 are the result of `deserialize`, where your input is a string and you get an arbitrary object; that’s an inherent type instability, and you’re already containing it by throwing `TypeError`s if it can’t be converted to an annotated type of another explicit variable used afterward. Without any knowledge of the large tensors you are actually `serialize`ing beforehand, we can’t suggest any type-stable alternatives for reading and writing data with any certainty.

31 and 40 are outputs of `calculate_ECCD(T2_old)`, which is potentially a problem with that method. 36 seems to lead to the `r2norm` part of `check_convergence`; it’s harder to follow because of the multiple assignment, but it appears that while the `conv` value is already inferred as `::Bool` so the `convert`/`typeassert` steps are unreachable, the `r2norm` value is `::Any`, again potentially a problem with that method. When you annotate left-hand variables, the `convert`/`typeassert` steps show up even if the right-hand expression is inferred to be the intended type already; that does complicate the lowered code but the compiler will remove redundant steps so it doesn’t happen at runtime.

The kernel function advice is good for isolating type instabilities, but you’re already doing that with annotated variables. I would first look into the functions you already have: `calculate_ECCD` and `check_convergence`. If you manage to get those type-stable, I would consider removing the variable annotations that do nothing unless they are inherently true e.g. `e_new` will always be a `Float64` no matter what methods or specializations you’ll ever run into, so you might as well annotate it and attempt the conversion if ever needed. Again, that wouldn’t change anything in runtime, but it does make lowered code easier to look at.

That means a variable is inferred to never have a value; if a return type of a method, it means the method is inferred to never return. That’s probably wrong.

1 Like

Either you have a bug in your code, causing it never to return, or you invoked `@code_warntype` wrong.