Add new symbolic bosonic operators and states

src/QSymbolicsBase/QSymbolicsBase.jl

@@ -43,8 +43,9 @@ export SymQObj,QObj,
        MixedState,IdentityOp,
        SApplyKet,SApplyBra,SMulOperator,SSuperOpApply,SCommutator,SAnticommutator,SBraKet,SOuterKetBra,
        HGate,XGate,YGate,ZGate,CPHASEGate,CNOTGate,
-       XBasisState,YBasisState,ZBasisState,FockState,CoherentState,SqueezedState,
+       XBasisState,YBasisState,ZBasisState,FockState,CoherentState,SqueezedState,EPRState,ThermalState,
        NumberOp,CreateOp,DestroyOp,PhaseShiftOp,DisplaceOp,SqueezeOp,
+       TwoSqueezeOp,BeamSplitterOp,
        XCXGate,XCYGate,XCZGate,YCXGate,YCYGate,YCZGate,ZCXGate,ZCYGate,ZCZGate,
        qsimplify,qsimplify_pauli,qsimplify_commutator,qsimplify_anticommutator,qsimplify_fock,
        qexpand,

src/QSymbolicsBase/predefined_fock.jl

@@ -26,6 +26,21 @@
 SqueezedState(z::Number) = SqueezedState(z, inf_fock_basis)
 symbollabel(x::SqueezedState) = "0,$(x.z)"
 
