Multivariate Multicycle codes using Koszul Complexes

ext/QuantumCliffordOscarExt/QuantumCliffordOscarExt.jl

@@ -18,11 +18,13 @@ import Oscar: free_group, small_group_identification, describe, order, FPGroupEl
     kronecker_product, FqMatrix, identity_matrix, iszero, FqPolyRingElem, laurent_polynomial_ring,
     hnf_with_transform, ideal, intersect, ==, is_coprime, quo, groebner_basis, length, FqMPolyRingElem,
     first, length, MPolyQuoRingElem, FqMPolyRingElem, modulus, ideal, monomials, terms, coeff, degree, mod
+import Oscar.koszul_complex
 import Oscar.Generic.MatSpaceElem
 import Oscar.Generic.DirectSumModule
 import Oscar.Generic.LaurentMPolyWrap
 import Oscar.Generic.exponent_vectors
 import Oscar.IdealGens
+import Combinatorics: combinations
 
 import QuantumClifford.ECC: two_block_group_algebra_codes, twobga_from_direct_product, twobga_from_fp_group,
     boundary_maps, max_xy_exponents
@@ -41,6 +43,7 @@ include("generalized_toric.jl")
 include("group_presentation.jl")
 include("d_dimensional_codes.jl")
 include("trivariate_tricycle.jl")
+include("multivariate_multicycle.jl")
 include("homological_product_codes.jl")
 
 end # module

ext/QuantumCliffordOscarExt/multivariate_multicycle.jl

