Weird behaviour of steadystate.eigenvector function returns any while it should return the type operator

Specific Domains Quantum
1
  • I was doing a simulation of calculating the population of a 3-state system interacting with a cavity mode.

Part of the code is given as follows

function pop(g, Ω23, Ω12, Δ23_list, Δ12_list, Δ, Γ31, Γ32, 
        κ, a, ad, σ13, σ23, σ31, σ32, σ33, σ11, σ12, σ21)
    m = length(Δ23_list)
    n = length(Δ12_list)
    popu = zeros(eltype(Float64),m, n)

    # Precompute constant part of V
    V_const = Ω12 * (σ12 + σ21)
    H0_val = g*(a * σ31 +  ad * σ13) + Ω23 * (σ32 + σ23) - Δ*ad*a
    J = [sqrt(Γ31) * σ13, sqrt(Γ32) * σ23, sqrt(κ) * a]

    for (l, Δ12) in enumerate(Δ12_list)
        for (k, Δ23) in enumerate(Δ23_list)
            # Compute H0
            H0delta = Δ23 * σ33 - Δ12 * σ11 + (Δ12 + Δ23) * ad * a
            # Compute H
            H_val = H0_val + V_const + H0delta
            H_val1 = DenseOperator(H_val)
            # Compute steady state eigenvector and population
            ρ = steadystate.eigenvector(H_val1, J)
            popu[k, l] += real(expect(σ33, ρ))
        end
    end

    return popu
end

when I typed @code_warntype I get \rho is of type Any but it should be of type Operator. Also the code is very slow allocates lots of memory. Possibly the type stability problem of my code in the QuantumOptics.jl package. Where am I getting wrong? The output of the @code_warntype is as follows:

@code_warntype pop(g, Ω23, Ω12, Δ23_list[1:2], Δ12_list[1:2], Δ, Γ31, Γ32, 
κ, a, ad, σ13, σ23, σ31, σ32, σ33, σ11, σ12, σ21)

out: MethodInstance for pop(::Float64, ::Float64, ::Float64, ::StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, ::StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, ::Float64, ::Float64, ::Float64, ::Float64, ::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}, ::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}, ::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}, ::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}, ::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}, ::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}, ::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}, ::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}, ::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}, ::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}})
  from pop(g, Ω23, Ω12, Δ23_list, Δ12_list, Δ, Γ31, Γ32, κ, a, ad, σ13, σ23, σ31, σ32, σ33, σ11, σ12, σ21) in Main at In[10]:1
Arguments
  #self#::Core.Const(pop)
  g::Float64
  Ω23::Float64
  Ω12::Float64
  Δ23_list::StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}
  Δ12_list::StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}
  Δ::Float64
  Γ31::Float64
  Γ32::Float64
  κ::Float64
  a::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
  ad::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
  σ13::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
  σ23::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
  σ31::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
  σ32::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
  σ33::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
  σ11::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
  σ12::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
  σ21::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
Locals
  @_21::Union{Nothing, Tuple{Tuple{Int64, Float64}, Tuple{Int64, Int64}}}
  J::Vector{Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}}
  H0_val::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
  V_const::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
  popu::Matrix{Float64}
  n::Int64
  m::Int64
  @_28::Union{Nothing, Tuple{Tuple{Int64, Float64}, Tuple{Int64, Int64}}}
  @_29::Int64
  Δ12::Float64
  l::Int64
  @_32::Int64
  Δ23::Float64
  k::Int64
  ρ::Any
  H_val1::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, Matrix{ComplexF64}}
  H_val::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
  H0delta::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
Body::Matrix{Float64}
1 ─       (m = Main.length(Δ23_list))
│         (n = Main.length(Δ12_list))
│   %3  = Main.eltype(Main.Float64)::Core.Const(Float64)
│   %4  = m::Int64
│         (popu = Main.zeros(%3, %4, n))
│   %6  = (σ12 + σ21)::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
│         (V_const = Ω12 * %6)
│   %8  = (a * σ31)::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
│   %9  = (ad * σ13)::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
│   %10 = (%8 + %9)::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
│   %11 = (g * %10)::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
│   %12 = (σ32 + σ23)::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
│   %13 = (Ω23 * %12)::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
│   %14 = (%11 + %13)::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
│   %15 = (Δ * ad * a)::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
│         (H0_val = %14 - %15)
│   %17 = Main.sqrt(Γ31)::Float64
│   %18 = (%17 * σ13)::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
│   %19 = Main.sqrt(Γ32)::Float64
│   %20 = (%19 * σ23)::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
│   %21 = Main.sqrt(κ)::Float64
│   %22 = (%21 * a)::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
│         (J = Base.vect(%18, %20, %22))
│   %24 = Main.enumerate(Δ12_list)::Base.Iterators.Enumerate{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}}
│         (@_21 = Base.iterate(%24))
│   %26 = (@_21 === nothing)::Bool
│   %27 = Base.not_int(%26)::Bool
└──       goto #7 if not %27
2 ┄ %29 = @_21::Tuple{Tuple{Int64, Float64}, Tuple{Int64, Int64}}
│   %30 = Core.getfield(%29, 1)::Tuple{Int64, Float64}

