fix runtests for HEOMSuperOp and spectra functions

src/HeomBase.jl

@@ -2,8 +2,8 @@ export AbstractHEOMLSMatrix
 
 abstract type AbstractHEOMLSMatrix{T} end
 
-QuantumToolbox._FType(::AbstractHEOMLSMatrix{<:AbstractArray{T}}) where {T<:Number} = _FType(T)
-QuantumToolbox._CType(::AbstractHEOMLSMatrix{<:AbstractArray{T}}) where {T<:Number} = _CType(T)
+QuantumToolbox._FType(M::AbstractHEOMLSMatrix) = _FType(eltype(M))
+QuantumToolbox._CType(M::AbstractHEOMLSMatrix) = _CType(eltype(M))
 
 _get_SciML_matrix_wrapper(M::AbstractArray) = QuantumToolbox.get_typename_wrapper(M){eltype(M)}
 _get_SciML_matrix_wrapper(M::MatrixOperator) = _get_SciML_matrix_wrapper(M.A)

src/density_of_states.jl

@@ -54,7 +54,7 @@ Calculate density of states for the fermionic system in frequency domain.
     # Handle d_op
     MType = _get_SciML_matrix_wrapper(M)
     _tr = transpose(_Tr(M))
-    Id_sys = I(prod(M.dims))
+    Id_sys = I(prod(d_op.dims))
     Id_HEOM = I(M.N)
     d_normal = HEOMSuperOp(spre(d_op, Id_sys), ODD, M; Id_cache = Id_HEOM)
     d_dagger = HEOMSuperOp(spre(d_op', Id_sys), ODD, M; Id_cache = Id_HEOM)

src/deprecated.jl

@@ -19,3 +19,16 @@ export evolution, SteadyState
 SteadyState(args...; kwargs...) = error("`SteadyState` has been deprecated, please use `steadystate` instead.")
 
 evolution(args...; kwargs...) = error("`evolution` has been deprecated, please use `HEOMsolve` instead.")
+
+_HEOMSuperOp_deprecated_method_error() = error("Specifying `mul_basis = \"L\", \"R\", or \"LR\"` to `HEOMSuperOp` has been deprecated, try to construct system SuperOperator manually by using `spre`, `spost`, or `sprepost`.")
+
+HEOMSuperOp(op, opParity::AbstractParity, dims::SVector, N::Int, mul_basis::AbstractString; Id_cache = I(N)) = _HEOMSuperOp_deprecated_method_error()
+HEOMSuperOp(
+    op,
+    opParity::AbstractParity,
+    refHEOMLS::AbstractHEOMLSMatrix,
+    mul_basis::AbstractString;
+    Id_cache = I(refHEOMLS.N),
+) = HEOMSuperOp(op, opParity, refHEOMLS.dims, refHEOMLS.N, mul_basis; Id_cache = Id_cache)
+HEOMSuperOp(op, opParity::AbstractParity, refADOs::ADOs, mul_basis::AbstractString; Id_cache = I(refADOs.N)) =
+    HEOMSuperOp(op, opParity, refADOs.dims, refADOs.N, mul_basis; Id_cache = Id_cache)

src/power_spectrum.jl

@@ -85,15 +85,14 @@ remember to set the parameters:
     end
     _check_sys_dim_and_ADOs_num(M, ados)
 
-    Id_sys = I(prod(M.dims))
     Id_HEOM = I(M.N)
 
     # Handle P_op
     if P_op isa HEOMSuperOp
         _check_sys_dim_and_ADOs_num(M, P_op)
         _P = P_op
     else
-        _P = HEOMSuperOp(spre(P_op, Id_sys), EVEN, M; Id_cache = Id_HEOM)
+        _P = HEOMSuperOp(spre(P_op), EVEN, M; Id_cache = Id_HEOM)
     end
     MType = _get_SciML_matrix_wrapper(M)
     _tr_P = transpose(_Tr(M)) * MType(_P).data
@@ -109,7 +108,7 @@ remember to set the parameters:
         else
             new_parity = !ados.parity
         end
-        _Q_ados = HEOMSuperOp(spre(Q_op, Id_sys), new_parity, M; Id_cache = Id_HEOM) * ados
+        _Q_ados = HEOMSuperOp(spre(Q_op), new_parity, M; Id_cache = Id_HEOM) * ados
     end
     b = _HandleVectorType(M, _Q_ados.data)
 

test/HEOMSuperOp.jl

@@ -45,15 +45,15 @@
     M_me = addFermionDissipator(M_me, jumpOPs)
 
     # create quantum jump superoperator (only the emission part)
-    J_me = HEOMSuperOp(γ_eR * d, ODD, M_me, "L") * HEOMSuperOp(γ_eR * d', ODD, M_me, "R")
+    J_me = HEOMSuperOp(γ_eR * spre(d), ODD, M_me) * HEOMSuperOp(γ_eR * spost(d'), ODD, M_me)
     @test size(M_me) == size(J_me)
     @test size(M_me, 1) == size(J_me, 1)
     @test eltype(M_me) == eltype(J_me)
     @test nnz(J_me.data) == 1
     @test J_me[4, 1] ≈ 0.0044570891355200214
     @test J_me.data * 2 ≈ (J_me - (-1) * J_me).data
     @test J_me.data * 2 ≈
-          (HEOMSuperOp(√(2) * γ_eR * d, ODD, M_me, "L") * HEOMSuperOp(√(2) * γ_eR * d', ODD, M_me, "R")).data
+          (HEOMSuperOp(√(2) * γ_eR * spre(d), ODD, M_me) * HEOMSuperOp(√(2) * γ_eR * spost(d'), ODD, M_me)).data
     @test show(devnull, MIME("text/plain"), J_me) == nothing
 
     # create HEOM Liouvuillian superoperator without quantum jump
@@ -77,7 +77,7 @@
     M_heom = addFermionDissipator(M_heom, jumpOPs)
 
     # create quantum jump superoperator (only the emission part)
-    J_heom = HEOMSuperOp(γ_eR * d, ODD, M_heom, "L") * HEOMSuperOp(γ_eR * d', ODD, M_heom, "R")
+    J_heom = HEOMSuperOp(γ_eR * spre(d), ODD, M_heom) * HEOMSuperOp(γ_eR * spost(d'), ODD, M_heom)
 
     # create HEOM Liouvuillian superoperator without quantum jump
     M1 = M_heom - J_heom
@@ -93,9 +93,9 @@
     end
 
     ## test for error
-    J_wrong1 = HEOMSuperOp(γ_eR * d, ODD, M_me, "L")
-    J_wrong2 = HEOMSuperOp(γ_eR * d, EVEN, M_heom, "L")
-    @test_throws ErrorException HEOMSuperOp(d, EVEN, M_heom, "LR")
+    J_wrong1 = HEOMSuperOp(γ_eR * spre(d), ODD, M_me)
+    J_wrong2 = HEOMSuperOp(γ_eR * spre(d), EVEN, M_heom)
+    @test_throws ErrorException HEOMSuperOp(sprepost(d), EVEN, M_heom)
     @test_throws ErrorException J_me * J_wrong2
     @test_throws ErrorException J_me + J_wrong1
     @test_throws ErrorException J_me + J_wrong2
@@ -105,4 +105,6 @@
     @test_throws ErrorException J_me - J_wrong2
     @test_throws ErrorException M_me - J_wrong1
     @test_throws ErrorException M_me - J_wrong2
+    @test_throws ErrorException HEOMSuperOp(d, ODD, M_heom, "L")
+    @test_throws ErrorException HEOMSuperOp(d, ODD, ados_s1, "L")
 end

test/density_of_states.jl

@@ -30,7 +30,7 @@
     if isfile("DOS.txt")
         rm("DOS.txt")
     end
-    d_up_normal = HEOMSuperOp(d_up, ODD, Le)
+    d_up_normal = HEOMSuperOp(spre(d_up), ODD, Le)
     dos1 = DensityOfStates(Lo, ados_s, d_up, ωlist; verbose = false, filename = "DOS")
     dos2 =
         PowerSpectrum(Lo, ados_s, d_up_normal, d_up', ωlist, true; verbose = false) .+

test/power_spectrum.jl

@@ -51,8 +51,8 @@
 
     mat = Qobj(spzeros(ComplexF64, 2, 2))
     mat2 = Qobj(spzeros(ComplexF64, 3, 3))
-    a_even = HEOMSuperOp(a, EVEN, L)
-    a_odd = HEOMSuperOp(a, ODD, L)
+    a_even = HEOMSuperOp(spre(a), EVEN, L)
+    a_odd = HEOMSuperOp(spre(a), ODD, L)
     bathf = Fermion_Lorentz_Pade(mat, 1, 1, 1, 1, 2)
     @test_throws ErrorException PowerSpectrum(L, ados_s, a, ωlist; verbose = false, filename = "PSD")
     @test_throws ErrorException PowerSpectrum(L, ados_s, a_even, a_odd, ωlist; verbose = false)

test/stationary_state.jl

@@ -37,5 +37,5 @@
     @test_throws ErrorException steadystate(L, mat2)
     @test_throws ErrorException steadystate(L, ADOs(zeros(8), 2))
     @test_throws ErrorException steadystate(L, ADOs(ados.data, ados.N, ODD))
-    @test_throws ErrorException steadystate(L, HEOMSuperOp(d, ODD, L) * ados)
+    @test_throws ErrorException steadystate(L, HEOMSuperOp(spre(d), ODD, L) * ados)
 end

test/time_evolution.jl

@@ -32,7 +32,7 @@
     ados_wrong1 = ADOs(zeros(8), 2)
     ados_wrong2 = ADOs(zeros(32), 2)
     ados_wrong3 = ADOs((ados_list[1]).data, (ados_list[1]).N, ODD)
-    ados_wrong4 = HEOMSuperOp(Q, ODD, ados_list[end]) * ados_list[end]
+    ados_wrong4 = HEOMSuperOp(spre(Q), ODD, ados_list[end]) * ados_list[end]
     ρ_list_p = getRho.(ados_list)
     @test show(devnull, MIME("text/plain"), sol_p) === nothing
     @test length(sol_p.ados) == 1