@@ -0,0 +1,198 @@
+"""
+We introduce a novel class of quantum CSS codes — *Multivariate Multicycle* codes — constructed
+from multivariate polynomial quotient ring formalism over finite fields over GF(2). Our discovery establishes
+that the boundary maps of these codes are governed by the combinatorial structure of *Koszul* complexes.
+Specifically, for a code defined by *t* polynomial relations, we demonstrate that the *k-th* boundary map
+is obtained by taking the **Koszul matrix** in degree *k* and replacing each variable entry with the
+corresponding circulant matrix derived from the code's defining relations. The **Koszul complex** provides
+the framework for the boundary map construction, ensuring the commutativity properties essential for the code
+construction. This correspondence reveals that multivariate multicycle codes possess the structure of
+**Koszul complexes** with circulant coefficients. This family of codes generalizes the bivariate bicycle,
+trivariate tricycle ([`TrivariateTricycleCode`](@ref)), and tetravariate tetracycle codes and it enables
+full single shot decoding in both X and Z directions.
+
+# Special Cases 
+
+## t = 2: Bivariate bicycle codes
+
+```jldoctest
+julia> using Oscar; using QuantumClifford.ECC;
+
+julia> l=21; m=18;
+
+julia> R, (x, y) = polynomial_ring(GF(2), [:x, :y]);
+
+julia> I = ideal(R, [x^l-1, y^m-1]);
+
+julia> S, _ = quo(R, I);
+
+julia> A = S(x^3 + y^10 + y^17);
+
+julia> B = S(y^5 + x^3  + x^19);
+
+julia> c = MultivariateMulticycleCode([l,m], [A,B]);
+
+julia> code_n(c), code_k(c)
+(756, 16)
+```
+
+## t = 3: Trivariate tricycle codes
+
+```jldoctest
+julia> using Oscar; using QuantumClifford.ECC;
+
+julia> l, m, p = 6, 6, 4;
+
+julia> R, (x, y, z) = polynomial_ring(GF(2), [:x, :y, :z]);
+
+julia> I = ideal(R, [x^l - 1, y^m - 1, z^p - 1]);
+
+julia> S, _ = quo(R, I);
+
+julia> A = S(1 + x*y*z^3 + x^3*y^4*z^2);
+
+julia> B = S(1 + x^3*y*z^2 + x^3*y^2*z^3);
+
+julia> C = S(1 + x^4*y^3*z^3 + x^5*z^2);
+
+julia> c = MultivariateMulticycleCode([l,m, p], [A, B, C]);
+
+julia> code_n(c), code_k(c)
+(432, 12)
+```
+
+## t = 4: Tetravariate Tetracycle Codes
+
+```jldoctest
+julia> using Oscar; using QuantumClifford.ECC;
+
+julia> l, m, p, r = 4, 3, 3, 3;
+
+julia> R, (w, x, y, z) = polynomial_ring(GF(2), [:w, :x, :y, :z]);
+
+julia> I = ideal(R, [w^l - 1, x^m - 1, y^p - 1, z^r - 1]);
+
+julia> S, _ = quo(R, I);
+
+julia> A = S((1 + x^2 )*(1 + w*x*y*z^2));
+
+julia> B = S((1 + x^2)*(1 + w*x^3*y^2*z));
+
+julia> C = S(1 + w^2*x^2*y^2*z^2);
+
+julia> D = S(1 + w^3*x^3*y^3*z^3);
+
+julia> c = MultivariateMulticycleCode([l, m, p, r], [A, B, C, D]);
+
+julia> code_n(c), code_k(c)
+(648, 18)
+```
+"""
+struct MultivariateMulticycleCode <: AbstractCSSCode
+    orders::Vector{Int}
+    polynomials::Vector{MPolyQuoRingElem{FqMPolyRingElem}}
+    function MultivariateMulticycleCode(orders::Vector{Int}, polys::Vector{<:MPolyQuoRingElem})
+        length(orders) == length(polys) || throw(ArgumentError("Mismatch orders/polys"))
+        all(x->x>0, orders) || throw(ArgumentError("All orders must be positive"))
+        length(orders) ≥ 2 || throw(ArgumentError("Need at least 2 variables to define a CSS code"))
+        new(orders, collect(polys))
+    end
+end
+
+_gf2_to_int(M) = [iszero(M[i,j]) ? 0 : 1 for i in 1:size(M,1), j in 1:size(M,2)]
+
+function _polynomial_to_circulant_matrix(f::MPolyQuoRingElem, orders::Vector{Int})
+    t = length(orders)
+    N = prod(orders)
+    M = zero_matrix(GF(2), N, N)
+    f_lift = lift(f)
+    for col_idx in 0:(N-1)
+        tmp = col_idx
+        idxs = Vector{Int}(undef, t)
+        for k in t:-1:1
+            idxs[k] = tmp % orders[k]
+            tmp ÷= orders[k]
+        end
+        for term in terms(f_lift)
+            c = coeff(term, 1)
+            iszero(c) && continue
+            exps = [degree(term, k) for k in 1:t]
+            row_idx = 0
+            for k in 1:t
+                row_idx = row_idx*orders[k]+mod(idxs[k]+exps[k], orders[k])
+            end
+            M[row_idx+1, col_idx+1] += c
+        end
+    end
+    return M
+end
+
+function boundary_maps(code::MultivariateMulticycleCode)
+    t = length(code.orders)
+    N = prod(code.orders)
+    circs = [_gf2_to_int(_polynomial_to_circulant_matrix(p, code.orders)) for p in code.polynomials]
+    maps = Vector{Matrix{Int}}(undef, t+1)
+    for k in 0:t
+        R, x = polynomial_ring(GF(2), ["x$i" for i in 1:t])
+        K = koszul_complex(x)
+        boundary_map = map(K, k)
+        KoszulMatrix = matrix(boundary_map)
+        nr, nc = size(KoszulMatrix)
+        M = zeros(Int, nr*N, nc*N)
+        for i in 1:nr, j in 1:nc
+            e = KoszulMatrix[i, j]
+            !iszero(e) || continue
+            for idx in 1:t
+                if e == x[idx]
+                    row_range = ((i-1)*N+1):(i*N)
+                    col_range = ((j-1)*N+1):(j*N)
+                    M[row_range, col_range] .= circs[idx]
+                    break
+                end
+            end
+        end
+        maps[k+1] = M
+    end
+    for k in 1:t
+        if !isempty(maps[k]) && !isempty(maps[k+1])
+            @assert iszero(mod.(maps[k]*maps[k-1],2)) "Exactness check failed: ∂$(k-1)∘∂$k ≠ 0 mod 2"
+        end
+    end
+    return maps
+end
+
+function parity_matrix_xz(code::MultivariateMulticycleCode)
+    maps = boundary_maps(code)
+    t = length(code.orders)
+    if t == 2
+        hx = maps[3]
+        hz = transpose(maps[2])
+    else
+        qd = t ÷ 2
+        hx = maps[qd+2]
+        hz = transpose(maps[qd+1])
+    end
+    return hx, hz
+end
+
+parity_matrix_x(c::MultivariateMulticycleCode) = parity_matrix_xz(c)[1]
+
+parity_matrix_z(c::MultivariateMulticycleCode) = parity_matrix_xz(c)[2]
+
+"""For t ≥ 4, provides metachecks for X-stabilizers enabling single-shot decoding."""
+function metacheck_matrix_x(code::MultivariateMulticycleCode)
+    maps = boundary_maps(code)
+    t = length(code.orders)
+    t ≥ 4 || throw(ArgumentError("X-metachecks require t ≥ 4 variables"))
+    qd = t÷2
+    return transpose(maps[qd+1])
+end
+
+"""For t ≥ 3, provides metachecks for Z-stabilizers enabling single-shot decoding."""
+function metacheck_matrix_z(code::MultivariateMulticycleCode)
+    maps = boundary_maps(code)
+    t = length(code.orders)
+    t ≥ 3 || throw(ArgumentError("Z-metachecks require t ≥ 3 variables"))
+    qd = t÷2
+    return maps[qd+2]
+end

