Error correction methods for creating Bicycle and Unicycle codes + belief decoder code evaluation

Project.toml

@@ -11,13 +11,16 @@ Graphs = "86223c79-3864-5bf0-83f7-82e725a168b6"
 HostCPUFeatures = "3e5b6fbb-0976-4d2c-9146-d79de83f2fb0"
 ILog2 = "2cd5bd5f-40a1-5050-9e10-fc8cdb6109f5"
 InteractiveUtils = "b77e0a4c-d291-57a0-90e8-8db25a27a240"
+LDPCDecoders = "3c486d74-64b9-4c60-8b1a-13a564e77efb"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 MacroTools = "1914dd2f-81c6-5fcd-8719-6d5c9610ff09"
 Nemo = "2edaba10-b0f1-5616-af89-8c11ac63239a"
 PrecompileTools = "aea7be01-6a6a-4083-8856-8a6e6704d82a"
 QuantumInterface = "5717a53b-5d69-4fa3-b976-0bf2f97ca1e5"
 Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
 SIMD = "fdea26ae-647d-5447-a871-4b548cad5224"
+SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
+Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
 SumTypes = "8e1ec7a9-0e02-4297-b0fe-6433085c89f2"
 
 [weakdeps]

src/ecc/ECC.jl

@@ -1,18 +1,23 @@
 module ECC
 
 using LinearAlgebra
-using QuantumClifford
-using QuantumClifford: AbstractOperation, AbstractStabilizer
+using QuantumClifford, SparseArrays, LDPCDecoders
+using QuantumClifford: AbstractOperation, AbstractStabilizer, Stabilizer
 import QuantumClifford: Stabilizer, MixedDestabilizer
 using DocStringExtensions
 using Combinatorics: combinations
+using Statistics: std
+using LinearAlgebra: rank
+using Nemo: ZZ, residue_ring, matrix
 
 abstract type AbstractECC end
 
 export Shor9, Steane7, Cleve8, Perfect5, Bitflip3,
     parity_checks, naive_syndrome_circuit, shor_syndrome_circuit, naive_encoding_circuit,
     code_n, code_s, code_k, rate, distance,
-    isdegenerate, faults_matrix
+    isdegenerate, faults_matrix,
+    Unicycle, Bicycle,
+    CSS
 
 """Parity check tableau of a code."""
 function parity_checks end
@@ -289,6 +294,8 @@ function isdegenerate(H::Stabilizer, d::Int=1)
 end
 
 include("circuits.jl")
+include("css.jl")
+include("simple_sparse_codes.jl")
 
 include("codes/bitflipcode.jl")
 include("codes/fivequbit.jl")

src/ecc/css.jl

@@ -0,0 +1,22 @@
+"""An arbitrary CSS error correcting code defined by its X and Z checks."""
+struct CSS <: AbstractECC
+    tab
+end
+
+"""Creates a CSS code using the two provided matrices where H contains the X checks and G contains the Z checks."""
+function CSS(H, G)
+    Hy, Hx = size(H)
+    Gy, Gx = size(G)
+    comp_matrix = falses(Hy + Gy, Hx + Gx)
+    comp_matrix[1:Hy, 1:Hx] = H
+    comp_matrix[Hy+1:end, Hx+1:end] = G
+    return CSS(comp_matrix)
+end
+
+"""Returns the stabilizer making up the parity check tableau."""
+function parity_checks(c::CSS)
+    Stabilizer(fill(0x0, size(c.tab, 2)), c.tab[:,1:end÷2], c.tab[:,end÷2+1:end])
+end
+
+"""Returns the block length of the code."""
+code_n(c::CSS) = size(c.tab,2)÷2

src/ecc/simple_sparse_codes.jl

