noncliff: improvements to apply! and dictvaltype

src/nonclifford.jl

@@ -89,16 +89,70 @@ function _stabmixdestab(mixeddestab, d)
     p
 end
 
+"""
+Updates the generalized stabilizer `τ = (χ, B(S, D))` by applying a Clifford gate `G`,
+where `χ` represents the density matrix and `B(S, D)` denotes the stabilizer basis in
+which `χ` is expressed.
+
+```jldoctest
+julia> sm = GeneralizedStabilizer(S"-X")
+A mixture ∑ ϕᵢⱼ Pᵢ ρ Pⱼ† where ρ is
+𝒟ℯ𝓈𝓉𝒶𝒷
++ Z
+𝒮𝓉𝒶𝒷
+- X
+with ϕᵢⱼ | Pᵢ | Pⱼ:
+ 1.0+0.0im | + _ | + _
+
+julia> apply!(sm, CliffordOperator(tHadamard))
+A mixture ∑ ϕᵢⱼ Pᵢ ρ Pⱼ† where ρ is
+𝒟ℯ𝓈𝓉𝒶𝒷
++ X
+𝒮𝓉𝒶𝒷
+- Z
+with ϕᵢⱼ | Pᵢ | Pⱼ:
+ 1.0+0.0im | + _ | + _
+```
+
+"""
 function apply!(state::GeneralizedStabilizer, gate::AbstractCliffordOperator) # TODO conjugate also the destabs
     apply!(state.stab, gate)
     state
 end
 
 """$(TYPEDSIGNATURES)
 
-Expectation value for the [PauliOperator](@ref) observable given the [`GeneralizedStabilizer`](@ref) state `s`."""
+Expectation value for the [PauliOperator](@ref) observable given the [`GeneralizedStabilizer`](@ref) state `s`.
+
+```jldoctest
+julia> sm = GeneralizedStabilizer(S"-X")
+A mixture ∑ ϕᵢⱼ Pᵢ ρ Pⱼ† where ρ is
+𝒟ℯ𝓈𝓉𝒶𝒷
++ Z
+𝒮𝓉𝒶𝒷
+- X
+with ϕᵢⱼ | Pᵢ | Pⱼ:
+ 1.0+0.0im | + _ | + _
+
+julia> apply!(sm, pcT)
+A mixture ∑ ϕᵢⱼ Pᵢ ρ Pⱼ† where ρ is
+𝒟ℯ𝓈𝓉𝒶𝒷
++ Z
+𝒮𝓉𝒶𝒷
+- X
+with ϕᵢⱼ | Pᵢ | Pⱼ:
+ 0.0+0.353553im | + _ | + Z
+ 0.0-0.353553im | + Z | + _
+ 0.853553+0.0im | + _ | + _
+ 0.146447+0.0im | + Z | + Z
+
+julia> expect(P"-X", sm)
+0.7071067811865475 + 0.0im
+```
+
+"""
 function expect(p::PauliOperator, s::GeneralizedStabilizer) # TODO optimize
-    χ′ = zero(_dictvaltype(s.destabweights))
+    χ′ = zero(valtype(s.destabweights))
     phase, b, c = rowdecompose(p, s.stab)
     for ((dᵢ,dⱼ), χ) in s.destabweights
         _allthreesumtozero(dᵢ,dⱼ,b) || continue
@@ -115,10 +169,6 @@ function _allthreesumtozero(a,b,c) # TODO consider using bitpacking and SIMD xor
     true
 end
 
-function _dictvaltype(dict)
-    return eltype(dict).parameters[2] # TODO there must be a cleaner way to do this
-end
-
 function project!(sm::GeneralizedStabilizer, p::PauliOperator)
     eval = expect(p, sm)
     prob₁ = (real(eval)+1)/2
@@ -173,13 +223,13 @@ end
 
 nqubits(pc::PauliChannel) = nqubits(pc.paulis[1][1])
 
-function apply!(state::GeneralizedStabilizer, gate::PauliChannel)
+function apply!(state::GeneralizedStabilizer, gate::AbstractPauliChannel; prune_threshold::Float64=1e-14)
     dict = state.destabweights
     stab = state.stab
-    dtype = _dictvaltype(dict)
+    dtype = valtype(dict)
     tzero = zero(dtype)
     tone = one(dtype)
-    newdict = typeof(dict)(tzero) # TODO jeez, this is ugly
+    newdict = typeof(dict)(tzero)
     for ((dᵢ,dⱼ), χ) in dict # the state
         for ((Pₗ,Pᵣ), w) in zip(gate.paulis,gate.weights) # the channel
             phaseₗ, dₗ, dₗˢᵗᵃᵇ = rowdecompose(Pₗ,stab)
@@ -188,14 +238,10 @@ function apply!(state::GeneralizedStabilizer, gate::PauliChannel)
             dᵢ′ = dₗ .⊻ dᵢ
             dⱼ′ = dᵣ .⊻ dⱼ
             χ′ = χ * w * (-tone)^c * (im)^(-phaseₗ+phaseᵣ+4)
-            newdict[(dᵢ′,dⱼ′)] += χ′
-        end
-    end
-    for (k,v) in newdict # TODO is it safe to modify a dict while iterating over it?
-        if abs(v) < 1e-14 # TODO parameterize this pruning parameter
-            delete!(newdict, k)
+            newdict[(dᵢ′,dⱼ′)] = get!(newdict,(dᵢ′,dⱼ′),0)+χ′
         end
     end
+    filter!(x -> abs(x[2]) > prune_threshold, newdict)
     state.destabweights = newdict
     state
 end
@@ -302,7 +348,7 @@ end
 
 nqubits(pc::UnitaryPauliChannel) = nqubits(pc.paulis[1])
 
-apply!(state::GeneralizedStabilizer, gate::UnitaryPauliChannel) = apply!(state, gate.paulichannel)
+apply!(state::GeneralizedStabilizer, gate::UnitaryPauliChannel; prune_threshold::Float64=1e-14) = apply!(state, gate.paulichannel; prune_threshold=prune_threshold)
 
 ##
 # Predefined Pauli Channels

test/test_nonclifford_quantumoptics.jl

@@ -41,12 +41,12 @@ qo_tgate.data[2,2] = exp(im*pi/4)
             @test qo_bigtgate ≈ Operator(bigtgate)
 
             for step in 1:10
-            # apply!(ket, qo_bigtgate) TODO implement this API
+                # apply!(ket, qo_bigtgate) TODO implement this API
                 ket = qo_bigtgate*ket
                 apply!(genstab, bigtgate)
                 @test dm(ket) ≈ Operator(genstab)
                 @test isapprox(expect(qo_pauli, ket), expect(pauli, genstab); atol=1e-5)
-           end
+            end
         end
     end
-end
+end