│   %31 = Base.indexed_iterate(%30, 1)::Core.PartialStruct(Tuple{Int64, Int64}, Any[Int64, Core.Const(2)])
│         (l = Core.getfield(%31, 1))
│         (@_29 = Core.getfield(%31, 2))
│   %34 = Base.indexed_iterate(%30, 2, @_29::Core.Const(2))::Core.PartialStruct(Tuple{Float64, Int64}, Any[Float64, Core.Const(3)])
│         (Δ12 = Core.getfield(%34, 1))
│   %36 = Core.getfield(%29, 2)::Tuple{Int64, Int64}
│   %37 = Main.enumerate(Δ23_list)::Base.Iterators.Enumerate{StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}}
│         (@_28 = Base.iterate(%37))
│   %39 = (@_28 === nothing)::Bool
│   %40 = Base.not_int(%39)::Bool
└──       goto #5 if not %40
3 ┄ %42 = @_28::Tuple{Tuple{Int64, Float64}, Tuple{Int64, Int64}}
│   %43 = Core.getfield(%42, 1)::Tuple{Int64, Float64}
│   %44 = Base.indexed_iterate(%43, 1)::Core.PartialStruct(Tuple{Int64, Int64}, Any[Int64, Core.Const(2)])
│         (k = Core.getfield(%44, 1))
│         (@_32 = Core.getfield(%44, 2))
│   %47 = Base.indexed_iterate(%43, 2, @_32::Core.Const(2))::Core.PartialStruct(Tuple{Float64, Int64}, Any[Float64, Core.Const(3)])
│         (Δ23 = Core.getfield(%47, 1))
│   %49 = Core.getfield(%42, 2)::Tuple{Int64, Int64}
│   %50 = (Δ23 * σ33)::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
│   %51 = (Δ12 * σ11)::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
│   %52 = (%50 - %51)::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
│   %53 = (Δ12 + Δ23)::Float64
│   %54 = (%53 * ad * a)::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, SparseArrays.SparseMatrixCSC{ComplexF64, Int64}}
│         (H0delta = %52 + %54)
│         (H_val = H0_val + V_const + H0delta)
│         (H_val1 = Main.DenseOperator(H_val))
│   %58 = QuantumOptics.steadystate.eigenvector::Core.Const(QuantumOptics.steadystate.eigenvector)
│   %59 = H_val1::Operator{CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, CompositeBasis{Vector{Int64}, Tuple{NLevelBasis{Int64}, FockBasis{Int64}}}, Matrix{ComplexF64}}
│         (ρ = (%58)(%59, J))
│   %61 = Base.getindex(popu, k, l)::Float64
│   %62 = Main.expect(σ33, ρ)::Any
│   %63 = Main.real(%62)::Any
│   %64 = (%61 + %63)::Any
│         Base.setindex!(popu, %64, k, l)
│         (@_28 = Base.iterate(%37, %49))
│   %67 = (@_28 === nothing)::Bool
│   %68 = Base.not_int(%67)::Bool
└──       goto #5 if not %68
4 ─       goto #3
5 ┄       (@_21 = Base.iterate(%24, %36))
│   %72 = (@_21 === nothing)::Bool
│   %73 = Base.not_int(%72)::Bool
└──       goto #7 if not %73
6 ─       goto #2
7 ┄       return popu

The only problem lies in \rho. Without using DenseOperator I get Arpack()
error. Will anyone please suggest me some solution of it? Thanks in advance.
The fockstate cutoff=4 and the other parameters are given as follows:

#params
κ=2*pi*1.0;
g=10.0*κ;
Ω23=3.0*κ;
Ω12=0.1*κ;
Γ31=0.5*κ;
Γ32=Γ31;
Δ23_list=2*pi*range(-30.0, stop=30.0, length=41)
Δ12_list=2*pi*range(-5.0, stop=5.0, length=41)
Δ=0.0;

The state and operators are defined as

b1=NLevelBasis(3)
b2=FockBasis(4)
b=b1 ⊗ b2;

σ31= transition(b1,3,1) ⊗ one(b2)
σ13= transition(b1,1,3) ⊗ one(b2)
σ32= transition(b1,3,2) ⊗ one(b2)
σ23= transition(b1,2,3) ⊗ one(b2)
σ12= transition(b1,1,2) ⊗ one(b2)
σ21= transition(b1,2,1) ⊗ one(b2)
σ33= transition(b1,3,3) ⊗ one(b2)
σ11= transition(b1,1,1) ⊗ one(b2)
a= one(b1) ⊗ destroy(b2)
ad= one(b1) ⊗ create(b2);