src/ecc/ECC.jl

@@ -41,7 +41,7 @@ export parity_checks, parity_matrix_x, parity_matrix_z, iscss,
     GeneralizedCirculantBivariateBicycle, GeneralizedHyperGraphProductCode,
     GeneralizedBicycleCode, ExtendedGeneralizedBicycleCode,
     HomologicalProductCode, DoubleHomologicalProductCode,
-    GeneralizedToricCode, TrivariateTricycleCode,
+    GeneralizedToricCode, TrivariateTricycleCode, MultivariateMulticycleCode,
     evaluate_decoder,
     CommutationCheckECCSetup, NaiveSyndromeECCSetup, ShorSyndromeECCSetup,
     TableDecoder,

src/ecc/codes/qeccs_using_oscar.jl

@@ -74,3 +74,13 @@ function TrivariateTricycleCode(args...; kwargs...)
     end
     return ext.TrivariateTricycleCode(args...; kwargs...)
 end
+
+"""Multivariate Multicycle codes.
+Implemented as a package extension with `Oscar`. Check the [QuantumClifford documentation](http://qc.quantumsavory.org/stable/ECC_API/) for more details on that extension."""
+function MultivariateMulticycleCode(args...; kwargs...)
+    ext = Base.get_extension(QuantumClifford, :QuantumCliffordOscarExt)
+    if isnothing(ext)
+        throw("The `MultivariateMulticycleCode` depends on the package `Oscar` but you have not installed or imported it yet. Immediately after you import `Oscar`, the `MultivariateMulticycleCode` will be available.")
+    end
+    return ext.MultivariateMulticycleCode(args...; kwargs...)
+end

test/test_ecc_base.jl

