Skip to content

Generalize the definition of liouvillian #541

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
wants to merge 6 commits into
base: main
Choose a base branch
from
Open

Generalize the definition of liouvillian #541

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
6 commits
Select commit Hold shift + click to select a range
ff1ac41
generalize the definition of `liouvillian`
ytdHuang Sep 3, 2025
dc40ca9
minor changes
ytdHuang Sep 3, 2025
affcfc9
Fix type instabilities
albertomercurio Sep 9, 2025
a17ce34
Fix dissipator init error
albertomercurio Sep 9, 2025
d8711a0
Fix errors
albertomercurio Sep 15, 2025
c035fc9
Fix wrong adjoint for ScaledOperator
albertomercurio Sep 15, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
3 changes: 3 additions & 0 deletions CHANGELOG.md
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -7,6 +7,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0

## [Unreleased](https://github.com/qutip/QuantumToolbox.jl/tree/main)

- Generalize the definition of `liouvillian`. It no longer expects the Hamiltonian to be Hermitian. ([#541])

## [v0.35.0]
Release date: 2025-09-03

Expand Down Expand Up @@ -315,4 +317,5 @@ Release date: 2024-11-13
[#536]: https://github.com/qutip/QuantumToolbox.jl/issues/536
[#537]: https://github.com/qutip/QuantumToolbox.jl/issues/537
[#539]: https://github.com/qutip/QuantumToolbox.jl/issues/539
[#541]: https://github.com/qutip/QuantumToolbox.jl/issues/541
[#546]: https://github.com/qutip/QuantumToolbox.jl/issues/546
24 changes: 5 additions & 19 deletions src/qobj/superoperators.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -41,9 +41,10 @@ end

## intrinsic liouvillian
_liouvillian(H::MT, Id::AbstractMatrix) where {MT<:Union{AbstractMatrix,AbstractSciMLOperator}} =
-1im * (_spre(H, Id) - _spost(H, Id))
-1im * (_spre(H, Id) - _spost(H', Id))
_liouvillian(H::MatrixOperator, Id::AbstractMatrix) = MatrixOperator(_liouvillian(H.A, Id))
_liouvillian(H::ScaledOperator, Id::AbstractMatrix) = ScaledOperator(H.λ, _liouvillian(H.L, Id))
Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We can remove the implementation for ScaledOperator, since it should just call the first method (for AbstractSciMLOperator)

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Actually by removing this we have issues with Enzyme.jl

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Ok we actually need it, otherwise we don't apply conj to the ScalarOperator. This creates some errors in Enzyme.jl

_liouvillian(H::ScaledOperator, Id::AbstractMatrix) =
-1im * (ScaledOperator(H.λ, _spre(H.L, Id)) - ScaledOperator(conj(H.λ), _spost(H.L', Id)))
_liouvillian(H::AddedOperator, Id::AbstractMatrix) = AddedOperator(map(op -> _liouvillian(op, Id), H.ops))

# intrinsic lindblad_dissipator
Expand Down Expand Up @@ -144,7 +145,7 @@ lindblad_dissipator(O::AbstractQuantumObject{SuperOperator}, Id_cache = nothing)
Construct the Liouvillian [`SuperOperator`](@ref) for a system Hamiltonian ``\hat{H}`` and a set of collapse operators ``\{\hat{C}_n\}_n``:

```math
\mathcal{L} [\cdot] = -i[\hat{H}, \cdot] + \sum_n \mathcal{D}(\hat{C}_n) [\cdot]
\mathcal{L} [\cdot] = -i\left(\hat{H}[\cdot] - [\cdot]\hat{H}^\dagger\right) + \sum_n \mathcal{D}(\hat{C}_n) [\cdot]
```

where
Expand Down Expand Up @@ -179,19 +180,4 @@ liouvillian(H::AbstractQuantumObject{Operator}, Id_cache::Diagonal = I(prod(H.di

liouvillian(H::AbstractQuantumObject{SuperOperator}, Id_cache::Diagonal) = H

function _sum_lindblad_dissipators(c_ops, Id_cache::Diagonal)
D = 0
# sum all the (time-independent) c_ops first
c_ops_ti = filter(op -> isa(op, QuantumObject), c_ops)
if !isempty(c_ops_ti)
D += mapreduce(op -> lindblad_dissipator(op, Id_cache), +, c_ops_ti)
end

# sum rest of the QobjEvo together
c_ops_td = filter(op -> isa(op, QuantumObjectEvolution), c_ops)
if !isempty(c_ops_td)
D += mapreduce(op -> lindblad_dissipator(op, Id_cache), +, c_ops_td)
end

return D
end
_sum_lindblad_dissipators(c_ops, Id_cache::Diagonal) = sum(op -> lindblad_dissipator(op, Id_cache), c_ops; init = 0)
9 changes: 4 additions & 5 deletions test/core-test/quantum_objects_evo.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -207,7 +207,7 @@
X = a * a'
c_op1 = QobjEvo(a', coef1)
c_op2 = QobjEvo(((a, coef2), (X, coef3)))
c_ops = [c_op1, c_op2]
c_ops = (c_op1, c_op2)
D1_ti = abs2(coef1(p, t)) * lindblad_dissipator(a')
D2_ti =
abs2(coef2(p, t)) * lindblad_dissipator(a) + # normal dissipator for first element in c_op2
Expand Down Expand Up @@ -237,10 +237,9 @@
@test_throws ArgumentError cache_operator(L_td, ψ)

@testset "Type Inference" begin
# we use destroy and create here because they somehow causes type instability before
H_td2 = H_td + QobjEvo(destroy(N) + create(N), coef3)
c_ops1 = (destroy(N), create(N))
c_ops2 = (destroy(N), QobjEvo(create(N), coef1))
H_td2 = H_td + QobjEvo(a + a', coef3)
c_ops1 = (a, a')
c_ops2 = (a, QobjEvo(a', coef1))

@inferred liouvillian(H_td, c_ops1)
@inferred liouvillian(H_td, c_ops2)
Expand Down
5 changes: 2 additions & 3 deletions test/core-test/steady_state.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -64,9 +64,8 @@
H_td = (H, (H_t, coeff))

sol_me = mesolve(H_td, psi0, t_l, c_ops, e_ops = e_ops, progress_bar = Val(false))
ρ_ss1 = steadystate_fourier(H, -1im * 0.5 * H_t, 1im * 0.5 * H_t, 1, c_ops, solver = SteadyStateLinearSolver())[1]
ρ_ss2 =
steadystate_fourier(H, -1im * 0.5 * H_t, 1im * 0.5 * H_t, 1, c_ops, solver = SSFloquetEffectiveLiouvillian())
ρ_ss1 = steadystate_fourier(H, 0.5 * H_t, 0.5 * H_t, 1, c_ops, solver = SteadyStateLinearSolver())[1]
ρ_ss2 = steadystate_fourier(H, 0.5 * H_t, 0.5 * H_t, 1, c_ops, solver = SSFloquetEffectiveLiouvillian())

@test abs(sum(sol_me.expect[1, (end-100):end]) / 101 - expect(e_ops[1], ρ_ss1)) < 1e-3
@test abs(sum(sol_me.expect[1, (end-100):end]) / 101 - expect(e_ops[1], ρ_ss2)) < 1e-3
Expand Down
Loading