Add identitysuperoperator and dagger for SuperOperator (#116)

* Add identitysuperoperator and dagger for SuperOperator

* import identitysuperoperator from QuantumInterface

* Add identitysuperoperator(b::Basis) and some tests

* whitelisth a pirate (identitysuperoperator)

---------

Co-authored-by: Stefan Krastanov <[email protected]>
Co-authored-by: Stefan Krastanov <[email protected]>

Author: akirakyle

src/QuantumOpticsBase.jl

@@ -3,7 +3,7 @@ module QuantumOpticsBase
 using SparseArrays, LinearAlgebra, LRUCache, Strided, UnsafeArrays, FillArrays
 import LinearAlgebra: mul!, rmul!
 
-import QuantumInterface: dagger, directsum, ⊕, dm, embed, expect, identityoperator,
+import QuantumInterface: dagger, directsum, ⊕, dm, embed, expect, identityoperator, identitysuperoperator,
         permutesystems, projector, ptrace, reduced, tensor, ⊗, variance, apply!, basis, AbstractSuperOperator
 
 # index helpers
@@ -34,6 +34,7 @@ export Basis, GenericBasis, CompositeBasis, basis,
         #superoperators
                 SuperOperator, DenseSuperOperator, DenseSuperOpType,
                 SparseSuperOperator, SparseSuperOpType, spre, spost, sprepost, liouvillian,
+                identitysuperoperator,
         #fock
                 FockBasis, number, destroy, create,
                 fockstate, coherentstate, coherentstate!,

src/superoperators.jl

@@ -97,6 +97,19 @@ end
 -(a::SuperOperator) = SuperOperator(a.basis_l, a.basis_r, -a.data)
 -(a::SuperOperator, b::SuperOperator) = throw(IncompatibleBases())
 
+identitysuperoperator(b::Basis) =
+    SuperOperator((b,b), (b,b), Eye{ComplexF64}(length(b)^2))
+
+identitysuperoperator(op::DenseSuperOpType) = 
+    SuperOperator(op.basis_l, op.basis_r, Matrix(one(eltype(op.data))I, size(op.data)))
+
+identitysuperoperator(op::SparseSuperOpType) = 
+    SuperOperator(op.basis_l, op.basis_r, sparse(one(eltype(op.data))I, size(op.data)))
+
+dagger(x::DenseSuperOpType) = SuperOperator(x.basis_r, x.basis_l, copy(adjoint(x.data)))
+dagger(x::SparseSuperOpType) = SuperOperator(x.basis_r, x.basis_l, sparse(adjoint(x.data)))
+
+
 """
     spre(op)
 

test/test_aqua.jl

@@ -9,6 +9,6 @@ using FillArrays
                   piracy=(broken=true,)
                   )
     # manual piracy check to exclude identityoperator
-    pirates = [pirate for pirate in Aqua.Piracy.hunt(QuantumOpticsBase) if pirate.name != :identityoperator]
+    pirates = [pirate for pirate in Aqua.Piracy.hunt(QuantumOpticsBase) if pirate.name ∉ [:identityoperator,:identitysuperoperator]]
     @test isempty(pirates)
 end # testset

test/test_superoperators.jl

@@ -120,7 +120,6 @@ x = -s1
 @test x.basis_l == (b1, b2)
 @test x.basis_r == (b3, b4)
 
-
 # TODO: Clean-up this part
 ωc = 1.2
 ωa = 0.9
@@ -153,6 +152,9 @@ J = [Ja, Jc]
 Ψ₀ = spinup(spinbasis) ⊗ fockstate(fockbasis, 5)
 ρ₀ = dm(Ψ₀)
 
+@test identitysuperoperator(spinbasis)*sx == sx
+@test identitysuperoperator(sparse(spre(sx)))*sx == sparse(sx)
+@test identitysuperoperator(dense(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])
@@ -161,6 +163,8 @@ op2 = DenseOperator(spinbasis, [0.2+0.1im 0.1+2.3im; 0.8+4.0im 0.3+1.4im])
 
 @test spre(sparse(op1))*op2 == op1*op2
 @test spost(sparse(op1))*op2 == op2*op1
+@test spre(sparse(dagger(op1)))*op2 == dagger(op1)*op2
+@test spre(dense(dagger(op1)))*op2 ≈ dagger(op1)*op2
 @test sprepost(sparse(op1), op1)*op2 ≈ op1*op2*op1
 
 @test spre(sparse(op1))*sparse(op2) == sparse(op1*op2)