Skip to content

Commit f2acb28

Browse files
akirakyleakirakyle
committed
More Superoperator WIP
1 parent 69c7e43 commit f2acb28

File tree

3 files changed

+83
-87
lines changed

3 files changed

+83
-87
lines changed

‎src/QuantumOpticsBase.jl‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -41,6 +41,7 @@ export Basis, GenericBasis, CompositeBasis, basis, basis_l, basis_r,
4141
#superoperators
4242
KetBraBasis, ChoiBasis, super, choi, kraus, vec,
4343
spre, spost, sprepost, liouvillian, identitysuperoperator,
44+
SuperOperatorType, DenseSuperOpType, SparseSuperOpType,
4445
#fock
4546
FockBasis, number, destroy, create,
4647
fockstate, coherentstate, coherentstate!,

‎src/superoperators.jl‎

Lines changed: 30 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -1,9 +1,25 @@
11
import QuantumInterface: KetBraBasis, ChoiBasis
22

3+
const SuperKetType = Ket{<:KetBraBasis}
4+
5+
const DenseSuperOpPureType{BL,BR} = Operator{<:KetBraBasis,<:KetBraBasis,<:Matrix}
6+
const DenseSuperOpAdjType{BL,BR} = Operator{<:KetBraBasis,<:KetBraBasis,<:Adjoint{<:Number,<:Matrix}}
7+
const DenseSuperOpType{BL,BR} = Union{DenseOpPureType{<:KetBraBasis,<:KetBraBasis},DenseOpAdjType{<:KetBraBasis,BR}}
8+
const SparseSuperOpPureType{BL,BR} = Operator{<:KetBraBasis,<:KetBraBasis,<:SparseMatrixCSC}
9+
const SparseSuperOpAdjType{BL,BR} = Operator{<:KetBraBasis,<:KetBraBasis,<:Adjoint{<:Number,<:SparseMatrixCSC}}
10+
const SparseSuperOpType{BL,BR} = Union{SparseOpPureType{<:KetBraBasis,<:KetBraBasis},SparseOpAdjType{<:KetBraBasis,BR}}
11+
312
const SuperOperatorType = Operator{<:KetBraBasis,<:KetBraBasis}
4-
const ChoiStateType = Operator{CompositeBasis{<:Integer,<:ChoiBasis},CompositeBasis{<:Integer,<:ChoiBasis}}
13+
14+
#const ChoiBasisType = Union{CompositeBasis{ChoiBasis,T} where {T}, CompositeBasis{ChoiBasis{B},T} where {B,T}}
15+
#const ChoiStateType = Operator{ChoiBasisType,ChoiBasisType}
16+
const ChoiStateType = Union{Operator{CompositeBasis{ChoiBasis,S},CompositeBasis{ChoiBasis,T}} where {S,T}, Operator{CompositeBasis{ChoiBasis{BL},S},CompositeBasis{ChoiBasis{BR},T}} where {BL,BR,S,T}}
517

618
vec(op::Operator) = Ket(KetBraBasis(basis_l(op), basis_r(op)), reshape(op.data, length(op.data)))
19+
function unvec(k::SuperKetType)
20+
bl, br = basis_l(basis(k)), basis_r(basis(k))
21+
return Operator(bl, br, reshape(k.data, length(bl), length(br)))
22+
end
723

824
function spre(op::Operator)
925
multiplicable(op, op) || throw(ArgumentError("It's not clear what spre of a non-square operator should be. See issue #113"))
@@ -12,10 +28,10 @@ end
1228

1329
function spost(op::Operator)
1430
multiplicable(op, op) || throw(ArgumentError("It's not clear what spost of a non-square operator should be. See issue #113"))
15-
SuperOperator(KetBraBasis(basis_r(op)), (op.basis_l, op.basis_l), kron(permutedims(op.data), identityoperator(op).data))
31+
Operator(KetBraBasis(basis_r(op)), KetBraBasis(basis_l(op)), kron(permutedims(op.data), identityoperator(op).data))
1632
end
1733

18-
sprepost(A::Operator, B::Operator) = SuperOperator((A.basis_l, B.basis_r), (A.basis_r, B.basis_l), kron(permutedims(B.data), A.data))
34+
sprepost(A::Operator, B::Operator) = Operator(KetBraBasis(A.basis_l, B.basis_r), KetBraBasis(A.basis_r, B.basis_l), kron(permutedims(B.data), A.data))
1935