@@ -0,0 +1,231 @@
+# Currently just has Bicycle and Unicycle codes, but open to all types of rudimentary sparse codes
+
+"""Takes a height and width of matrix and generates a bicycle code to the specified height and width.
+
+Parameters:
+- n: width of array, should be >= 2
+- m: height of array, should be >= 2 and a multiple of 2"""
+function Bicycle(n::Int, m::Int)
+    if m%2 == 1
+        throw(DomainError(m, " M should be a multiple for 2 for bicycle codes."))
+    end
+    if m < 2
+        throw(DomainError(m, " M 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))
+    while size(bsc)[1] > m/2
+        bsc = reduce_bicycle(bsc)
+    end
+    return CSS(bsc, bsc)
+end
+
+"""Takes a height and width of matrix and generates a bicycle code to the specified height and width.
+
+Parameters:
+- n: width of array, should be >= 1
+- set: array of indices that are 'active' checks in the circulant code"""
+function Unicycle(n::Int, set)
+    usc = circ_to_unicycle_h0(set, n)
+    rusc = reduce_unicycle(usc)
+    return CSS(rusc, rusc)
+end
+
+"""Takes an untrimmed bicycle matrix and removes the row which keeps the spread of the column weights minimal.
+
+Required before the bicycle code can be used.
+
+Typical usage:
+ReduceBicycle(Circ2BicycleH0(array_indices, (block length / 2) ) )"""
+function reduce_bicycle(H0::Matrix{Bool})
+    m, n = size(H0)
+    r_i = 0
+    std_min = Inf
+    for i in 1:m
+        t_H0 = vcat(H0[1:i-1, :], H0[i+1:end, :])
+        std_temp = std(convert(Array, sum(t_H0, dims = 1)))
+        if std_temp < std_min
+            std_min = std_temp
+            r_i = i
+        end
+    end
+    return vcat(H0[1:r_i-1, :], H0[r_i+1:end, :])
+end
+
+"""Takes a list of indices and creates the base of the bicycle matrix.
+
+For example:
+Circ2BicycleH0([1, 2, 4], 7)
+
+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)
+    circ_matrix = Matrix{Bool}(undef, n, n)
+    comp_matrix = Matrix{Bool}(undef, n, 2*n)
+    for i = 1:n
+        if Int(i-1) in circ_indices
+            circ_arr[i] = true
+        else
+            circ_arr[i] = false
+        end
+    end
+    for i = 1:n
+        circ_matrix[i,1:n] = circ_arr
+        li = circ_arr[end]
+        circ_arr[2:end] = circ_arr[1:end-1]
+        circ_arr[1] = li
+    end
+    comp_matrix[1:n,1:n] = circ_matrix
+    comp_matrix[1:n,n+1:2*n] = transpose(circ_matrix)
+    return comp_matrix
+end
+
+"""Takes an untrimmed unicycle matrix and removes linearly dependent rows.
+
+Required before the unicycle code can be used.
+
+Typical usage:
+ReduceUnicycle(Circ2UnicycleH0(array_indices, block length) )"""
+function reduce_unicycle(m::Matrix{Bool})
+    rrzz = residue_ring(ZZ, 2)
+    nm = matrix(rrzz, m)
+    r = LinearAlgebra.rank(nm)
+    for i in 1:size(m)[1]
+        tm = vcat(m[1:i-1,:], m[i+1:end,:])
+        tr = LinearAlgebra.rank(matrix(rrzz, tm))
+        if(tr == r)
+            m = tm
+            i -= 1
+            if(size(m)[1] == r)
+                break
+            end
+        end
+    end
+    return m
+end
+
+"""Takes a list of indices and creates the base of the unicycle matrix.
+
+For example:
+Circ2UnicycleH0([1, 2, 4], 7)
+
+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)
+    one_col = transpose(fill(true, n))
+    circ_matrix = Matrix{Bool}(undef, n, n)
+    comp_matrix = Matrix{Bool}(undef, n, n+1)
+    for i = 1:n
+        if i in circ_indices
+            circ_arr[i] = true
+        else
+            circ_arr[i] = false
+        end
+    end
+    for i = 1:n
+        circ_matrix[i,1:n] = circ_arr
+        li = circ_arr[end]
+        circ_arr[2:end] = circ_arr[1:end-1]
+        circ_arr[1] = li
+    end
+    comp_matrix[1:n,1:n] = circ_matrix
+    comp_matrix[1:n,n+1] = one_col
+    return comp_matrix
+end
+
+"""Attempts to generate a list of indices to be used in a bicycle code using a search method"""
+function bicycle_set_gen(N::Int)
+    circ_arr::Array{Int} = [0]
+    diff_arr::Array{Int} = []
+    circ_arr[1] = 0
+    # test new elements
+    for add_i = (circ_arr[end] + 1):N - 1
+        valid = true
+        temp_circ_arr = copy(circ_arr)
+        temp_diff_arr::Array{Int} = []
+        push!(temp_circ_arr, add_i)
+        for j = 1:size(temp_circ_arr)[1]
+            temp_arr = copy(temp_circ_arr)
+            # add lesser elements + N to temp_arr
+            for k = 1:size(temp_circ_arr)[1]
+                if k < j
+                    push!(temp_arr, temp_circ_arr[k] + N)
+                else
+                    break
+                end
+            end
+            # test if new index is valid
+            for k = 1:(size(temp_circ_arr)[1] - 2)
+                t_diff = (temp_arr[j + k] - temp_arr[j]) % N
+                if ((t_diff) in temp_diff_arr)
+                    valid = false
+                    break
+                else
+                    push!(temp_diff_arr, t_diff)
+                end
+            end
+            if !valid
+                break
+            end
+        end
+        if valid
+            circ_arr = copy(temp_circ_arr)
+            diff_arr = copy(temp_diff_arr)
+        end
+    end
+    return circ_arr
+end
+
+"""Attempts to generate a list of indices to be used in a bicycle code using a randomized check method
+
+Note: This is very slow for large N"""
+function bicycle_set_gen_rand(N::Int, d::Int)
+    circ_arr::Array{Int} = [0]
+    diff_arr::Array{Int} = []
+    atmp_add::Array{Int} = [0]
+    circ_arr[1] = 0
+    # test new elements
+    for i = (circ_arr[end] + 1):(N^2)
+        valid = true
+        temp_circ_arr = copy(circ_arr)
+        temp_diff_arr::Array{Int} = []
+        add_i = rand(1: N-1)
+        atmp_add = push!(atmp_add, add_i)
+        if add_i in circ_arr
+            continue
+        end
+        push!(temp_circ_arr, add_i)
+        for j = 1:size(temp_circ_arr)[1]
+            temp_arr = copy(temp_circ_arr)
+            # add lesser elements + N to temp_arr
+            for k = 1:size(temp_circ_arr)[1]
+                if k < j
+                    push!(temp_arr, temp_circ_arr[k] + N)
+                else
+                    break
+                end
+            end
+            # test if new index is valid
+            for k = 1:(size(temp_circ_arr)[1] - 2)
+                t_diff = (temp_arr[j + k] - temp_arr[j]) % N
+                if ((t_diff) in temp_diff_arr)
+                    valid = false
+                    break
+                else
+                    push!(temp_diff_arr, t_diff)
+                end
+            end
+            if !valid
+                break
+            end
+        end
+        if valid
+            circ_arr = copy(temp_circ_arr)
+            diff_arr = copy(temp_diff_arr)
+            if (size(atmp_add)[1] == N) || (size(circ_arr)[1] == d)
+                break
+            end
+        end
+    end
+    return circ_arr
+end