Skip to content
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.

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

Create ChoiState type and conversions to/from SuperOperator #115

Merged
merged 4 commits into from
Nov 17, 2024
Merged
Show file tree
Hide file tree
Changes from 1 commit
Commits
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
4 changes: 2 additions & 2 deletions src/QuantumOpticsBase.jl
Original file line number Diff line number Diff line change
Expand Up @@ -36,8 +36,8 @@ export Basis, GenericBasis, CompositeBasis, basis,
current_time, time_shift, time_stretch, time_restrict, static_operator,
#superoperators
SuperOperator, DenseSuperOperator, DenseSuperOpType,
SparseSuperOperator, SparseSuperOpType, spre, spost, sprepost, liouvillian,
identitysuperoperator,
SparseSuperOperator, SparseSuperOpType, ChoiState,
Copy link
Collaborator

Choose a reason for hiding this comment

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

Let's keep it unexported until more of the related functionality is completed and it has been tested a bit more. I want to avoid a breaking release.

Copy link
Member Author

Choose a reason for hiding this comment

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

Sounds good to me, especially since #177 may rework things a bit!

spre, spost, sprepost, liouvillian, identitysuperoperator,
#fock
FockBasis, number, destroy, create,
fockstate, coherentstate, coherentstate!,
Expand Down
64 changes: 62 additions & 2 deletions src/superoperators.jl
Original file line number Diff line number Diff line change
Expand Up @@ -11,8 +11,9 @@
basis_r::B2
data::T
function SuperOperator{BL,BR,T}(basis_l::BL, basis_r::BR, data::T) where {BL,BR,T}
if length(basis_l[1])*length(basis_l[2]) != size(data, 1) ||
length(basis_r[1])*length(basis_r[2]) != size(data, 2)
if (length(basis_l) != 2 || length(basis_r) != 2 ||
length(basis_l[1])*length(basis_l[2]) != size(data, 1) ||
length(basis_r[1])*length(basis_r[2]) != size(data, 2))
throw(DimensionMismatch("Tried to assign data of size $(size(data)) to Hilbert spaces of sizes $(length.(basis_l)), $(length.(basis_r))"))
end
new(basis_l, basis_r, data)
Expand Down Expand Up @@ -316,3 +317,62 @@
}
throw(IncompatibleBases())
end

"""
ChoiState <: AbstractSuperOperator

Superoperator represented as a choi state.
"""
mutable struct ChoiState{B1,B2,T} <: AbstractSuperOperator{B1,B2}
basis_l::B1
basis_r::B2
data::T
function ChoiState{BL,BR,T}(basis_l::BL, basis_r::BR, data::T) where {BL,BR,T}
if (length(basis_l) != 2 || length(basis_r) != 2 ||
length(basis_l[1])*length(basis_l[2]) != size(data, 1) ||
length(basis_r[1])*length(basis_r[2]) != size(data, 2))
throw(DimensionMismatch("Tried to assign data of size $(size(data)) to Hilbert spaces of sizes $(length.(basis_l)), $(length.(basis_r))"))
end
new(basis_l, basis_r, data)
end
end
ChoiState{BL,BR}(b1::BL, b2::BR, data::T) where {BL,BR,T} = ChoiState{BL,BR,T}(b1, b2, data)

Check warning on line 339 in src/superoperators.jl

View check run for this annotation

Codecov / codecov/patch

src/superoperators.jl#L339

Added line #L339 was not covered by tests
Copy link
Collaborator

Choose a reason for hiding this comment

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

just to make sure everything is really fine, could you add a test for this constructor as well

Copy link
Member Author

Choose a reason for hiding this comment

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

Ah this is left over from my copy-paste of the SuperOperator constructors. I don't actually use this so I just went ahead and deleted it. It can easily be added back if it's needed at any point in the future.

ChoiState(b1::BL, b2::BR, data::T) where {BL,BR,T} = ChoiState{BL,BR,T}(b1, b2, data)

dense(a::ChoiState) = ChoiState(a.basis_l, a.basis_r, Matrix(a.data))
sparse(a::ChoiState) = ChoiState(a.basis_l, a.basis_r, sparse(a.data))
dagger(a::ChoiState) = ChoiState(dagger(SuperOperator(a)))
*(a::ChoiState, b::ChoiState) = ChoiState(SuperOperator(a)*SuperOperator(b))
*(a::ChoiState, b::Operator) = SuperOperator(a)*b
==(a::ChoiState, b::ChoiState) = (SuperOperator(a) == SuperOperator(b))