2036
# reshape swaps within systems due to colum major ordering
2137
# https://docs.qojulia.org/quantumobjects/operators/#tensor_order
@@ -51,15 +67,24 @@ end
5167

5268
function super(op::ChoiStateType)
5369
bl, br = basis_l(op), basis_r(op)
54-
(l1, l2), (r1, r2) = (basis_l(bl), basis_r(bl)), (basis_l(br), basis_r(br))
70+
(l1, l2), (r1, r2) = (basis(bl[1]), basis(bl[2])), (basis(br[1]), basis(br[2]))
5571
(l1, l2), (r1, r2), data = _super_choi((l1, l2), (r1, r2), op.data)
5672
return Operator(KetBraBasis(l1,l2), KetBraBasis(r1,r2), data)
5773
end
5874

5975
dagger(a::ChoiStateType) = choi(dagger(super(a)))
6076

61-
*(a::SuperOperatorType, b::Operator) = a*vec(b)
77+
*(a::SuperOperatorType, b::SuperOperatorType) = (check_multiplicable(a,b); Operator(a.basis_l, b.basis_r, a.data*b.data))
78+
*(a::SuperOperatorType, b::Operator) = unvec(a*vec(b))
6279
*(a::ChoiStateType, b::SuperOperatorType) = super(a)*b
6380
*(a::SuperOperatorType, b::ChoiStateType) = a*super(b)
6481
*(a::ChoiStateType, b::ChoiStateType) = choi(super(a)*super(b))
6582
*(a::ChoiStateType, b::Operator) = super(a)*b
83+
84+
identitysuperoperator(b::Basis) = Operator(KetBraBasis(b,b), KetBraBasis(b,b), Eye{ComplexF64}(length(b)^2))
85+
86+
identitysuperoperator(op::DenseSuperOpType) =
87+
Operator(op.basis_l, op.basis_r, Matrix(one(eltype(op.data))I, size(op.data)))
88+
89+
identitysuperoperator(op::SparseSuperOpType) =
90+
Operator(op.basis_l, op.basis_r, sparse(one(eltype(op.data))I, size(op.data)))

‎test/test_superoperators.jl‎

Lines changed: 52 additions & 82 deletions
Original file line numberDiff line numberDiff line change
@@ -1,38 +1,9 @@
11
using Test
22
using QuantumOpticsBase
33
using SparseArrays, LinearAlgebra
4-
import QuantumOpticsBase: ChoiState # Remove when ChoiState is publicly exported
54

65
@testset "superoperators" begin
76

8-
# Test creation
9-
b = GenericBasis(3)
10-
@test_throws DimensionMismatch DenseSuperOperator((b, b), (b, b), zeros(ComplexF64, 3, 3))
11-
@test_throws DimensionMismatch SparseSuperOperator((b, b), (b, b), spzeros(ComplexF64, 3, 3))
12-
@test_throws DimensionMismatch ChoiState((b, b), (b, b), zeros(ComplexF64, 3, 3))
13-
14-
# Test copy, sparse and dense
15-
b1 = GenericBasis(2)
16-
b2 = GenericBasis(7)
17-
b3 = GenericBasis(5)
18-
b4 = GenericBasis(3)
19-
20-
s = DenseSuperOperator((b1, b2), (b3, b4))
21-
s_ = dense(s)
22-
s_.data[1,1] = 1
23-
@test s.data[1,1] == 0
24-
s_sparse = sparse(s_)
25-
@test isa(s_sparse, SparseSuperOpType)
26-
@test s_sparse.data[1,1] == 1
27-
28-
s = SparseSuperOperator((b1, b2), (b3, b4))
29-
s_ = sparse(s)
30-
s_.data[1,1] = 1
31-
@test s.data[1,1] == 0
32-
s_dense = dense(s_)
33-
@test isa(s_dense, DenseSuperOpType)
34-
@test s_dense.data[1,1] == 1
35-
367
# Test length
378
b1 = GenericBasis(3)
389
op = DenseOperator(b1, b1)
@@ -51,76 +22,76 @@ b2 = GenericBasis(5)
5122
b3 = GenericBasis(5)
5223
b4 = GenericBasis(3)
5324

