build: move to SQA v0.10

Project.toml

@@ -29,10 +29,10 @@ LinearAlgebra = "1.10"
 ModelingToolkitBase = "1.42.2"
 NumericalIntegration = "0.3"
 OrdinaryDiffEq = "6"
-QuantumCumulants = "0.5"
+QuantumCumulants = "0.7"
 QuantumOptics = "1"
 QuantumOpticsBase = "0.5"
-SecondQuantizedAlgebra = "0.9"
+SecondQuantizedAlgebra = "0.10"
 SpecialFunctions = "2"
 StaticArrays = "1"
 SymbolicUtils = "4"

docs/Project.toml

@@ -24,8 +24,9 @@ MacroTools = "0.5"
 ModelingToolkitBase = "1"
 OrdinaryDiffEq = "6"
 Plots = "1"
-QuantumCumulants = "0.5"
+QuantumCumulants = "0.7"
 QuantumOptics = "1"
 QuantumOpticsBase = "0.5"
+SecondQuantizedAlgebra = "0.10"
 SymbolicUtils = "4"
 Symbolics = "7"

examples/05-1_N-QDs_bidirectional-waveguide_coherent-pulse.jl

@@ -8,7 +8,7 @@ using QuantumInputOutput
 using SecondQuantizedAlgebra
 using Symbolics: Symbolics
 using QuantumOptics
-using QuantumOpticsBase: dagger
+using QuantumOpticsBase: dagger, static_operator
 using Plots
 using LaTeXStrings
 
@@ -159,9 +159,12 @@ lT = length(T)
 G2 = zeros(lT, lT) # transmission
 G2_ref = zeros(lT, lT) # reflection
 
-L0(t) = L_R_QO(t)
+# `L_R_QO(t)`/`L_L_QO(t)` return a lazy `TimeDependentSum`; materialize it to a concrete
+# operator at each time so the quantum-regression products below yield a plain operator
+# (a `TimeDependentSum` product cannot serve as the solver's initial state).
+L0(t) = dense(static_operator(L_R_QO(t)))
 L0_dag(t) = dagger(L0(t))
-L0_ref(t) = L_L_QO(t)
+L0_ref(t) = dense(static_operator(L_L_QO(t)))
 L0_ref_dag(t) = dagger(L0_ref(t))
 
 for it1 = 1:(lT-1)

examples/08-1_pulse-delay__simple.jl

@@ -190,6 +190,9 @@ p_t_sym = [g_u, g_in, g_out, g_v, M_ls...]
 p_t_num = [gu_, gin_, gout_, gv_, M_t_ls...]
 dict_p_t_int = Dict(p_t_sym .=> p_t_num)
 
+# The interaction picture eliminates the delay cavity `d`, so the numeric operators live
+# on the two-mode basis `bu ⊗ bv` (with `ad` mapped to the identity). Build on that basis.
+b_int = bu ⊗ bv
 au_int = destroy(bu) ⊗ one(bv)
 ## This is amazing! # hide
 ad_int = one(bu ⊗ bv)
@@ -198,13 +201,13 @@ operators = Dict(
     [au, au', ad, ad', av, av'] .=> [au_int, au_int', ad_int, ad_int', av_int, av_int'],
 )
 
-H_int_QO = to_numeric(H_int_sym, b; time_parameter = dict_p_t_int, operators)
+H_int_QO = to_numeric(H_int_sym, b_int; time_parameter = dict_p_t_int, operators)
 L_int_QO = [
-    to_numeric(L_int_sym[i], b; time_parameter = dict_p_t_int, operators) for
+    to_numeric(L_int_sym[i], b_int; time_parameter = dict_p_t_int, operators) for
     i = 1:length(L_int_sym)
 ]