# reshape swaps within systems due to colum major ordering
# https://docs.qojulia.org/quantumobjects/operators/#tensor_order
function _super_choi((l1, l2), (r1, r2), data)
data = reshape(data, map(length, (l2, l1, r2, r1)))
(l1, l2), (r1, r2) = (r2, l2), (r1, l1)
data = permutedims(data, (1, 3, 2, 4))
data = reshape(data, map(length, (l1⊗l2, r1⊗r2)))
return (l1, l2), (r1, r2), data
end

function _super_choi((r2, l2), (r1, l1), data::SparseMatrixCSC)
data = _permutedims(data, map(length, (l2, r2, l1, r1)), (1, 3, 2, 4))
data = reshape(data, map(length, (l1⊗l2, r1⊗r2)))
# sparse(data) is necessary since reshape of a sparse array returns a
# ReshapedSparseArray which is not a subtype of AbstractArray and so
# _permutedims fails to acces the ".m" field
# https://github.com/qojulia/QuantumOpticsBase.jl/pull/83
# https://github.com/JuliaSparse/SparseArrays.jl/issues/24
# permutedims in SparseArrays.jl only implements perm (2,1) and so
# _permutedims should probably be upstreamed
# https://github.com/JuliaLang/julia/issues/26534
return (l1, l2), (r1, r2), sparse(data)
end

ChoiState(op::SuperOperator) = ChoiState(_super_choi(op.basis_l, op.basis_r, op.data)...)
SuperOperator(op::ChoiState) = SuperOperator(_super_choi(op.basis_l, op.basis_r, op.data)...)
Copy link
Collaborator

Choose a reason for hiding this comment

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

does it make sense to also add a Ket(::ChoiState) constructor, so that we can actually extract it as a state?

Copy link
Member Author

Choose a reason for hiding this comment

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

I don't think it makes sense since in general the ChoiState is mixed iff the corresponding channel is non-unitary, so it is akin to asking for a Ket(::Operator) method which has multiple ways in which it could make sense, e.g. closest pure state under a given operator norms.

Copy link
Collaborator

Choose a reason for hiding this comment

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

Ah, noted. Then maybe a dm(::ChoiState) would make sense, but no need to block the PR on this.

Copy link
Member Author

Choose a reason for hiding this comment

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

I have something like that added as part of #177. The tricky part is ensuring the input and ancilla bases make sense and can be appropriately manipulated!


*(a::ChoiState, b::SuperOperator) = SuperOperator(a)*b
*(a::SuperOperator, b::ChoiState) = a*SuperOperator(b)

35 changes: 35 additions & 0 deletions test/test_superoperators.jl
Original file line number Diff line number Diff line change
Expand Up @@ -8,6 +8,7 @@ using SparseArrays, LinearAlgebra
b = GenericBasis(3)
@test_throws DimensionMismatch DenseSuperOperator((b, b), (b, b), zeros(ComplexF64, 3, 3))
@test_throws DimensionMismatch SparseSuperOperator((b, b), (b, b), spzeros(ComplexF64, 3, 3))
@test_throws DimensionMismatch ChoiState((b, b), (b, b), zeros(ComplexF64, 3, 3))

# Test copy, sparse and dense
b1 = GenericBasis(2)
Expand Down Expand Up @@ -153,29 +154,45 @@ J = [Ja, Jc]
ρ₀ = dm(Ψ₀)

@test identitysuperoperator(spinbasis)*sx == sx
@test ChoiState(identitysuperoperator(spinbasis))*sx == sx
@test identitysuperoperator(sparse(spre(sx)))*sx == sparse(sx)
@test ChoiState(identitysuperoperator(sparse(spre(sx))))*sx == sparse(sx)
@test sparse(ChoiState(identitysuperoperator(spre(sx))))*sx == sparse(sx)
@test identitysuperoperator(dense(spre(sx)))*sx == dense(sx)
@test ChoiState(identitysuperoperator(dense(spre(sx))))*sx == dense(sx)
@test dense(ChoiState(identitysuperoperator(spre(sx))))*sx == dense(sx)

op1 = DenseOperator(spinbasis, [1.2+0.3im 0.7+1.2im;0.3+0.1im 0.8+3.2im])
op2 = DenseOperator(spinbasis, [0.2+0.1im 0.1+2.3im; 0.8+4.0im 0.3+1.4im])
@test tracedistance(spre(op1)*op2, op1*op2) < 1e-12
@test tracedistance(ChoiState(spre(op1))*op2, op1*op2) < 1e-12
@test tracedistance(spost(op1)*op2, op2*op1) < 1e-12
@test tracedistance(ChoiState(spost(op1))*op2, op2*op1) < 1e-12