@@ -523,13 +523,53 @@ const code_instance_args = Dict(
     α1₆ = (0,  3)
     α2₆ = (19, 1)
 
+    # Multivariate Multicycle Codes
+    # t = 2; Bivariate Bicycle codes
+    # [[72, 12, 6]]
+    l_mm₁ =6; m_mm₁ = 6
+    R_mm₁, (x, y) = polynomial_ring(GF(2), [:x, :y])
+    I_mm₁ = ideal(R_mm₁, [x^l_mm₁-1, y^m_mm₁-1])
+    S_mm₁, _ = quo(R_mm₁, I_mm₁)
+    A_mm₁ = S_mm₁(x^3 + y + y^2)
+    B_mm₁ = S_mm₁(y^3 + x + x^2)
+
+    # [[90, 8, 10]]
+    l_mm₂ =15; m_mm₂ = 3
+    R_mm₂, (x, y) = polynomial_ring(GF(2), [:x, :y])
+    I_mm₂ = ideal(R_mm₂, [x^l_mm₂-1, y^m_mm₂-1])
+    S_mm₂, _ = quo(R_mm₂, I_mm₂)
+    A_mm₂ = S_mm₂(x^9 + y   + y^2)
+    B_mm₂ = S_mm₂(1   + x^2 + x^7)
+
+    # t = 3; Trivariate Tricycle codes
+    # [[60, 3, 4]]
+    ℓ_mm₃, m_mm₃, p_mm₃ = 5, 2, 2
+    F2_mm₃ = GF(2)
+    R_mm₃, (x, y, z) = polynomial_ring(F2_mm₃, [:x, :y, :z])
+    I_mm₃ = ideal(R_mm₃, [x^ℓ_mm₃ - 1, y^m_mm₃ - 1, z^p_mm₃ - 1])
+    S_mm₃, _ = quo(R_mm₃, I_mm₃)
+    A_mm₃ = S_mm₃(1 + x*z)
+    B_mm₃ = S_mm₃(1 + x*y)
+    C_mm₃ = S_mm₃(1 + x*y*z)   
+
+    # [[90, 3, 5]]
+    ℓ_mm₄, m_mm₄, p_mm₄ = 5, 3, 2
+    F2_mm₄ = GF(2)
+    R_mm₄, (x, y, z) = polynomial_ring(F2_mm₄, [:x, :y, :z])
+    I_mm₄ = ideal(R_mm₄, [x^ℓ_mm₄ - 1, y^m_mm₄ - 1, z^p_mm₄ - 1])
+    S_mm₄, _ = quo(R_mm₄, I_mm₄)
+    A_mm₄ = S_mm₄(1 + x)
+    B_mm₄ = S_mm₄(1 + x*y)
+    C_mm₄ = S_mm₄(1 + x^2*y^2*z)
+
     oscar_code_instance_args = Dict(
         :DDimensionalSurfaceCode => [(2, 3), (3, 2), (3, 3), (4, 2)],
         :DDimensionalToricCode => [(2, 3), (3, 2), (3, 3), (4, 2)],
         :GeneralizedToricCode => [(f₁, g₁, α1₁, α2₁), (f₂, g₂, α1₂, α2₂), (f₃, g₃, α1₃, α2₃), (f₄, g₄, α1₄, α2₄), (f₅, g₅, α1₅, α2₅), (f₆, g₆, α1₆, α2₆)],
         :HomologicalProductCode => [([H₁, transpose(H₁)], l₁), ([H₂, transpose(H₂)], l₂), ([H₃, transpose(H₃)],), ([δ₄, δ₄, δ₄],)],
         :DoubleHomologicalProductCode => [(δ₁,), (δ₂,)],
-        :TrivariateTricycleCode => [(ℓ₁, m₁, p₁, A₁, B₁, C₁), (ℓ₂, m₂, p₂, A₂, B₂, C₂), (ℓ₃, m₃, p₃, A₃, B₃, C₃), (ℓ₄, m₄, p₄, A₄, B₄, C₄)]
+        :TrivariateTricycleCode => [(ℓ₁, m₁, p₁, A₁, B₁, C₁), (ℓ₂, m₂, p₂, A₂, B₂, C₂), (ℓ₃, m₃, p₃, A₃, B₃, C₃), (ℓ₄, m₄, p₄, A₄, B₄, C₄)],
+        :MultivariateMulticycleCode =>[([l_mm₁,m_mm₁], [A_mm₁, B_mm₁]), ([l_mm₂,m_mm₂], [A_mm₂, B_mm₂]), ([ℓ_mm₃, m_mm₃, p_mm₃], [A_mm₃, B_mm₃, C_mm₃]), ([ℓ_mm₄, m_mm₄, p_mm₄], [A_mm₄, B_mm₄, C_mm₄])]
     )
     merge!(code_instance_args, oscar_code_instance_args)
   end

test/test_ecc_trivariate_tricycle.jl

@@ -78,6 +78,16 @@
                 @test code_k(c) == k == code_k(stab)
                 @test stab_looks_good(stab, remove_redundant_rows=true) == true
                 @test iszero(mod.(metacheck_matrix_z(c)*parity_matrix_z(c), 2))
+                c = MultivariateMulticycleCode([ℓ, m, p], [A, B, C])
+                stab = parity_checks(c)
+                mat = matrix(GF(2), stab_to_gf2(stab))
+                computed_rank = rank(mat)
+                @test computed_rank == code_n(c) - code_k(c)
+                # A TT code is defined on n = 3*ℓ*m*p data qubits.
+                @test code_n(c) == n == code_n(stab) == 3*ℓ*m*p
+                @test code_k(c) == k == code_k(stab)
+                @test stab_looks_good(stab, remove_redundant_rows=true) == true
+                @test iszero(mod.(metacheck_matrix_z(c)*parity_matrix_z(c)', 2))
             end
         end
     end