-## H_int_QO = to_numeric(H_int_sym, b; time_parameter=dict_p_t_int) # hide 
-## L_int_QO = [to_numeric(L_int_sym[i], b; time_parameter=dict_p_t_int) for i=1:length(L_int_sym)] # hide
+## H_int_QO = to_numeric(H_int_sym, b_int; time_parameter=dict_p_t_int) # hide
+## L_int_QO = [to_numeric(L_int_sym[i], b_int; time_parameter=dict_p_t_int) for i=1:length(L_int_sym)] # hide
 nothing # hide
 
 #

examples/drafts/02-2_traveling-cat-state__PRA2020_102-023717_fig3b.jl

@@ -13,6 +13,7 @@
 using QuantumInputOutput
 using SecondQuantizedAlgebra
 using QuantumOptics
+using QuantumOpticsBase: dagger
 using Plots
 using LinearAlgebra
 using DataInterpolations
@@ -23,7 +24,7 @@ using DataInterpolations
 hc = FockSpace(:c)
 
 ## symbolic operator
-a = Destroy(hc, :a, 1)
+a = Destroy(hc, :a)
 
 ## symbolic parameters
 @variables γ::Real K::Real Δ::Real p::Number

examples/drafts/08-2_pulse-propagation-delay__PRA2026_113-013730_fig4.jl

@@ -11,6 +11,7 @@ using QuantumInputOutput
 using SecondQuantizedAlgebra
 using Symbolics: Symbolics
 using QuantumOptics
+using QuantumOpticsBase: dagger
 using SymbolicUtils
 using LinearAlgebra
 using Plots