+"""Two-mode squeezed vacuum state, or EPR state, in defined Fock basis."""
+@withmetadata struct EPRState <: SpecialKet
+    z::Number
+    basis::CompositeBasis
+    function EPRState(z::Number, basis::CompositeBasis)
+        bases = basis.bases
+        length(bases) == 2 && all(x -> x isa FockBasis, bases) || 
+                throw(ArgumentError(lazy"The underlying basis for an EPR state must be a
+                                        tensor product of two bases for single-mode quantum systems."))
+        return new(z, basis)
+    end
+end
+EPRState(z::Number) = EPRState(z, inf_fock_basis^2)
+symbollabel(x::EPRState) = "0,$(x.z)"
+
 const inf_fock_basis = FockBasis(Inf,0.)
 """Vacuum basis state of n"""
 const vac = const F₀ = const F0 = FockState(0)
@@ -37,9 +52,13 @@
 ##
 
 abstract type AbstractSingleBosonOp <: Symbolic{AbstractOperator} end
+abstract type AbstractTwoBosonOp <: Symbolic{AbstractOperator} end
 abstract type AbstractSingleBosonGate <: AbstractSingleBosonOp end # TODO maybe an IsUnitaryTrait is a better choice
+abstract type AbstractTwoBosonGate <: AbstractTwoBosonOp end
 isexpr(::AbstractSingleBosonGate) = false
 basis(x::AbstractSingleBosonOp) = inf_fock_basis
+isexpr(::AbstractTwoBosonGate) = false
+basis(x::AbstractTwoBosonOp) = inf_fock_basis^2
 
 """Number operator.
 
@@ -161,4 +180,42 @@
     basis::FockBasis
 end
 SqueezeOp(z::Number) = SqueezeOp(z, inf_fock_basis)
-symbollabel(x::SqueezeOp) = "S($(x.z))"
+symbollabel(x::SqueezeOp) = "S($(x.z))"
+
+"""Thermal bosonic state in defined Fock basis."""
+@withmetadata struct ThermalState <: AbstractSingleBosonOp
+    photons::Int
+    basis::FockBasis
+end
+ThermalState(photons::Int) = ThermalState(photons, inf_fock_basis)
+symbollabel(x::ThermalState) = "ρₜₕ($(x.photons))"
+
+"""Two-mode squeezing operator in defined Fock basis."""
+@withmetadata struct TwoSqueezeOp <: AbstractTwoBosonOp
+    z::Number
+    basis::CompositeBasis
+    function TwoSqueezeOp(z::Number, basis::CompositeBasis)
+        bases = basis.bases
+        length(bases) == 2 && all(x -> x isa FockBasis, bases) || 
+                throw(ArgumentError(lazy"The underlying basis for a two-mode squeeze operator must be a
+                                        tensor product of two bases for single-mode quantum systems."))
+        return new(z, basis)
+    end
+end
+TwoSqueezeOp(z::Number) = TwoSqueezeOp(z, inf_fock_basis^2)
+symbollabel(x::TwoSqueezeOp) = "S₂($(x.z))"
+
+"""Two-mode beamsplitter operator in defined Fock basis."""
+@withmetadata struct BeamSplitterOp <: AbstractTwoBosonOp
+    transmit::Number
+    basis::CompositeBasis
+    function BeamSplitterOp(transmit::Number, basis::CompositeBasis)
+        bases = basis.bases
+        length(bases) == 2 && all(x -> x isa FockBasis, bases) || 
+                throw(ArgumentError(lazy"The underlying basis for a beam splitter operator must be a
+                                        tensor product of two bases for single-mode quantum systems."))
+        return new(transmit, basis)
+    end
+end
+BeamSplitterOp(transmit::Number) = BeamSplitterOp(transmit, inf_fock_basis^2)
+symbollabel(x::BeamSplitterOp) = "B($(x.transmit))"

src/QSymbolicsBase/rules.jl

@@ -11,7 +11,7 @@ function hasscalings(xs)
     end
 end
 _isa(T) = x->isa(x,T)
-_isequal(obj) = x->(x==obj)
+_isequal(obj) = x->isequal(x, obj)
 _vecisa(T) = x->all(_isa(T), x)
 
 ##
@@ -106,7 +106,8 @@ RULES_FOCK = [
     @rule(dagger(~o1::_isa(DisplaceOp)) * ~o2::_isa(DestroyOp) * ~o1 => (~o2) + (~o1).alpha*IdentityOp((~o2).basis)),
     @rule(dagger(~o1::_isa(DisplaceOp)) * ~o2::_isa(CreateOp) * ~o1 => (~o2) + conj((~o1).alpha)*IdentityOp((~o2).basis)),
     @rule(~o::_isa(DisplaceOp) * ~k::((x->(isa(x,FockState) && x.idx == 0))) => CoherentState((~o).alpha)),
-    @rule(~o::_isa(SqueezeOp) * ~k::_isequal(vac) => SqueezedState((~o).z, (~o).basis))
+    @rule(~o::_isa(SqueezeOp) * ~k::_isequal(vac) => SqueezedState((~o).z, (~o).basis)),
+    @rule(~o::_isa(TwoSqueezeOp) * ~k::_isequal(vac ⊗ vac) => EPRState((~o).z, (~o).basis))
 ]
 
 RULES_SIMPLIFY = [RULES_PAULI; RULES_COMMUTATOR; RULES_ANTICOMMUTATOR; RULES_FOCK]

test/test_fock.jl

@@ -7,7 +7,9 @@
     phase2 = PhaseShiftOp(pi)
     displace = DisplaceOp(im)
     squeezeop = SqueezeOp(pi)
+    twosqueezeop = TwoSqueezeOp(pi)
     sstate = SqueezedState(pi)
+    tsstate = EPRState(pi)
 
     @testset "ladder and number operators" begin
         @test isequal(qsimplify(Destroy*vac, rewriter=qsimplify_fock), SZeroKet())
@@ -28,5 +30,6 @@
 
     @testset "Squeeze operators" begin
         @test isequal(qsimplify(squeezeop*vac, rewriter=qsimplify_fock), sstate)
+        @test isequal(qsimplify(twosqueezeop*(vac ⊗ vac), rewriter=qsimplify_fock), tsstate)
     end
 end

test/test_throws.jl

@@ -4,6 +4,7 @@
 
     @op A SpinBasis(1//2) ⊗ SpinBasis(1//2); @op B; @op C;
     @ket k; @bra b; @ket l SpinBasis(1//2) ⊗ SpinBasis(1//2);
+    fb = FockBasis(10)
 
     @test_throws IncompatibleBases A*B
     @test_throws IncompatibleBases commutator(A, B)
@@ -14,4 +15,7 @@
     @test_throws IncompatibleBases b*l
     @test_throws ArgumentError ptrace(B, 2)
     @test_throws ArgumentError ptrace(B+C, 2)
+    @test_throws ArgumentError EPRState(pi, fb^3)
+    @test_throws ArgumentError TwoSqueezeOp(pi, fb^3)
+    @test_throws ArgumentError BeamSplitterOp(pi, fb^3)
 end