Quantum Pulse Bi-Directional Waveguide

This example mirrors the example Bi-Directional Waveguide but uses the numeric SLH struct SLHqo, to circumvent the symbolic derivation part. Furthermore it drives the system with a quantum single-photon pulse via a virtual cavity.

As usual, we start by loading the packages and defining the operators and parameters of the system.

using QuantumInputOutput
using QuantumOptics
using Plots

N = 2 # number of quantum dots

# numeric Hilbert space
bu = FockBasis(6) # virtual input cavity
ba = NLevelBasis(2)
b_qds = tensor([ba for _ = 1:N]...)
b = bu ⊗ b_qds

# QD operators
σ(i, j, k) = embed(b, i+1, transition(ba, j, k))

# input mode operator
a_u = destroy(bu) ⊗ one(b_qds)

# parameters
γ = 1.0
β = 0.9
γRn = fill(γ * β / 2, N)
γLn = fill(γ * β / 2, N)
γ_add = fill(γ * (1 - β), N)
Δn = fill(0.0, N)
ϕn = fill(π/10, max(N - 1, 0))

# quantum pulse u(t) and virtual cavity coupling g_u(t)
σt = 0.8 # pulse with
t0 = 4σt  # pulse peak
Tend = 3t0
T = [0:0.005:1;]*Tend
ΔT = T[2]-T[1]
u1(t) = 1/(sqrt(σt)*π^(1/4)) * exp(-(t - t0)^2 / (2*σt^2))

gu_t = u_to_gu(u1, T)

We use the SLHqo to directly use QuantumOptics.jl operators and functions to the model the system.

G_u = SLHqo(1, t -> gu_t(t) * a_u, 0*one(b))
G_ϕ(i, j) = SLHqo(exp(1im * ϕn[i]), 0*one(b), 0*one(b))
G_R(i) = SLHqo(1, √(γRn[i]) * σ(i, 1, 2), -Δn[i] * σ(i, 2, 2))
G_L(i) = SLHqo(1, √(γLn[i]) * σ(i, 1, 2), 0*one(b))

# Cascade right-moving channel
G_R_t = G_u ▷ G_R(1) ▷ G_ϕ(1, 2) ▷ G_R(2)

# Cascade left-moving channel (reverse order)
G_L_t = G_L(2) ▷ G_ϕ(1, 2) ▷ G_L(1)

G_t = G_R_t ⊞ G_L_t

The full Hamiltonian and Lindblad terms are extracted from the final SLH element. Note that as soon as one time-dependent function is involved in a cascade or concatenate, the returned $H$ and $L$ will also be time-dependent.

H = get_hamiltonian(G_t)
L = get_lindblad(G_t)
L_R = L[1]
L_L = L[2]

J_add = [√(γ_add[i]) * σ(i, 1, 2) for i = 1:N]

function input_output(t, ρ)
    Ht = H(t)
    J = [L_R(t), L_L(t), J_add...]
    return Ht, J, dagger.(J)
end

# time evolution
α0 = √(0.1) # √ of total photon number
ψ0 = coherentstate(bu, α0) ⊗ tensor([nlevelstate(ba, 1) for _ = 1:N]...)
t, ρt = timeevolution.master_dynamic(T, ψ0, input_output)

# transmitted and reflected intensity
I_R = zeros(length(t))
I_L = zeros(length(t))
for (i, ti) in enumerate(t)
    LR = L_R(ti)
    LL = L_L(ti)
    I_R[i] = real(expect(LR' * LR, ρt[i]))
    I_L[i] = real(expect(LL' * LL, ρt[i]))
end

p = plot(t, I_R; label = "Transmission")
plot!(p, t, I_L; label = "Reflection")
plot!(p, t, abs2.(α0*u1.(t)); color = :grey, ls = :dash, label = "Input")
plot!(
    p;
    xlabel = "time",
    ylabel = "intensity",
    legend = :best,
    grid = true,
    size = (500, 320),
)
p

Package versions

These results were obtained using the following versions:

using InteractiveUtils
versioninfo()

using Pkg
Pkg.status(["QuantumInputOutput", "QuantumOptics", "Plots"], mode = PKGMODE_MANIFEST)

Julia Version 1.12.5
Commit 5fe89b8ddc1 (2026-02-09 16:05 UTC)
Build Info:
  Official https://julialang.org release
Platform Info:
  OS: Linux (x86_64-linux-gnu)
  CPU: 4 × AMD EPYC 9V74 80-Core Processor
  WORD_SIZE: 64
  LLVM: libLLVM-18.1.7 (ORCJIT, znver4)
  GC: Built with stock GC
Threads: 4 default, 1 interactive, 4 GC (on 4 virtual cores)
Environment:
  JULIA_PKG_SERVER_REGISTRY_PREFERENCE = eager
  JULIA_DEBUG = Documenter
  JULIA_NUM_THREADS = auto
Status `~/work/QuantumInputOutput.jl/QuantumInputOutput.jl/docs/Manifest.toml`
  [91a5bcdd] Plots v1.41.6
  [18f9eda6] QuantumInputOutput v0.1.0 `~/work/QuantumInputOutput.jl/QuantumInputOutput.jl`
  [6e0679c1] QuantumOptics v1.2.4

This page was generated using Literate.jl.