@@ -181,9 +182,9 @@ t_, ρt = timeevolution.master_dynamic(T, ψ0_fock, input_output)
 
 expect(au_QO'au_QO, ρt)
 
-p = plot(T, expect(au_QO' * au_QO, ρt); label = "")
-plot!(p, T, expect(ad1_QO' * ad1_QO, ρt); label = "")
-plot!(p, T, expect(ad2_QO' * ad2_QO, ρt); label = "")
+p = plot(T, real.(expect(au_QO' * au_QO, ρt)); label = "")
+plot!(p, T, real.(expect(ad1_QO' * ad1_QO, ρt)); label = "")
+plot!(p, T, real.(expect(ad2_QO' * ad2_QO, ρt)); label = "")
 p
 
 
@@ -197,8 +198,6 @@ expect(σ22_QO, ρt[end])
 pop_fock = zeros(length(Δ_ls))
 pop_coh = zeros(length(Δ_ls))
 idx_t1 = findmin(abs.(T .- t1))[2]
-using ProgressMeter
-prog = Progress(length(Δ_ls))
 
 for (it, Δn) in enumerate(Δ_ls)
     Δn = Δ_ls[it]
@@ -219,7 +218,6 @@ for (it, Δn) in enumerate(Δ_ls)
 
     # t_, ρt = timeevolution.master_dynamic(T, ψ0_coh, input_output)
     # pop_coh[it] = real(expect(σ22_QO, ρt[idx_t1]))
-    next!(prog)
 end
 nothing # hide
 

src/pulses.jl

@@ -6,6 +6,15 @@ const _tol_div = 1e-10
 const _extrapolate = ExtrapolationType.Extension
 const _ϵ = 1e-10
 
+# Wrap a sampled interpolation in an explicit single-argument closure `t -> g(t)`.
+# SecondQuantizedAlgebra v0.10 rejects `time_parameter` values whose callable has
+# ambiguous arity (a `LinearInterpolation` responds to both `g(t)` and `g(t, order)`),
+# so the coupling functions hand back an unambiguous `t -> value` closure that plugs
+# straight into `to_numeric`'s `time_parameter` while still being callable as `g(t)`.
+_as_time_fn(itp) = let g = itp
+    t -> g(t)
+end
+
 _mode_interp(mode::AbstractVector, T::AbstractVector) =
     LinearInterpolation(mode, T; extrapolation = _extrapolate)
 _mode_interp(mode, T::AbstractVector) = mode
@@ -54,7 +63,7 @@ function _compute_coupling(mode::Vector, T::Vector, denom_fn)
             g[i] = mode[i]' / sqrt(d)
         end
     end
-    return LinearInterpolation(g, T; extrapolation = _extrapolate)
+    return _as_time_fn(LinearInterpolation(g, T; extrapolation = _extrapolate))
 end
 
 # ──────────────────────────────────────────────
@@ -65,7 +74,8 @@ end
     coupling_input(u, T)
 
 Compute the virtual-cavity input coupling ``g_u(t)`` from an input mode `u(t)`
-sampled on time grid `T`. Returns the `LinearInterpolation` directly (callable as `g(t)`).
+sampled on time grid `T`. Returns a single-argument closure `g(t)` (an interpolation of
+the sampled coupling) that plugs directly into a `to_numeric` `time_parameter`.
 """
 coupling_input(u::Vector, T::Vector) = _compute_coupling(u, T, x -> abs(1 - x) + _ϵ)
 coupling_input(u::Function, T::Vector) = coupling_input(u.(T), T)
@@ -85,7 +95,8 @@ end
     coupling_output(v, T)
 
 Compute the virtual-cavity output coupling ``g_v(t)`` from an output mode `v(t)`
-sampled on time grid `T`. Returns the `LinearInterpolation` directly (callable as `g(t)`).
+sampled on time grid `T`. Returns a single-argument closure `g(t)` (an interpolation of
+the sampled coupling) that plugs directly into a `to_numeric` `time_parameter`.
 """
 coupling_output(v::Vector, T::Vector) = _compute_coupling(-v, T, x -> x + _ϵ)
 coupling_output(v::Function, T::Vector) = coupling_output(v.(T), T)
@@ -290,14 +301,14 @@ function _compute_coupling_delay(num_mode::Vector, u::Vector, v::Vector, T::Vect
             g[i] = num_mode[i]' / sqrt(d + _ϵ)
         end
     end
-    return LinearInterpolation(g, T; extrapolation = _extrapolate)
+    return _as_time_fn(LinearInterpolation(g, T; extrapolation = _extrapolate))
 end
 
 """
     coupling_delay_out(u, v, T)
 
 Compute the out-coupling strength for a delay cavity.
-Returns `LinearInterpolation` directly.
+Returns a single-argument closure `g(t)`.
 """
 coupling_delay_out(u::Vector, v::Vector, T::Vector) = _compute_coupling_delay(u, u, v, T)
 coupling_delay_out(u::Function, v::Function, T::Vector) =
@@ -309,7 +320,7 @@ coupling_delay_out(u::LinearInterpolation, v::LinearInterpolation, T::Vector) =
     coupling_delay_in(u, v, T)
 
 Compute the in-coupling strength for a delay cavity.
-Returns `LinearInterpolation` directly.
+Returns a single-argument closure `g(t)`.
 """
 coupling_delay_in(u::Vector, v::Vector, T::Vector) = _compute_coupling_delay(-v, u, v, T)
 coupling_delay_in(u::Function, v::Function, T::Vector) = coupling_delay_in(u.(T), v.(T), T)

test/test_translate.jl

@@ -2,9 +2,15 @@ using QuantumInputOutput
 using SecondQuantizedAlgebra
 using SymbolicUtils
 using QuantumOpticsBase
-using QuantumOpticsBase: dagger
+using QuantumOpticsBase: dagger, static_operator, TimeDependentSum
 using Test
 
+# Since SecondQuantizedAlgebra v0.10, `to_numeric` returns a `TimeDependentSum` for
+# time-dependent inputs instead of a plain callable. Calling it at a time `t` sets the
+# current time and still returns a `TimeDependentSum`; `mat` materializes it to a dense
+# `Operator` for numeric comparison.
+mat(F, t) = dense(static_operator(F(t)))
+
 @testset "translate" begin
     @variables κ_L::Real κ_R::Real Δ::Real g::Real γ::Real
     @variables E::Complex E1::Complex E2::Complex
@@ -56,16 +62,16 @@ using Test
             time_parameter = dict_p_t2,
             operators = ops_dict,
         )
-        @test isa(F1, Function)
-        @test isequal(F1(0.1), a_QO*3)
+        @test F1 isa TimeDependentSum
+        @test sum(abs.((mat(F1, 0.1) - dense(a_QO*3)).data)) < 1e-12
         F2 = to_numeric(
             a*E,
             bc1;
             parameter = dict_p1,
             time_parameter = dict_p_t2,
             operators = ops_dict,
         )
-        @test isequal(F2(0.1), a_QO*E_t(0.1))
+        @test sum(abs.((mat(F2, 0.1) - dense(a_QO*E_t(0.1))).data)) < 1e-12
         # association of the scalar prefactor differs from `Δn*A'*A`, compare numerically
         @test sum(
             abs.(
@@ -101,38 +107,37 @@ using Test
         abs.((dense(to_numeric(Δ, b; parameter = dict_p1)) - dense(one(b)*Δn)).data),
     ) < 1e-12
     F3 = to_numeric(Δ, b; parameter = dict_p1, time_parameter = dict_p_t2)
-    @test sum(abs.((dense(F3(4)) - dense(one(b)*Δn)).data)) < 1e-8
+    @test sum(abs.((mat(F3, 4) - dense(one(b)*Δn)).data)) < 1e-8
     F4 = to_numeric(
         a*3*conj(E) + Δ*σ(2, 2),
         b;
         parameter = dict_p1,
         time_parameter = dict_p_t2,
     )
-    @test sum(abs.((F4(0.2) - dense(a_QO2*3*E_t_c(0.2) + Δn*σ_QO(2, 2))).data)) < 1e-8
+    @test sum(abs.((mat(F4, 0.2) - dense(a_QO2*3*E_t_c(0.2) + Δn*σ_QO(2, 2))).data)) < 1e-8
     F5 = to_numeric(
         a*conj(E) + Δ*σ(2, 2),
         b;
         parameter = dict_p1,
         time_parameter = dict_p_t2,
     )
-    @test sum(abs.((F5(0.2) - dense(a_QO2*E_t_c(0.2) + Δn*σ_QO(2, 2))).data)) < 1e-8
-    @inferred F5(0.2)  # QAdd time-dependent path should return concrete type
+    @test sum(abs.((mat(F5, 0.2) - dense(a_QO2*E_t_c(0.2) + Δn*σ_QO(2, 2))).data)) < 1e-8
     F5_ = to_numeric(a*conj(E), b; parameter = dict_p1, time_parameter = dict_p_t2)
-    @test sum(abs.((dense(F5_(0.2)) - dense(a_QO2*E_t_c(0.2))).data)) < 1e-8
+    @test sum(abs.((mat(F5_, 0.2) - dense(a_QO2*E_t_c(0.2))).data)) < 1e-8
     F6 = to_numeric(conj(E), b; parameter = dict_p1, time_parameter = dict_p_t2)
-    @test sum(abs.((dense(F6(0.2)) - dense(E_t_c(0.2)*one(b))).data)) < 1e-8
+    @test sum(abs.((mat(F6, 0.2) - dense(E_t_c(0.2)*one(b))).data)) < 1e-8
     F7 = to_numeric(conj(E) + Δ*σ(2, 2), b; parameter = dict_p1, time_parameter = dict_p_t2)
-    @test sum(abs.((F7(0.2) - dense(E_t_c(0.2)*one(b) + Δn*σ_QO(2, 2))).data)) < 1e-8
+    @test sum(abs.((mat(F7, 0.2) - dense(E_t_c(0.2)*one(b) + Δn*σ_QO(2, 2))).data)) < 1e-8
     # a bare symbolic scalar without a numeric/time value cannot be translated
     @test_throws ArgumentError to_numeric(conj(E), b; parameter = dict_p1)
     F8 = to_numeric(E^2, b; parameter = dict_p1, time_parameter = dict_p_t2)
-    @test sum(abs.((dense(F8(0.2)) - dense(E_t(0.2)^2*one(b))).data)) < 1e-8
+    @test sum(abs.((mat(F8, 0.2) - dense(E_t(0.2)^2*one(b))).data)) < 1e-8
 
 
     @testset "time_parameter_normalization" begin
         dict_p_t_num = Dict([E] .=> [2.5])
         F_num = to_numeric(a*E, b; parameter = dict_p1, time_parameter = dict_p_t_num)
-        @test sum(abs.((dense(F_num(0.2)) - dense(a_QO2 * 2.5)).data)) < 1e-8
+        @test sum(abs.((mat(F_num, 0.2) - dense(a_QO2 * 2.5)).data)) < 1e-8
     end
 
     @testset "multiple_time_prefactors" begin
@@ -147,7 +152,7 @@ using Test
             time_parameter = dict_p_t_multi,
         )
         expected = dense(a_QO2 * E1_t(0.4) * E2_t_c(0.4))
-        @test sum(abs.((dense(F_multi(0.4)) - expected).data)) < 1e-8
+        @test sum(abs.((mat(F_multi, 0.4) - expected).data)) < 1e-8
     end
 
     @testset "to_numeric vector overload" begin
@@ -160,10 +165,10 @@ using Test
             operators = ops_dict,
         )
         @test length(translated_ops) == length(ops)
-        @test all(op -> op isa Function, translated_ops)
-        @test isequal(translated_ops[1](0.3), a_QO)
-        @test isequal(translated_ops[2](0.3), one(bc1) * Δn)
-        @test isequal(translated_ops[3](0.3), a_QO * E_t(0.3))
+        @test all(op -> op isa TimeDependentSum, translated_ops)
+        @test sum(abs.((mat(translated_ops[1], 0.3) - dense(a_QO)).data)) < 1e-12
+        @test sum(abs.((mat(translated_ops[2], 0.3) - dense(one(bc1) * Δn)).data)) < 1e-12
+        @test sum(abs.((mat(translated_ops[3], 0.3) - dense(a_QO * E_t(0.3))).data)) < 1e-12
     end
 
     @testset "to_numeric SLH overload" begin
@@ -176,11 +181,11 @@ using Test
             operators = ops_dict,
         )
 
-        @test H_QO isa Function
+        @test H_QO isa TimeDependentSum
         @test length(L_QO) == 1
-        @test L_QO[1] isa Function
-        @test sum(abs.((H_QO(0.4) - dense(Δn * dagger(a_QO) * a_QO)).data)) < 1e-8
-        @test sum(abs.((L_QO[1](0.4) - dense(sqrt(κ_Rn) * a_QO * E_t(0.4))).data)) < 1e-8
+        @test L_QO[1] isa TimeDependentSum
+        @test sum(abs.((mat(H_QO, 0.4) - dense(Δn * dagger(a_QO) * a_QO)).data)) < 1e-8
+        @test sum(abs.((mat(L_QO[1], 0.4) - dense(sqrt(κ_Rn) * a_QO * E_t(0.4))).data)) < 1e-8
     end
 
     @testset "substitute operators qmul" begin
@@ -231,7 +236,7 @@ using Test
         @variables gR::Real
         gR_t(t) = 1.0 + 2.0im
         F = to_numeric(im * gR * a3, b3; time_parameter = Dict(gR => gR_t))
-        @test F isa Function
-        @test sum(abs.((F(0.0) - dense((im * (1.0 + 2.0im)) * a3_QO)).data)) < 1e-8
+        @test F isa TimeDependentSum
+        @test sum(abs.((mat(F, 0.0) - dense((im * (1.0 + 2.0im)) * a3_QO)).data)) < 1e-8
     end
 end