define BellSwap, NoisyBellSwap and T1NoiseOp (#23)

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

Author: 3 people

.gitignore

@@ -4,3 +4,5 @@ LocalPreferences.toml
 */.*swp
 scratch/
 *.cov
+.history
+.vscode

CHANGELOG.md

@@ -0,0 +1,6 @@
+## v1.2.0 - 2024-11-04
+
+- Implementing convenience constructors for T1 and T2 noise (`T1NoiseOp` and `T2NoiseOp`)
+- Implementing a bilateral swap gate `BellSwap`
+
+## Older - before 2024-11-04 unrecorded

Project.toml

@@ -1,7 +1,7 @@
 name = "BPGates"
 uuid = "2c1a90cd-bec4-4415-ae35-245016909a8f"
-authors = ["Shu Ge <[email protected]>", "Stefan Krastanov <[email protected]>"]
-version = "1.1.0"
+authors = ["Shu Ge <[email protected]>", "Stefan Krastanov <[email protected]>", "BPGates contributors"]
+version = "1.2.0"
 
 [deps]
 QuantumClifford = "0525e862-1e90-11e9-3e4d-1b39d7109de1"
@@ -15,6 +15,6 @@ BPGatesQuantikzExt = "Quantikz"
 
 [compat]
 Quantikz = "1.3.1"
-QuantumClifford = "0.8, 0.9"
+QuantumClifford = "0.9"
 Random = "1"
 julia = "1.9"

src/BPGates.jl

@@ -9,17 +9,20 @@ export BellState,
     BellSinglePermutation, BellDoublePermutation, BellPauliPermutation,
     BellMeasure, bellmeasure!,
     BellGate, CNOTPerm, GoodSingleQubitPerm,
-    PauliNoiseOp, PauliNoiseBellGate, NoisyBellMeasure, NoisyBellMeasureNoisyReset
+    PauliNoiseOp, PauliNoiseBellGate, NoisyBellMeasure, NoisyBellMeasureNoisyReset,
+    BellSwap, NoisyBellSwap
 
-function int_to_bit(int,digits)
+const IT = Union{Int8,Int16,Int32,Int64,UInt8,UInt16,UInt32,UInt64}
+
+function int_to_bit(int::IT,digits)
     int = int - 1 # -1 so that we use julia indexing conventions
     Bool[int>>shift&0x1 for shift in 0:digits-1]
 end
-@inline function int_to_bit(int,::Val{2}) # faster if we know we need only two bits
+@inline function int_to_bit(int::IT,::Val{2}) # faster if we know we need only two bits
     int = int - 1 # -1 so that we use julia indexing conventions
     int & 0x1, int>>1 & 0x1
 end
-@inline function int_to_bit(int,::Val{4}) # faster if we know we need only two bits
+@inline function int_to_bit(int::IT,::Val{4}) # faster if we know we need only two bits
     int = int - 1 # -1 so that we use julia indexing conventions
     int & 0x1, int>>1 & 0x1, int>>2 & 0x1, int>>3 & 0x1
 end
@@ -367,6 +370,27 @@ function QuantumClifford.apply!(state::BellState, g::BellGate)
     return state
 end
 
+##############################
+# SWAP gate
+##############################
+
+"""SWAP gate"""
+struct BellSwap <: BellOp
+    idx1::Int
+    idx2::Int
+end
+
+function QuantumClifford.apply!(state::BellState, op::BellSwap)
+    phase = state.phases
+    @inbounds temp1 = phase[op.idx1*2-1]
+    @inbounds temp2 = phase[op.idx1*2]
+    @inbounds phase[op.idx1*2-1] = phase[op.idx2*2-1]
+    @inbounds phase[op.idx1*2] = phase[op.idx2*2]
+    @inbounds phase[op.idx2*2-1] = temp1
+    @inbounds phase[op.idx2*2] = temp2
+    return state
+end
+
 ##############################
 # Typically good operations
 ##############################
@@ -518,6 +542,112 @@ function QuantumClifford.apply!(state::BellState, g::PauliNoiseOp)
     return state
 end
 
+"""Simulates bilateral twirled T1 noise with per-qubit Kraus ops `|0⟩⟨0| + √(1-λ) |1⟩⟨1|` and `√λ |0⟩⟨1|`"""
+struct T1NoiseOp <: BellOp
+    idx::Int
+    λ₁::Float64
+end
+
+function QuantumClifford.apply!(state::BellState, g::T1NoiseOp)
+    phase = state.phases
+    input_state = bit_to_int(phase[g.idx*2-1],phase[g.idx*2]) # this is my initial state
+
+    r = rand()
+    λ₁ = g.λ₁
+
+    output_state = if input_state==1
+        if     r < 0.5*λ₁^2 - λ₁ + 1
+            1
+        elseif r < 0.5*λ₁^2 - λ₁ + 1  +  0.5*λ₁^2
+            2 # XXX
+        elseif r < 0.5*λ₁^2 - λ₁ + 1  +  0.5*λ₁^2 +  0.5*λ₁*(1-λ₁)
+            3 # XXX
+        else # r < 1 = 0.5*λ₁^2 - λ₁ + 1  +  0.5*λ₁*(1-λ₁)  +  0.5*λ₁^2  +   0.5*λ₁*(1-λ₁)
+            4
+        end
+    elseif input_state==2
+        if     r < 0.5*λ₁^2
+            1
+        elseif r < 0.5*λ₁^2 + 0.5*λ₁^2 -λ₁ + 1
+            2 # XXX
+        elseif r < 0.5*λ₁^2 + 0.5*λ₁^2 -λ₁ + 1 + 0.5*λ₁*(1-λ₁)
+            3 # XXX
+        else # r < 1 = 0.5*λ₁^2 + 0.5*λ₁^2 -λ₁ + 1 + 0.5*λ₁*(1-λ₁) + 0.5*λ₁*(1-λ₁)
+            4
+        end
+    elseif input_state==3
+        if     r < 0.5*λ₁
+            1
+        elseif r < 0.5*λ₁  +  0.5*λ₁
+            2 # XXX
+        elseif r < 0.5*λ₁  +  0.5*λ₁  +  (1-λ₁)
+            3 # XXX
+        else # r < 1 = 0.5*λ₁  +  0.5*λ₁  +  (1-λ₁)  +  0
+            4
+        end
+    else # input_state==4
+        if     r < 0.5*λ₁
+            1
+        elseif r < 0.5*λ₁ + 0.5*λ₁
+            2 # XXX
+        # elseif r < 0.5*λ₁ + 0.5*λ₁ + 0             # output_state 3 is never reached
+        #     3 # XXX
+        else # r < 1 = 0.5*λ₁ + 0.5*λ₁ + 0 + 1-λ₁
+            4
+        end
+    end
+
+    bit1, bit2 = int_to_bit(output_state, Val(2))
+    @inbounds phase[g.idx*2-1] = bit1
+    @inbounds phase[g.idx*2] = bit2
+    return state
+end
+
+"""Simulates bilateral T2 noise with per-qubit Kraus ops `√(1-λ/2) I` and `√(λ/2) Z`"""
+struct T2NoiseOp <: BellOp
+    idx::Int
+    λ₂::Float64
+end
+
+function QuantumClifford.apply!(state::BellState, g::T2NoiseOp)
+    phase = state.phases
+    input_state = bit_to_int(phase[g.idx*2-1],phase[g.idx*2]) # this is my initial state
+
+    r = rand()
+    λ₂ = g.λ₂
+
+    output_state = if input_state==1
+        if     r < 0.5*λ₂^2 - λ₂ + 1
+            1
+        else # r < 1 = 0.5*λ₂^2 - λ₂ + 1  +  0.5*λ₂*(2-λ₂)
+            2
+        end
+    elseif input_state==2
+        if     r < 0.5*λ₂^2 - λ₂ + 1
+            2
+        else # r < 1 = 0.5*λ₂^2 - λ₂ + 1  +  0.5*λ₂*(2-λ₂)
+            1
+        end
+    elseif input_state==3
+        if     r < 0.5*λ₂^2 - λ₂ + 1
+            3
+        else # r < 1 = 0.5*λ₂^2 - λ₂ + 1  +  0.5*λ₂*(2-λ₂)
+            4
+        end
+    else # input_state==4
+        if     r < 0.5*λ₂^2 - λ₂ + 1
+            4
+        else # r < 1 = 0.5*λ₂^2 - λ₂ + 1  +  0.5*λ₂*(2-λ₂)
+            3
+        end
+    end
+
+    bit1, bit2 = int_to_bit(output_state, Val(2))
+    @inbounds phase[g.idx*2-1] = bit1
+    @inbounds phase[g.idx*2] = bit2
+    return state
+end
+
 """A wrapper for [`BellMeasure`](@ref) that implements measurement noise."""
 struct NoisyBellMeasure <: BellOp # TODO make it work with the QuantumClifford noise ops
     m::BellMeasure
@@ -545,6 +675,11 @@ function QuantumClifford.applywstatus!(state::BellState, op::NoisyBellMeasureNoi
     state, cont ? continue_stat : failure_stat
 end
 
+function NoisyBellSwap(idx1,idx2,px,py,pz)
+    return PauliNoiseBellGate(BellSwap(idx1,idx2), px,py,pz)
+end
+
+
 ##############################
 # Random
 ##############################
@@ -703,12 +838,12 @@ function toBPpermutation1(circ) # this one works only on single pair mappings th
 end
 
 function BellState(s::Stabilizer)
-    r, c = size(s)
+    r, c = size(s)::Tuple{Int, Int}
     r==c || throw(ArgumentError("Conversion to `BellState` failed. The stabilizer state has to be square in order to be convertible to a `BellState`."))
-    s = canonicalize_rref!(copy(s))[1][end:-1:1]
+    s = canonicalize_rref!(copy(s))[1][end:-1:1]::Stabilizer
     bits = Bool[]
     for i in 1:r÷2
-        j = 2i-1
+        j = (2i-1)
         s[j  ,j] == s[j  ,j+1] == (true, false) && all(==((false,false)), (s[j  ,k] for k in 1:c if k!=j && k!=j+1)) || throw(ArgumentError(lazy"Conversion to `BellState` failed. Row $(j  ) of the stabilizer state has to be of the form `..._XX_...` for it to be a valid `BellState`"))
         s[j+1,j] == s[j+1,j+1] == (false, true) && all(==((false,false)), (s[j+1,k] for k in 1:c if k!=j && k!=j+1)) || throw(ArgumentError(lazy"Conversion to `BellState` failed. Row $(j+1) row of the stabilizer state has to be of the form `..._ZZ_...` for it to be a valid `BellState`"))
         push!(bits, s.tab.phases[j]÷2)

test/Project.toml

@@ -8,8 +8,9 @@ JET = "c3a54625-cd67-489e-a8e7-0a5a0ff4e31b"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 Quantikz = "b0d11df0-eea3-4d79-b4a5-421488cbf74b"
 QuantumClifford = "0525e862-1e90-11e9-3e4d-1b39d7109de1"
+QuantumOpticsBase = "4f57444f-1401-5e15-980d-4471b28d5678"
 Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
-SafeTestsets = "1bc83da4-3b8d-516f-aca4-4fe02f6d838f"
 StableRNGs = "860ef19b-820b-49d6-a774-d7a799459cd3"
 Static = "aedffcd0-7271-4cad-89d0-dc628f76c6d3"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
+TestItemRunner = "f8b46487-2199-4994-9208-9a1283c18c0a"

test/runtests.jl

@@ -1,34 +1,20 @@
-using SafeTestsets
 using BPGates
+using TestItemRunner
 
-function doset(descr)
-    if length(ARGS) == 0
-        return true
-    end
-    for a in ARGS
-        if occursin(lowercase(a), lowercase(descr))
-            return true
-        end
-    end
-    return false
-end
+# filter for the test
+testfilter = ti -> begin
+  exclude = Symbol[]
+  if get(ENV,"JET_TEST","")!="true"
+    push!(exclude, :jet)
+  end
+  if !(VERSION >= v"1.10")
+    push!(exclude, :doctests)
+    push!(exclude, :aqua)
+  end
 
-macro doset(descr)
-    quote
-        if doset($descr)
-            @safetestset $descr begin include("test_"*$descr*".jl") end
-        end
-    end
+  return all(!in(exclude), ti.tags)
 end
 
 println("Starting tests with $(Threads.nthreads()) threads out of `Sys.CPU_THREADS = $(Sys.CPU_THREADS)`...")
 
-@doset "quantumclifford"
-@doset "quantikz"
-get(ENV,"JET_TEST","")=="true" && @doset "jet"
-@doset "doctests"
-
-using Aqua
-doset("aqua") && begin
-    Aqua.test_all(BPGates)
-end
+@run_package_tests filter=testfilter

test/test_aqua.jl

@@ -0,0 +1,7 @@
+@testitem "Aqua analysis" begin
+
+using Aqua, BPGates
+
+Aqua.test_all(BPGates)
+
+end

test/test_doctests.jl

@@ -1,5 +1,9 @@
+@testitem "doctests" begin
+
 using Documenter
 using BPGates
 
 DocMeta.setdocmeta!(BPGates, :DocTestSetup, :(using QuantumClifford, BPGates); recursive=true)
 doctest(BPGates)
+
+end

test/test_jet.jl

@@ -1,3 +1,5 @@
+@testitem "JET analysis" tags=[:jet] begin
+
 using JET
 using Test
 using BPGates
@@ -9,3 +11,5 @@ rep = report_package("BPGates";
 )
 @show rep
 @test length(JET.get_reports(rep)) == 0
+
+end

test/test_quantikz.jl

@@ -1,3 +1,5 @@
+@testitem "Quantikz circuit plotting" begin
+
 using BPGates
 using Quantikz
 using Test
@@ -6,3 +8,5 @@ using Test
     quantikz_circ = Quantikz.QuantikzOp.([BellSinglePermutation(2,2), BellDoublePermutation(1,2,3), BellPauliPermutation(4,1), BellMeasure(3,1), BellGate(1,1,1,1,1,1,2), CNOTPerm(1,2,3,4), GoodSingleQubitPerm(1,2), PauliNoiseOp(1,0.0,0.0,0.0), PauliNoiseBellGate(CNOTPerm(1,2,3,4),0.0,0.0,0.0), NoisyBellMeasure(BellMeasure(1,2),0.0), NoisyBellMeasureNoisyReset(BellMeasure(1,2),0.0,0.0,0.0,0.0)])
     img = Quantikz.circuit2image(quantikz_circ)
 end
+
+end