54-
d1 = DenseSuperOperator((b1, b2), (b3, b4))
55-
d2 = DenseSuperOperator((b3, b4), (b1, b2))
56-
s1 = SparseSuperOperator((b1, b2), (b3, b4))
57-
s2 = SparseSuperOperator((b3, b4), (b1, b2))
25+
d1 = DenseOperator(KetBraBasis(b1, b2), KetBraBasis(b3, b4))
26+
d2 = DenseOperator(KetBraBasis(b3, b4), KetBraBasis(b1, b2))
27+
s1 = SparseOperator(KetBraBasis(b1, b2), KetBraBasis(b3, b4))
28+
s2 = SparseOperator(KetBraBasis(b3, b4), KetBraBasis(b1, b2))
5829

5930
x = d1*d2
6031
@test isa(x, DenseSuperOpType)
61-
@test x.basis_l == x.basis_r == (b1, b2)
32+
@test x.basis_l == x.basis_r == KetBraBasis(b1, b2)
6233

6334
x = s1*s2
6435
@test isa(x, SparseSuperOpType)
65-
@test x.basis_l == x.basis_r == (b1, b2)
36+
@test x.basis_l == x.basis_r == KetBraBasis(b1, b2)
6637

6738
x = d1*s2
6839
@test isa(x, DenseSuperOpType)
69-
@test x.basis_l == x.basis_r == (b1, b2)
40+
@test x.basis_l == x.basis_r == KetBraBasis(b1, b2)
7041

7142
x = s1*d2
7243
@test isa(x, DenseSuperOpType)
73-
@test x.basis_l == x.basis_r == (b1, b2)
44+
@test x.basis_l == x.basis_r == KetBraBasis(b1, b2)
7445

7546
x = d1*3
7647
@test isa(x, DenseSuperOpType)
77-
@test x.basis_l == (b1, b2)
78-
@test x.basis_r == (b3, b4)
48+
@test x.basis_l == KetBraBasis(b1, b2)
49+
@test x.basis_r == KetBraBasis(b3, b4)
7950

8051
x = 2.5*s1
8152
@test isa(x, SparseSuperOpType)
82-
@test x.basis_l == (b1, b2)
83-
@test x.basis_r == (b3, b4)
53+
@test x.basis_l == KetBraBasis(b1, b2)
54+
@test x.basis_r == KetBraBasis(b3, b4)
8455

8556
x = d1 + d1
8657
@test isa(x, DenseSuperOpType)
87-
@test x.basis_l == (b1, b2)
88-
@test x.basis_r == (b3, b4)
58+
@test x.basis_l == KetBraBasis(b1, b2)
59+
@test x.basis_r == KetBraBasis(b3, b4)
8960

9061
x = s1 + s1
9162
@test isa(x, SparseSuperOpType)
92-
@test x.basis_l == (b1, b2)
93-
@test x.basis_r == (b3, b4)
63+
@test x.basis_l == KetBraBasis(b1, b2)
64+
@test x.basis_r == KetBraBasis(b3, b4)
9465

9566
x = d1 + s1
9667
@test isa(x, DenseSuperOpType)
97-
@test x.basis_l == (b1, b2)
98-
@test x.basis_r == (b3, b4)
68+
@test x.basis_l == KetBraBasis(b1, b2)
69+
@test x.basis_r == KetBraBasis(b3, b4)
9970

10071
x = d1 - d1
10172
@test isa(x, DenseSuperOpType)
102-
@test x.basis_l == (b1, b2)
103-
@test x.basis_r == (b3, b4)
73+
@test x.basis_l == KetBraBasis(b1, b2)
74+
@test x.basis_r == KetBraBasis(b3, b4)
10475

10576
x = s1 - s1
10677
@test isa(x, SparseSuperOpType)
107-
@test x.basis_l == (b1, b2)
108-
@test x.basis_r == (b3, b4)
78+
@test x.basis_l == KetBraBasis(b1, b2)
79+
@test x.basis_r == KetBraBasis(b3, b4)
10980

11081
x = d1 - s1
11182
@test isa(x, DenseSuperOpType)
112-
@test x.basis_l == (b1, b2)
113-
@test x.basis_r == (b3, b4)
83+
@test x.basis_l == KetBraBasis(b1, b2)
84+
@test x.basis_r == KetBraBasis(b3, b4)
11485

11586
x = -d1
11687
@test isa(x, DenseSuperOpType)
117-
@test x.basis_l == (b1, b2)
118-
@test x.basis_r == (b3, b4)
88+
@test x.basis_l == KetBraBasis(b1, b2)
89+
@test x.basis_r == KetBraBasis(b3, b4)
11990