@test spre(sparse(op1))*op2 == op1*op2
@test ChoiState(spre(sparse(op1)))*op2 == op1*op2
@test spost(sparse(op1))*op2 == op2*op1
@test ChoiState(spost(sparse(op1)))*op2 == op2*op1
@test spre(sparse(dagger(op1)))*op2 == dagger(op1)*op2
@test ChoiState(spre(sparse(dagger(op1))))*op2 == dagger(op1)*op2
@test spre(dense(dagger(op1)))*op2 ≈ dagger(op1)*op2
@test ChoiState(spre(dense(dagger(op1))))*op2 ≈ dagger(op1)*op2
@test sprepost(sparse(op1), op1)*op2 ≈ op1*op2*op1
@test ChoiState(sprepost(sparse(op1), op1))*op2 ≈ op1*op2*op1

@test spre(sparse(op1))*sparse(op2) == sparse(op1*op2)
@test ChoiState(spre(sparse(op1)))*sparse(op2) == sparse(op1*op2)
@test spost(sparse(op1))*sparse(op2) == sparse(op2*op1)
@test ChoiState(spost(sparse(op1)))*sparse(op2) == sparse(op2*op1)
@test sprepost(sparse(op1), sparse(op1))*sparse(op2) ≈ sparse(op1*op2*op1)
@test ChoiState(sprepost(sparse(op1), sparse(op1)))*sparse(op2) ≈ sparse(op1*op2*op1)

@test sprepost(op1, op2) ≈ spre(op1)*spost(op2)
b1 = FockBasis(1)
b2 = FockBasis(5)
op = fockstate(b1, 0) ⊗ dagger(fockstate(b2, 0))
@test sprepost(dagger(op), op)*dm(fockstate(b1, 0)) == dm(fockstate(b2, 0))
@test ChoiState(sprepost(dagger(op), op))*dm(fockstate(b1, 0)) == dm(fockstate(b2, 0))
@test_throws ArgumentError spre(op)
@test_throws ArgumentError spost(op)

Expand All @@ -185,12 +202,14 @@ for j=J
ρ .+= j*ρ₀*dagger(j) - 0.5*dagger(j)*j*ρ₀ - 0.5*ρ₀*dagger(j)*j
end
@test tracedistance(L*ρ₀, ρ) < 1e-10
@test tracedistance(ChoiState(L)*ρ₀, ρ) < 1e-10

# tout, ρt = timeevolution.master([0.,1.], ρ₀, H, J; reltol=1e-7)
# @test tracedistance(exp(dense(L))*ρ₀, ρt[end]) < 1e-6
# @test tracedistance(exp(sparse(L))*ρ₀, ρt[end]) < 1e-6

@test dense(spre(op1)) == spre(op1)
@test dense(ChoiState(spre(op1))) == ChoiState(spre(op1))

@test L/2.0 == 0.5*L == L*0.5
@test -L == SparseSuperOperator(L.basis_l, L.basis_r, -L.data)
Expand Down Expand Up @@ -228,4 +247,20 @@ N = exp(log(2) * sparse(L)) # 50% loss channel
@test (0.5 - real(tr(ρ^2))) < 1e-10 # one photon state becomes maximally mixed
@test tracedistance(ρ, normalize(dm(fockstate(b, 0)) + dm(fockstate(b, 1)))) < 1e-10

# Testing 0-2-4 binomial code encoder
Copy link
Collaborator

Choose a reason for hiding this comment

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

this is a really neat test!

b_logical = SpinBasis(1//2)
b_fock = FockBasis(5)
z_l = normalize(fockstate(b_fock, 0) + fockstate(b_fock, 4))
o_l = fockstate(b_fock, 2)
encoder_kraus = z_l ⊗ dagger(spinup(b_logical)) + o_l ⊗ dagger(spindown(b_logical))
encoder_sup = sprepost(encoder_kraus, dagger(encoder_kraus))
decoder_sup = sprepost(dagger(encoder_kraus), encoder_kraus)
@test SuperOperator(ChoiState(encoder_sup)).data == encoder_sup.data
@test decoder_sup == dagger(encoder_sup)
@test ChoiState(decoder_sup) == dagger(ChoiState(encoder_sup))
@test decoder_sup*encoder_sup ≈ dense(identitysuperoperator(b_logical))
@test decoder_sup*ChoiState(encoder_sup) ≈ dense(identitysuperoperator(b_logical))
@test ChoiState(decoder_sup)*encoder_sup ≈ dense(identitysuperoperator(b_logical))
@test SuperOperator(ChoiState(decoder_sup)*ChoiState(encoder_sup)) ≈ dense(identitysuperoperator(b_logical))

end # testset
Loading