From 38f69722640b80ec99719546b7ce28c4cd7c8dba Mon Sep 17 00:00:00 2001 From: Benzillaist Date: Sat, 6 Jan 2024 02:04:57 -0500 Subject: [PATCH] Added bicycle and unicycle tests, minor consistency changes to CSS struct -Added tests for each of the methods used in creating bicycle or unicycle ECCs. - Fixed the final return step of some CSS methods, specifically the constructor and boolean_tableau - Changed name of `code_m` to `code_s` --- .gitignore | 3 +- src/ecc/codes/css.jl | 16 ++++-- src/ecc/codes/simple_sparse_codes.jl | 83 ++++++++++++++++++++++++---- 3 files changed, 85 insertions(+), 17 deletions(-) diff --git a/.gitignore b/.gitignore index ebcf4e4b5..5d180000f 100644 --- a/.gitignore +++ b/.gitignore @@ -3,4 +3,5 @@ Manifest.toml LocalPreferences.toml */.*swp scratch/ -*.cov \ No newline at end of file +*.cov +.vscode/settings.json \ No newline at end of file diff --git a/src/ecc/codes/css.jl b/src/ecc/codes/css.jl index 4c1a92fa3..bd8ce4cf0 100644 --- a/src/ecc/codes/css.jl +++ b/src/ecc/codes/css.jl @@ -1,13 +1,17 @@ -"""An arbitrary CSS error correcting code defined by its X and Z checks.""" +"""An arbitrary CSS error correcting code defined by its X and Z checks. + +Hx: A boolean matrix describing the X checks +Hz: A boolean matrix describing the Z checks +""" struct CSS <: AbstractECC - Hx - Hz + Hx::Matrix{Bool} + Hz::Matrix{Bool} """Creates a CSS code using the two provided matrices where Hx contains the X checks and Hz contains the Z checks.""" function CSS(Hx, Hz) n = size(Hx, 2) if n != size(Hz, 2) error("When constructing a CSS quantum code, the two classical codes are required to have the same block size") end if size(Hx,1)+size(Hz,1) >= n error("When constructing a CSS quantum code, the total number of checks (rows) in the parity checks of the two classical codes have to be lower than the block size (the number of columns).") end - return new(Hx, Hz) + new(Hx, Hz) end end @@ -17,7 +21,7 @@ function boolean_tableau(c::CSS) checks_matrix = falses(Hx_height + Hz_height, Hx_width + Hz_width) checks_matrix[1:Hx_height, 1:Hx_width] = c.Hx checks_matrix[Hx_height+1:end, Hx_width+1:end] = c.Hz - return CSS(checks_matrix) + return checks_matrix end """Returns the stabilizer making up the parity check tableau.""" @@ -31,7 +35,7 @@ end code_n(c::CSS) = size(c.Hx,2) """Returns the depth of the parity check matrix""" -code_m(c::CSS) = size(c.Hx, 1) + size(c.Hz, 1) +code_s(c::CSS) = size(c.Hx, 1) + size(c.Hz, 1) """Returns the number of encoded qubits""" code_k(c::CSS) = (2 * size(c.Hx,2)) - code_m(c) diff --git a/src/ecc/codes/simple_sparse_codes.jl b/src/ecc/codes/simple_sparse_codes.jl index ff1ecd3ab..47bbad63d 100644 --- a/src/ecc/codes/simple_sparse_codes.jl +++ b/src/ecc/codes/simple_sparse_codes.jl @@ -4,13 +4,35 @@ Parameters: - n: width of array, should be >= 2 -- m: height of array, should be >= 2 and a multiple of 2""" +- m: height of array, should be >= 2 and a multiple of 2 + + +``` jldoctest Bicycle +julia> using QuantumClifford.ECC +julia> parity_checks(Bicycle(6, 4)) ++ XX_X_X ++ X_X_XX ++ ZZ_Z_Z ++ Z_Z_ZZ + +julia> Bicycle(6,4).Hx +2×6 Matrix{Bool}: + 1 1 0 1 0 1 + 1 0 1 0 1 1 + + julia> typeof(Bicycle(6, 4)) + CSS +``` +""" function Bicycle(n::Int, m::Int) + if n%2 == 1 + throw(DomainError(m, "The number of qubits should be a multiple for 2 for bicycle codes.")) + end if m%2 == 1 - throw(DomainError(m, " M should be a multiple for 2 for bicycle codes.")) + throw(DomainError(m, "Code depth should be a multiple for 2 for CSS codes.")) end if m < 2 - throw(DomainError(m, " M is too small, make it greater than 1.")) + throw(DomainError(m, "Code depth is too small, make it greater than 1.")) end bs = bicycle_set_gen(Int(n/2)) bsc = circ_to_bicycle_h0(bs, Int(n/2)) @@ -24,11 +46,21 @@ end Parameters: - n: width of array, should be >= 1 -- set: array of indices that are 'active' checks in the circulant code""" +- set: array of indices that are 'active' checks in the circulant code + +``` jldoctest Unicycle +julia> Unicycle(7, [1, 2, 4]) +4×8 Matrix{Bool}: + 0 1 1 0 1 0 0 1 + 0 0 0 1 1 0 1 1 + 0 1 0 0 0 1 1 1 + 1 0 1 0 0 0 1 1 +``` +""" function Unicycle(n::Int, set) usc = circ_to_unicycle_h0(set, n) rusc = reduce_unicycle(usc) - return CSS(rusc, rusc) + return rusc end """Takes an untrimmed bicycle matrix and removes the row which keeps the spread of the column weights minimal. @@ -36,7 +68,15 @@ end Required before the bicycle code can be used. Typical usage: -ReduceBicycle(Circ2BicycleH0(array_indices, (block length / 2) ) )""" + +``` jldoctest reduce_bicycle +julia> reduce_bicycle(Bool[1 1 0 1 0 1; 0 1 1 1 1 0; 1 0 1 0 1 1]) +Bool[1 1 0 1 0 1; 1 0 1 0 1 1] + +julia> reduce_bicycle(Bool[1 1 0 0 1 0 0 1; 0 1 1 0 1 1 0 0; 0 0 1 1 0 1 1 0; 1 0 0 1 0 0 1 1]) +Bool[0 1 1 0 1 1 0 0; 0 0 1 1 0 1 1 0; 1 0 0 1 0 0 1 1] +``` +""" function reduce_bicycle(H0::Matrix{Bool}) m, n = size(H0) r_i = 0 @@ -55,8 +95,13 @@ end """Takes a list of indices and creates the base of the bicycle matrix. For example: -Circ2BicycleH0([1, 2, 4], 7) +``` jldoctest circ_to_bicycle_h0 +julia> circ_to_bicycle_h0([0, 1], 3) +Bool[1 1 0 1 0 1; 0 1 1 1 1 0; 1 0 1 0 1 1] +julia> circ_to_bicycle_h0([0, 1], 4) +Bool[1 1 0 0 1 0 0 1; 0 1 1 0 1 1 0 0; 0 0 1 1 0 1 1 0; 1 0 0 1 0 0 1 1] +``` See https://arxiv.org/abs/quant-ph/0304161 for more details""" function circ_to_bicycle_h0(circ_indices::Array{Int}, n::Int) circ_arr = Array{Bool}(undef, n) @@ -85,7 +130,12 @@ end Required before the unicycle code can be used. Typical usage: -`ReduceUnicycle(Circ2UnicycleH0(array_indices, block length) )`""" +`reduce_unicycle(circ_to_unicycle_h0(array_indices, block length) )` + +``` jldoctest reduce_unicycle +julia> reduce_unicycle(Bool[1 1 0 1 0 0 0 1; 0 1 1 0 1 0 0 1; 0 0 1 1 0 1 0 1; 0 0 0 1 1 0 1 1; 1 0 0 0 1 1 0 1; 0 1 0 0 0 1 1 1; 1 0 1 0 0 0 1 1]) +Bool[0 1 1 0 1 0 0 1; 0 0 0 1 1 0 1 1; 0 1 0 0 0 1 1 1; 1 0 1 0 0 0 1 1] +```""" function reduce_unicycle(m::Matrix{Bool}) rrzz = residue_ring(ZZ, 2) nm = matrix(rrzz, m) @@ -107,8 +157,12 @@ end """Takes a list of indices and creates the base of the unicycle matrix. For example: -`Circ2UnicycleH0([1, 2, 4], 7)` +`circ_to_unicycle_h0([1, 2, 4], 7)` +``` jldoctest circ_to_unicycle +julia> circ_to_unicycle([1, 2, 4], 7) +Bool[1 1 0 1 0 0 0 1; 0 1 1 0 1 0 0 1; 0 0 1 1 0 1 0 1; 0 0 0 1 1 0 1 1; 1 0 0 0 1 1 0 1; 0 1 0 0 0 1 1 1; 1 0 1 0 0 0 1 1] +``` See https://arxiv.org/abs/quant-ph/0304161 for more details""" function circ_to_unicycle_h0(circ_indices::Array{Int}, n::Int) circ_arr = fill(false, n) @@ -133,7 +187,16 @@ function circ_to_unicycle_h0(circ_indices::Array{Int}, n::Int) return comp_matrix end -"""Attempts to generate a list of indices to be used in a bicycle code using a search method""" +"""Attempts to generate a list of indices to be used in a bicycle code using a search method + +``` jldoctest bicycle_set_gen +julia> bicycle_set_gen(3) +[0, 1] + +julia> bicycle_set_gen(1000) +[0, 1, 3, 7, 12, 20, 30, 44, 65, 80, 96, 122, 147, 181, 203, 289] +``` +""" function bicycle_set_gen(N::Int) circ_arr::Array{Int} = [0] diff_arr::Array{Int} = []