12091
x = -s1
12192
@test isa(x, SparseSuperOpType)
122-
@test x.basis_l == (b1, b2)
123-
@test x.basis_r == (b3, b4)
93+
@test x.basis_l == KetBraBasis(b1, b2)
94+
@test x.basis_r == KetBraBasis(b3, b4)
12495

12596
# TODO: Clean-up this part
12697
ωc = 1.2
@@ -131,7 +102,6 @@ g = 1.0
131102

132103
T = Float64[0.,1.]
133104

134-
135105
fockbasis = FockBasis(7)
136106
spinbasis = SpinBasis(1//2)
137107
basis = tensor(spinbasis, fockbasis)
@@ -155,45 +125,45 @@ J = [Ja, Jc]
155125
ρ₀ = dm(Ψ₀)
156126

157127
@test identitysuperoperator(spinbasis)*sx == sx
158-
@test ChoiState(identitysuperoperator(spinbasis))*sx == sx
128+
@test choi(identitysuperoperator(spinbasis))*sx == sx
159129
@test identitysuperoperator(sparse(spre(sx)))*sx == sparse(sx)
160-
@test ChoiState(identitysuperoperator(sparse(spre(sx))))*sx == sparse(sx)
161-
@test sparse(ChoiState(identitysuperoperator(spre(sx))))*sx == sparse(sx)
130+
@test choi(identitysuperoperator(sparse(spre(sx))))*sx == sparse(sx)
131+
@test sparse(choi(identitysuperoperator(spre(sx))))*sx == sparse(sx)
162132
@test identitysuperoperator(dense(spre(sx)))*sx == dense(sx)
163-
@test ChoiState(identitysuperoperator(dense(spre(sx))))*sx == dense(sx)
164-
@test dense(ChoiState(identitysuperoperator(spre(sx))))*sx == dense(sx)
133+
@test choi(identitysuperoperator(dense(spre(sx))))*sx == dense(sx)
134+
@test dense(choi(identitysuperoperator(spre(sx))))*sx == dense(sx)
165135

166136
op1 = DenseOperator(spinbasis, [1.2+0.3im 0.7+1.2im;0.3+0.1im 0.8+3.2im])
167137
op2 = DenseOperator(spinbasis, [0.2+0.1im 0.1+2.3im; 0.8+4.0im 0.3+1.4im])
168138
@test tracedistance(spre(op1)*op2, op1*op2) < 1e-12
169-
@test tracedistance(ChoiState(spre(op1))*op2, op1*op2) < 1e-12
139+
@test tracedistance(choi(spre(op1))*op2, op1*op2) < 1e-12
170140
@test tracedistance(spost(op1)*op2, op2*op1) < 1e-12
171-
@test tracedistance(ChoiState(spost(op1))*op2, op2*op1) < 1e-12
141+
@test tracedistance(choi(spost(op1))*op2, op2*op1) < 1e-12
172142

173143
@test spre(sparse(op1))*op2 == op1*op2
174-
@test ChoiState(spre(sparse(op1)))*op2 == op1*op2
144+
@test choi(spre(sparse(op1)))*op2 == op1*op2
175145
@test spost(sparse(op1))*op2 == op2*op1
176-
@test ChoiState(spost(sparse(op1)))*op2 == op2*op1
146+
@test choi(spost(sparse(op1)))*op2 == op2*op1
177147
@test spre(sparse(dagger(op1)))*op2 == dagger(op1)*op2
178-
@test ChoiState(spre(sparse(dagger(op1))))*op2 == dagger(op1)*op2
148+
@test choi(spre(sparse(dagger(op1))))*op2 == dagger(op1)*op2
179149
@test spre(dense(dagger(op1)))*op2 ≈ dagger(op1)*op2
180-
@test ChoiState(spre(dense(dagger(op1))))*op2 ≈ dagger(op1)*op2
150+
@test choi(spre(dense(dagger(op1))))*op2 ≈ dagger(op1)*op2
181151
@test sprepost(sparse(op1), op1)*op2 ≈ op1*op2*op1
182-
@test ChoiState(sprepost(sparse(op1), op1))*op2 ≈ op1*op2*op1
152+
@test choi(sprepost(sparse(op1), op1))*op2 ≈ op1*op2*op1
183153

184154
@test spre(sparse(op1))*sparse(op2) == sparse(op1*op2)
185-
@test ChoiState(spre(sparse(op1)))*sparse(op2) == sparse(op1*op2)
155+
@test choi(spre(sparse(op1)))*sparse(op2) == sparse(op1*op2)
186156
@test spost(sparse(op1))*sparse(op2) == sparse(op2*op1)
187-
@test ChoiState(spost(sparse(op1)))*sparse(op2) == sparse(op2*op1)
157+
@test choi(spost(sparse(op1)))*sparse(op2) == sparse(op2*op1)
188158
@test sprepost(sparse(op1), sparse(op1))*sparse(op2) ≈ sparse(op1*op2*op1)
189-
@test ChoiState(sprepost(sparse(op1), sparse(op1)))*sparse(op2) ≈ sparse(op1*op2*op1)
159+
@test choi(sprepost(sparse(op1), sparse(op1)))*sparse(op2) ≈ sparse(op1*op2*op1)
190160

191161
@test sprepost(op1, op2) ≈ spre(op1)*spost(op2)
192162
b1 = FockBasis(1)
193163
b2 = FockBasis(5)
194164
op = fockstate(b1, 0) ⊗ dagger(fockstate(b2, 0))
195165
@test sprepost(dagger(op), op)*dm(fockstate(b1, 0)) == dm(fockstate(b2, 0))
196-
@test ChoiState(sprepost(dagger(op), op))*dm(fockstate(b1, 0)) == dm(fockstate(b2, 0))
166+
@test choi(sprepost(dagger(op), op))*dm(fockstate(b1, 0)) == dm(fockstate(b2, 0))
197167
@test_throws ArgumentError spre(op)
198168
@test_throws ArgumentError spost(op)
199169

@@ -203,17 +173,17 @@ for j=J
203173
ρ .+= j*ρ₀*dagger(j) - 0.5*dagger(j)*j*ρ₀ - 0.5*ρ₀*dagger(j)*j
204174
end
205175
@test tracedistance(L*ρ₀, ρ) < 1e-10
206-
@test tracedistance(ChoiState(L)*ρ₀, ρ) < 1e-10
176+
@test tracedistance(choi(L)*ρ₀, ρ) < 1e-10
207177

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

212182
@test dense(spre(op1)) == spre(op1)
213-
@test dense(ChoiState(spre(op1))) == ChoiState(spre(op1))
183+
@test dense(choi(spre(op1))) == choi(spre(op1))
214184

215185
@test L/2.0 == 0.5*L == L*0.5
216-
@test -L == SparseSuperOperator(L.basis_l, L.basis_r, -L.data)
186+
@test -L == SparseOperator(L.basis_l, L.basis_r, -L.data)
217187

218188
@test_throws AssertionError liouvillian(H, J; rates=zeros(4, 4))
219189

@@ -256,12 +226,12 @@ o_l = fockstate(b_fock, 2)
256226
encoder_kraus = z_l ⊗ dagger(spinup(b_logical)) + o_l ⊗ dagger(spindown(b_logical))
257227
encoder_sup = sprepost(encoder_kraus, dagger(encoder_kraus))
258228
decoder_sup = sprepost(dagger(encoder_kraus), encoder_kraus)
259-
@test SuperOperator(ChoiState(encoder_sup)).data == encoder_sup.data
229+
@test super(choi(encoder_sup)).data == encoder_sup.data
260230
@test decoder_sup == dagger(encoder_sup)
261-
@test ChoiState(decoder_sup) == dagger(ChoiState(encoder_sup))
231+
@test choi(decoder_sup) == dagger(choi(encoder_sup))
262232
@test decoder_sup*encoder_sup ≈ dense(identitysuperoperator(b_logical))
263-
@test decoder_sup*ChoiState(encoder_sup) ≈ dense(identitysuperoperator(b_logical))
264-
@test ChoiState(decoder_sup)*encoder_sup ≈ dense(identitysuperoperator(b_logical))
265-
@test SuperOperator(ChoiState(decoder_sup)*ChoiState(encoder_sup)) ≈ dense(identitysuperoperator(b_logical))
233+
@test decoder_sup*choi(encoder_sup) ≈ dense(identitysuperoperator(b_logical))
234+
@test choi(decoder_sup)*encoder_sup ≈ dense(identitysuperoperator(b_logical))
235+
@test super(choi(decoder_sup)*choi(encoder_sup)) ≈ dense(identitysuperoperator(b_logical))
266236

267237
end # testset

0 commit comments

Comments
 (0)