Full ECC API (autogenerated)
QuantumClifford.ECC.CircuitCode
— TypeCircuitCode
is defined by a given encoding circuit circ
.
n
: qubit numbercirc
: the encoding circuitencode_qubits
: the qubits to be encoded
See also: random_all_to_all_circuit_code
, random_brickwork_circuit_code
QuantumClifford.ECC.CommutationCheckECCSetup
— TypeConfiguration for ECC evaluator that does not simulate any ECC circuits, rather it simply checks the commutation of the parity check and the Pauli error.
This is much faster than any other simulation method, but it is incapable of noisy-circuit simulations and thus useless for fault-tolerance studies.
See also: NaiveSyndromeECCSetup
, ShorSyndromeECCSetup
QuantumClifford.ECC.Concat
— TypeConcat(c₁, c₂)
is a code concatenation of two quantum codes (Knill and Laflamme, 1996).
The inner code c₁ and the outer code c₂. The construction is the following: replace each qubit in code c₂ with logical qubits encoded by code c₁. The resulting code will have n = n₁ × n₂
qubits and k = k₁ × k₂
logical qubits.
QuantumClifford.ECC.DistanceMIPAlgorithm
— Typestruct DistanceMIPAlgorithm <: QECCore.AbstractDistanceAlg
A Mixed Integer Programming (MIP) method for computing the code distance of CSS stabilizer codes by finding the minimum-weight non-trivial logical PauliOperator
(either X
-type or Z
-type). Used with distance
to select MIP as the method of finding the distance of a code.
- Requires a
JuMP
-compatible MIP solver (e.g.,HiGHS
,SCIP
). X
-type andZ
-type logical operators yield identical code distance results.- For stabilizer codes, the
X
-distance andZ
-distance are equal.
upper_bound
: iftrue
(default=false
), uses the provided value as an upper bound for the code distancelogical_qubit
: index of the logical qubit to compute code distance for (nothing means compute for all logical qubits)logical_operator_type
: type of logical operator to consider (:X or :Z, defaults to :X) - both types yield identical distance results for CSS stabilizer codes.solver
:JuMP
-compatible MIP solver (e.g.,HiGHS
,SCIP
)opt_summary
: whentrue
(default=false
), prints the MIP solver's solution summarytime_limit
: time limit (in seconds) for the MIP solver's execution (default=60.0)
QuantumClifford.ECC.Gottesman
— TypeThe family of [[2ʲ, 2ʲ - j - 2, 3]]
Gottesman codes, also known as quantum Hamming codes, as described in Gottesman's 1997 PhD thesis and in (Gottesman, 1996).
You might be interested in consulting (Yu et al., 2013) and (Chao and Reichardt, 2017) as well.
The ECC Zoo has an entry for this family
QuantumClifford.ECC.NaiveSyndromeECCSetup
— TypeConfiguration for ECC evaluator that runs the simplest syndrome measurement circuit.
The circuit is being simulated (as opposed to doing only a quick commutation check). This circuit would give poor performance if there is non-zero gate noise.
See also: CommutationCheckECCSetup
, ShorSyndromeECCSetup
QuantumClifford.ECC.QuantumReedMuller
— TypeThe family of [[2ᵐ - 1, 1, 3]]
CSS Quantum-Reed-Muller codes, as discovered by Steane in his 1999 paper (Steane, 1999).
Quantum codes are constructed from shortened Reed-Muller codes RM(1, m)
, by removing the first row and column of the generator matrix Gₘ
. Similarly, we can define truncated dual codes RM(m - 2, m)
using the generator matrix Hₘ
(Anderson et al., 2014). The quantum Reed-Muller codes QRM(m)
derived from RM(1, m)
are CSS codes.
Given that the stabilizers of the quantum code are defined through the generator matrix of the classical code, the minimum distance of the quantum code corresponds to the minimum distance of the dual classical code, which is d = 3
, thus it can correct any single qubit error. Since one stabilizer from the original and one from the dual code are removed in the truncation process, the code parameters are [[2ᵐ - 1, 1, 3]]
.
You might be interested in consulting (Anderson et al., 2014) and (Campbell et al., 2012) as well.
The ECC Zoo has an entry for this family.
QuantumClifford.ECC.ShorSyndromeECCSetup
— TypeConfiguration for ECC evaluators that simulate the Shor-style syndrome measurement (without a flag qubit).
The simulated circuit includes:
- perfect noiseless encoding (encoding and its fault tolerance are not being studied here)
- one round of "memory noise" after the encoding but before the syndrome measurement
- perfect preparation of entangled ancillary qubits
- noisy Shor-style syndrome measurement (only two-qubit gate noise)
- noiseless "logical state measurement" (providing the comparison data when evaluating the decoder)
See also: CommutationCheckECCSetup
, NaiveSyndromeECCSetup
QuantumClifford.ECC.TableDecoder
— TypeA simple look-up table decoder for error correcting codes.
The lookup table contains only weight=1 errors, thus it is small, but at best it provides only for distance=3 decoding.
The size of the lookup table would grow exponentially quickly for higher distances.
QuantumClifford.ECC.Triangular488
— TypeTriangular code following the 4.8.8 tiling. Constructor take a distance d as input.
Example
Here is [[17,1, 5]]
color code following the 4.8.8
tiling:
julia> import HiGHS; import JuMP; # hide
julia> using QuantumClifford.ECC: Triangular666, DistanceMIPAlgorithm; # hide
julia> c = Triangular488(5);
julia> code = Stabilizer(c)
+ XXXX_____________
+ X_X_XX___________
+ __XX_XX__XX__XX__
+ ____XX__XX_______
+ ______XX__XX_____
+ _______X___X___XX
+ ________XX__XX___
+ __________XX__XX_
+ ZZZZ_____________
+ Z_Z_ZZ___________
+ __ZZ_ZZ__ZZ__ZZ__
+ ____ZZ__ZZ_______
+ ______ZZ__ZZ_____
+ _______Z___Z___ZZ
+ ________ZZ__ZZ___
+ __________ZZ__ZZ_
julia> distance(c, DistanceMIPAlgorithm(solver=HiGHS))
5
More information can be seen in (Landahl et al., 2011)
QuantumClifford.ECC.Triangular666
— TypeTriangular code following the 6.6.6 tiling. Constructor take a distance d as input.
Example
Here is [[19,1, 5]]
color code following the 6.6.6
tiling:
julia> import HiGHS; import JuMP; # hide
julia> using QuantumClifford.ECC: Triangular666, DistanceMIPAlgorithm; # hide
julia> c = Triangular666(5);
julia> code = Stabilizer(c)
+ XXXX_______________
+ _X_X_XX____________
+ __XXXX_XX__________
+ ____X__X__XX_______
+ _________X___X___XX
+ _____XX_XX__XX_____
+ _______XX__XX__XX__
+ __________XX__XX___
+ ____________XX__XX_
+ ZZZZ_______________
+ _Z_Z_ZZ____________
+ __ZZZZ_ZZ__________
+ ____Z__Z__ZZ_______
+ _________Z___Z___ZZ
+ _____ZZ_ZZ__ZZ_____
+ _______ZZ__ZZ__ZZ__
+ __________ZZ__ZZ___
+ ____________ZZ__ZZ_
julia> distance(c, DistanceMIPAlgorithm(solver=HiGHS))
5
More information can be seen in (Landahl et al., 2011)
QECCore.code_k
— MethodThe number of logical qubits in a code.
Note that when redundant rows exist in the parity check matrix, the number of logical qubits code_k(c)
will be greater than code_n(c) - code_s(c)
, where the difference equals the redundancy.
QECCore.parity_matrix_x
— MethodParity check boolean matrix of a code (only the X entries in the tableau, i.e. the checks for Z errors).
Only CSS codes have this method.
See also: parity_checks
QECCore.parity_matrix_z
— MethodParity check boolean matrix of a code (only the Z entries in the tableau, i.e. the checks for X errors).
Only CSS codes have this method.
See also: parity_checks
QuantumClifford.ECC.BeliefPropDecoder
— MethodA simple Belief Propagation decoder built around tools from LDPCDecoders.jl
.
QuantumClifford.ECC.BitFlipDecoder
— MethodAn Iterative Bitflip decoder built around tools from LDPCDecoders.jl
.
QuantumClifford.ECC.LPCode
— MethodLifted product codes ((Panteleev and Kalachev, 2021), (Panteleev and Kalachev, Jun 2022))
Implemented as a package extension with Hecke. Check the QuantumClifford documentation for more details on that extension.
QuantumClifford.ECC.LiftedCode
— MethodClassical codes lifted over a group algebra, used for lifted product code construction ((Panteleev and Kalachev, 2021), (Panteleev and Kalachev, Jun 2022))
Implemented as a package extension with Hecke. Check the QuantumClifford documentation for more details on that extension.
QuantumClifford.ECC.PyBeliefPropDecoder
— MethodA Belief Propagation decoder built around tools from the python package ldpc
available from the julia package PyQDecoders.jl
.
QuantumClifford.ECC.PyBeliefPropOSDecoder
— MethodA Belief Propagation decoder with ordered statistics decoding, built around tools from the python package ldpc
available from the julia package PyQDecoders.jl
.
QuantumClifford.ECC.PyMatchingDecoder
— MethodA perfect matching decoder built around tools from the python package pymatching
available from the julia package PyQDecoders.jl
.
QuantumClifford.ECC.bicycle_codes
— FunctionImplemented in a package extension with Hecke. Check the docs for the Hecke extension
QuantumClifford.ECC.evaluate_decoder
— MethodEvaluate the performance of a given decoder (e.g. TableDecoder
) and a given style of running an ECC code (e.g. ShorSyndromeECCSetup
)
QuantumClifford.ECC.evaluate_decoder
— MethodEvaluate the performance of an error-correcting circuit.
This method requires you give the circuit that performs both syndrome measurements and (probably noiseless) logical state measurements. The faults matrix that translates an error vector into corresponding logical errors is necessary as well.
This is a relatively barebones method that assumes the user prepares necessary circuits, etc. It is a method that is used internally by more user-frienly methods providing automatic conversion of codes and noise models to the necessary noisy circuits.
QuantumClifford.ECC.faults_matrix
— MethodError-to-logical-observable map (a.k.a. fault matrix) of a code.
For a code with n physical qubits and k logical qubits this function returns a 2k × 2n binary matrix O such that O[i,j]
is true if the logical observable of index i
is flipped by the single physical qubit error of index j
. Indexing is such that:
O[1:k,:]
is the error-to-logical-X-observable map (logical X observable, i.e. triggered by logical Z errors)O[k+1:2k,:]
is the error-to-logical-Z-observable mapO[:,1:n]
is the X-physical-error-to-logical-observable mapO[n+1:2n,:]
is the Z-physical-error-to-logical-observable map
E.g. for k=1
, n=10
, then if O[2,5]
is true, then the logical Z observable is flipped by a X₅ error; and if O[1,12]
is true, then the logical X observable is flipped by a Z₂ error.
Of note is that there is a lot of freedom in choosing the logical operations! A logical operator multiplied by a stabilizer operator is still a logical operator. Similarly there is a different fault matrix for each choice of logical operators. But once the logical operators are picked, the fault matrix is fixed.
Below we show an example that uses the Shor code. While it is not the smallest code, it is a convenient choice to showcase the importance of the fault matrix when dealing with degenerate codes where a correction operation and an error do not need to be the same.
First, consider a single-qubit error, potential correction operations, and their effect on the Shor code:
julia> using QuantumClifford.ECC: faults_matrix, Shor9
julia> state = MixedDestabilizer(Shor9())
𝒟ℯ𝓈𝓉𝒶𝒷━━━━━
+ Z________
+ ___Z_____
+ _X_______
+ __X______
+ ____X____
+ _____X___
+ ______X__
+ _______X_
𝒳ₗ━━━━━━━━━
+ ______XXX
𝒮𝓉𝒶𝒷━━━━━━━
+ XXX___XXX
+ ___XXXXXX
+ ZZ_______
+ Z_Z______
+ ___ZZ____
+ ___Z_Z___
+ ______Z_Z
+ _______ZZ
𝒵ₗ━━━━━━━━━
+ Z__Z____Z
julia> err_Z₁ = single_z(9,1) # the error we will simulate
+ Z________
julia> cor_Z₂ = single_z(9,2) # the correction operation we will perform
+ _Z_______
julia> err_Z₁ * state # observe that one of the syndrome bits is now flipped
𝒟ℯ𝓈𝓉𝒶𝒷━━━━━
+ Z________
+ ___Z_____
+ _X_______
+ __X______
+ ____X____
+ _____X___
+ ______X__
+ _______X_
𝒳ₗ━━━━━━━━━
+ ______XXX
𝒮𝓉𝒶𝒷━━━━━━━
- XXX___XXX
+ ___XXXXXX
+ ZZ_______
+ Z_Z______
+ ___ZZ____
+ ___Z_Z___
+ ______Z_Z
+ _______ZZ
𝒵ₗ━━━━━━━━━
+ Z__Z____Z
julia> cor_Z₂ * err_Z₁ * state # we are back to a good code state
𝒟ℯ𝓈𝓉𝒶𝒷━━━━━
+ Z________
+ ___Z_____
- _X_______
+ __X______
+ ____X____
+ _____X___
+ ______X__
+ _______X_
𝒳ₗ━━━━━━━━━
+ ______XXX
𝒮𝓉𝒶𝒷━━━━━━━
+ XXX___XXX
+ ___XXXXXX
+ ZZ_______
+ Z_Z______
+ ___ZZ____
+ ___Z_Z___
+ ______Z_Z
+ _______ZZ
𝒵ₗ━━━━━━━━━
+ Z__Z____Z
julia> bad_Z₆Z₉ = single_z(9,6) * single_z(9,9) # a different "correction" operation
+ _____Z__Z
julia> bad_Z₆Z₉ * err_Z₁ * state # the syndrome is trivial, but now we have a logical error
𝒟ℯ𝓈𝓉𝒶𝒷━━━━━
+ Z________
+ ___Z_____
+ _X_______
+ __X______
+ ____X____
- _____X___
+ ______X__
+ _______X_
𝒳ₗ━━━━━━━━━
- ______XXX
𝒮𝓉𝒶𝒷━━━━━━━
+ XXX___XXX
+ ___XXXXXX
+ ZZ_______
+ Z_Z______
+ ___ZZ____
+ ___Z_Z___
+ ______Z_Z
+ _______ZZ
𝒵ₗ━━━━━━━━━
+ Z__Z____Z
The success of cor_Z₂
and the failure of bad_Z₆Z₉
can be immediately seen through the fault matrix, as the wrong "correction" does not result in the same logical flips ad the error:
julia> O = faults_matrix(Shor9())
2×18 BitMatrix:
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1
1 0 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 0
julia> O * stab_to_gf2(err_Z₁)
2-element Vector{Int64}:
0
0
julia> O * stab_to_gf2(cor_Z₂)
2-element Vector{Int64}:
0
0
julia> O * stab_to_gf2(bad_Z₆Z₉)
2-element Vector{Int64}:
1
0
While its use in this situation is rather contrived, the fault matrix is incredibly useful when running large scale simulations in which we want a separate fast error sampling process, (e.g. with Pauli frames) and a syndrome decoding process, without coupling between them. We just gather all our syndrome measurement and logical observables from the Pauli frame simulations, and then use them with the fault matrix in the syndrome decoding simulation.
QuantumClifford.ECC.generalized_bicycle_codes
— FunctionImplemented in a package extension with Hecke. Check the docs for the Hecke extension
QuantumClifford.ECC.haah_cubic_codes
— FunctionImplemented in a package extension with Hecke. Check the docs for the Hecke extension
QuantumClifford.ECC.iscss
— MethodCheck if the code is CSS.
Return nothing
if unknown from the type.
QuantumClifford.ECC.isdegenerate
— FunctionCheck if the code is degenerate with respect to a given set of error or with respect to all "up to d physical-qubit" errors (defaulting to d=1).
julia> using QuantumClifford.ECC
julia> isdegenerate(Shor9(), [single_z(9,1), single_z(9,2)])
true
julia> isdegenerate(Shor9(), [single_z(9,1), single_x(9,1)])
false
julia> isdegenerate(Steane7(), 1)
false
julia> isdegenerate(Steane7(), 2)
true
QuantumClifford.ECC.naive_encoding_circuit
— MethodEncoding physical qubits into a larger logical code.
The initial physical qubits to be encoded have to be at indices n-k+1:n
.
Encoding circuits are not fault-tolerant, and thus should not be used in practice. Instead, you should measure the stabilizers of the code and the logical observables, thus projecting into the code space (which can be fault-tolerant).
The canonicalization operation performed on the code may permute the qubits (see canonicalize_gott!
). That permutation is corrected for with SWAP gates by default (controlled by the undoperm
keyword argument).
Based on (Cleve and Gottesman, 1997) and (Gottesman, 1997), however it seems the published algorithm has some errors. Consult the erratum, as well as the more recent (Grassl, 2002) and (Grassl, 2011), and be aware that this implementation also uses H instead of Z gates.
QuantumClifford.ECC.naive_syndrome_circuit
— FunctionGenerate the non-fault-tolerant stabilizer measurement cicuit for a given code instance or parity check tableau.
Use the ancillary_index
and bit_index
arguments to offset where the corresponding part the circuit starts.
Returns the circuit, the number of ancillary qubits that were added, and a list of bit indices that will store the measurement results.
See also: shor_syndrome_circuit
QuantumClifford.ECC.parity_checks
— FunctionParity check tableau of a code.
See also: parity_matrix_x
and parity_matrix_z
QuantumClifford.ECC.random_all_to_all_circuit_code
— FunctionRandom all-to-all Clifford circuit code (Brown and Fawzi, Jul 2013).
The code of n
qubits is generated by an all-to-all random Clifford circuit of ngates
gates that encodes a subset of qubits encode_qubits
into logical qubits.
Because of the random picking, the size of encode_qubits
is the only thing that matters for the code, referred to as k
.
See also: random_all_to_all_clifford_circuit
, CircuitCode
QuantumClifford.ECC.random_brickwork_circuit_code
— FunctionRandom brickwork Clifford circuit code (Brown and Fawzi, Jul 2013).
The code is generated by a brickwork random Clifford circuit of nlayers
layers that encodes a subset of qubits encode_qubits
into logical qubits.
See also: random_brickwork_clifford_circuit
, CircuitCode
QuantumClifford.ECC.shor_syndrome_circuit
— FunctionGenerate the Shor fault-tolerant stabilizer measurement cicuit for a given code instance or parity check tableau.
Use the ancillary_index
and bit_index
arguments to offset where the corresponding part the circuit starts. Ancillary qubits
Returns:
- The ancillary cat state preparation circuit.
- The Shor syndrome measurement circuit.
- The number of ancillary qubits that were added.
- The list of bit indices that store the final measurement results.
See also: naive_syndrome_circuit
QuantumClifford.ECC.two_block_group_algebra_codes
— FunctionImplemented in a package extension with Hecke. Check the docs for the Hecke extension
QuantumClifford.ECC.twobga_from_direct_product
— FunctionImplemented in a package extension with Oscar. Check the docs for the Oscar extension
QuantumClifford.ECC.twobga_from_fp_group
— FunctionImplemented in a package extension with Oscar. Check the docs for the Oscar extension
Implemented in an extension requiring Hecke.jl
QuantumCliffordHeckeExt.LPCode
— Typestruct LPCode <: QECCore.AbstractECC
Lifted product codes ((Panteleev and Kalachev, 2021), (Panteleev and Kalachev, Jun 2022))
A lifted product code is defined by the hypergraph product of a base matrices A
and the conjugate of another base matrix B'
. Here, the hypergraph product is taken over a group algebra, of which the base matrices are consisting.
The binary parity check matrix is obtained by applying repr
to each element of the matrix resulted from the hypergraph product, which is mathematically a linear map from each group algebra element to a binary matrix.
Constructors
Multiple constructors are available:
Two base matrices of group algebra elements.
Two lifted codes, whose base matrices are for quantum code construction.
Two base matrices of group elements, where each group element will be considered as a group algebra element by assigning a unit coefficient.
Two base matrices of integers, where each integer represent the shift of a cyclic permutation. The order of the cyclic permutation should be specified.
Below is a list of all constructors:
LPCode(A, B; GA, repr)
defined at /home/runner/work/QuantumClifford.jl/QuantumClifford.jl/ext/QuantumCliffordHeckeExt/lifted_product.jl:105
.
LPCode(c₁, c₂; GA, repr)
defined at /home/runner/work/QuantumClifford.jl/QuantumClifford.jl/ext/QuantumCliffordHeckeExt/lifted_product.jl:110
.
LPCode(A, B; GA)
defined at /home/runner/work/QuantumClifford.jl/QuantumClifford.jl/ext/QuantumCliffordHeckeExt/lifted_product.jl:118
.
LPCode(group_elem_array1, group_elem_array2; GA)
defined at /home/runner/work/QuantumClifford.jl/QuantumClifford.jl/ext/QuantumCliffordHeckeExt/lifted_product.jl:123
.
LPCode(shift_array1, shift_array2, l; GA)
defined at /home/runner/work/QuantumClifford.jl/QuantumClifford.jl/ext/QuantumCliffordHeckeExt/lifted_product.jl:128
.
Examples
A [[882, 24, d ≤ 24]] code from Appendix B of (Roffe et al., 2023). We use the 1st constructor to generate the code and check its length and dimension. During the construction, we do arithmetic operations to get the group algebra elements in base matrices A
and B
. Here x
is the generator of the group algebra, i.e., offset-1 cyclic permutation, and GA(1)
is the unit element.
julia> import Hecke: group_algebra, GF, abelian_group, gens; import LinearAlgebra: diagind; using QuantumClifford.ECC;
julia> l = 63; GA = group_algebra(GF(2), abelian_group(l)); x = gens(GA)[];
julia> A = zeros(GA, 7, 7);
julia> A[diagind(A)] .= x^27;
julia> A[diagind(A, -1)] .= x^54;
julia> A[diagind(A, 6)] .= x^54;
julia> A[diagind(A, -2)] .= GA(1);
julia> A[diagind(A, 5)] .= GA(1);
julia> B = reshape([1 + x + x^6], (1, 1));
julia> c1 = LPCode(A, B);
julia> code_n(c1), code_k(c1)
(882, 24)
A [[175, 19, d ≤ 10]] code from Eq. (18) in Appendix A of (Raveendran et al., 2022), following the 4th constructor.
julia> import Hecke; using QuantumClifford.ECC;
julia> base_matrix = [0 0 0 0; 0 1 2 5; 0 6 3 1]; l = 7;
julia> c2 = LPCode(base_matrix, l .- base_matrix', l);
julia> code_n(c2), code_k(c2)
(175, 19)
Code subfamilies and convenience constructors for them
- When the base matrices of the
LPCode
are 1×1, the code is called a two-block group-algebra codetwo_block_group_algebra_codes
. - When the base matrices of the
LPCode
are 1×1 and their elements are sums of cyclic permutations, the code is called a generalized bicycle codegeneralized_bicycle_codes
. - When the two matrices are adjoint to each other, the code is called a bicycle code
bicycle_codes
.
The representation function
We use the default representation function Hecke.representation_matrix
to convert a GF(2)
-group algebra element to a binary matrix. The default representation, provided by Hecke
, is the permutation representation.
We also accept a custom representation function as detailed in LiftedCode
.
See also: LiftedCode
, two_block_group_algebra_codes
, generalized_bicycle_codes
, bicycle_codes
, haah_cubic_codes
.
All fields:
A::Union{LinearAlgebra.Adjoint{<:Hecke.GroupAlgebraElem, <:Matrix{<:Hecke.GroupAlgebraElem}}, Matrix{<:Hecke.GroupAlgebraElem}}
: the first base matrix of the code, whose elements are in a group algebra.B::Union{LinearAlgebra.Adjoint{<:Hecke.GroupAlgebraElem, <:Matrix{<:Hecke.GroupAlgebraElem}}, Matrix{<:Hecke.GroupAlgebraElem}}
: the second base matrix of the code, whose elements are in the same group algebra asA
.GA::Hecke.GroupAlgebra
: the group algebra for which elements inA
andB
are from.repr::Function
: a function that converts a group algebra element to a binary matrix; default to be the permutation representation for GF(2)-algebra.
QuantumCliffordHeckeExt.LiftedCode
— Typestruct LiftedCode <: QECCore.AbstractCECC
Classical codes lifted over a group algebra, used for lifted product code construction ((Panteleev and Kalachev, 2021), (Panteleev and Kalachev, Jun 2022))
The parity-check matrix is constructed by applying repr
to each element of A
, which is mathematically a linear map from a group algebra element to a binary matrix. The size of the parity check matrix will enlarged with each element of A
being inflated into a matrix. The procedure is called a lift (Panteleev and Kalachev, Jun 2022).
Constructors
A lifted code can be constructed via the following approaches:
A matrix of group algebra elements.
A matrix of group elements, where a group element will be considered as a group algebra element by assigning a unit coefficient.
A matrix of integers, where each integer represent the shift of a cyclic permutation. The order of the cyclic permutation should be specified.
The default GA
is the group algebra of A[1, 1]
, the default representation repr
is the permutation representation.
Below is a list of all constructors:
LiftedCode(A; GA, repr)
defined at /home/runner/work/QuantumClifford.jl/QuantumClifford.jl/ext/QuantumCliffordHeckeExt/lifted.jl:58
.
LiftedCode(A; GA)
defined at /home/runner/work/QuantumClifford.jl/QuantumClifford.jl/ext/QuantumCliffordHeckeExt/lifted.jl:68
.
LiftedCode(group_elem_array; GA, repr)
defined at /home/runner/work/QuantumClifford.jl/QuantumClifford.jl/ext/QuantumCliffordHeckeExt/lifted.jl:74
.
LiftedCode(shift_array, l; GA)
defined at /home/runner/work/QuantumClifford.jl/QuantumClifford.jl/ext/QuantumCliffordHeckeExt/lifted.jl:87
.
The representation function repr
We use the default representation function Hecke.representation_matrix
to convert a GF(2)
-group algebra element to a binary matrix. The default representation, provided by Hecke
, is the permutation representation.
We also accept a custom representation function (the repr
field of the constructor). Whatever the representation, the matrix elements need to be convertible to Integers (e.g. permit lift(ZZ, ...)
). Such a customization would be useful to reduce the number of bits required by the code construction.
For example, if we use a D4 group for lifting, our default representation will be 8×8
permutation matrices, where 8 is the group's order. However, we can find a 4×4
matrix representation for the group, e.g. by using the typical 2×2
representation and converting it into binary representation by replacing "1" with the Pauli I, and "-1" with the Pauli X matrix.
See also: QuantumClifford.ECC.LPCode
.
All fields:
A::Union{LinearAlgebra.Adjoint{<:Hecke.GroupAlgebraElem, <:Matrix{<:Hecke.GroupAlgebraElem}}, Matrix{<:Hecke.GroupAlgebraElem}}
: the base matrix of the code, whose elements are in a group algebra.GA::Hecke.GroupAlgebra
: the group algebra for which elements inA
are from.repr::Function
: a function that converts a group algebra element to a binary matrix; default to be the permutation representation for GF(2)-algebra.
QuantumClifford.ECC.bicycle_codes
— Methodbicycle_codes(a_shifts::Array{Int64}, l::Int64) -> Any
Bicycle codes are a special case of generalized bicycle codes, where a
and b
are conjugate to each other. The order of the cyclic group is l
, and the shifts a_shifts
and b_shifts
are reverse to each other. Thus you need to provide only the a_shifts
and the rest of the conversions and conjugations are taken care of.
See also: two_block_group_algebra_codes
, generalized_bicycle_codes
, haah_cubic_codes
.
QuantumClifford.ECC.generalized_bicycle_codes
— Methodgeneralized_bicycle_codes(
a_shifts::Array{Int64},
b_shifts::Array{Int64},
l::Int64
) -> Any
Generalized bicycle codes, which are a special case of abelian 2GBA codes (and therefore of lifted product codes). Here the group is chosen as the cyclic group of order l
, and the base matrices a
and b
are the sum of the group algebra elements corresponding to the shifts a_shifts
and b_shifts
.
Behind the scenes, the shifts are converted to the corresponding group algebra elements and piped to two_block_group_algebra_codes
.
See also: two_block_group_algebra_codes
, bicycle_codes
.
Examples
A [[254, 28, 14 ≤ d ≤ 20]] code from (A1) in Appendix B of (Panteleev and Kalachev, 2021).
julia> import Hecke; using QuantumClifford.ECC
julia> c = generalized_bicycle_codes([0, 15, 20, 28, 66], [0, 58, 59, 100, 121], 127);
julia> code_n(c), code_k(c)
(254, 28)
An [[70, 8, 10]] abelian 2BGA code from Table 1 of (Lin and Pryadko, 2024), with cyclic group of order l = 35
, illustrates that abelian 2BGA codes can be viewed as GB codes.
julia> import Hecke; using QuantumClifford.ECC
julia> l = 35;
julia> c1 = generalized_bicycle_codes([0, 15, 16, 18], [0, 1, 24, 27], l);
julia> code_n(c1), code_k(c1)
(70, 8)
QuantumClifford.ECC.haah_cubic_codes
— Methodhaah_cubic_codes(
a_shifts::Array{Int64},
b_shifts::Array{Int64},
l::Int64
) -> Any
Haah’s cubic codes (Haah, 2011) can be viewed as generalized bicycle (GB) codes with the group G = Cₗ × Cₗ × Cₗ
, where l
denotes the lattice size. In particular, a GB code with the group G = ℤ₃ˣ³
corresponds to a cubic code.
Behind the scenes, this function is just a simple shortcut for preparing the group G
, before piping the arguments to generalized_bicycle_codes
.
The ECC Zoo has an entry for this family.
See also: bicycle_codes
, generalized_bicycle_codes
, two_block_group_algebra_codes
.
Examples
julia> import Hecke; using QuantumClifford.ECC;
julia> c = haah_cubic_codes([0, 15, 20, 28, 66], [0, 58, 59, 100, 121], 6);
julia> code_n(c), code_k(c)
(432, 8)
QuantumClifford.ECC.two_block_group_algebra_codes
— Methodtwo_block_group_algebra_codes(
a::Hecke.GroupAlgebraElem,
b::Hecke.GroupAlgebraElem
) -> Any
Two-block group algebra (2BGA) codes, which are a special case of lifted product codes from two group algebra elements a
and b
, used as 1×1
base matrices. To build them, you pick a group and specific generators for that group, then you pick two polynomials made of the group generators, and then, behind the scenes, these two polynomials a
and b
are piped to the lifted product code constructor as the elements of 1×1
matrices.
See also: QuantumClifford.ECC.LPCode
, generalized_bicycle_codes
, bicycle_codes
, haah_cubic_codes
.
Examples of 2BGA code subfamilies
C₄ x C₂
Here is an example of a [[56, 28, 2]] 2BGA code from Table 2 of (Lin and Pryadko, 2024) build out of polymonials of generators of the direct product C₄ × C₂
.
julia> import Hecke: group_algebra, GF, abelian_group, gens; using QuantumClifford.ECC;
julia> GA = group_algebra(GF(2), abelian_group([14,2]));
julia> x, s = gens(GA);
julia> A = 1 + x^7;
julia> B = 1 + x^7 + s + x^8 + s*x^7 + x;
julia> c = two_block_group_algebra_codes(A,B);
julia> import HiGHS;
julia> code_n(c), code_k(c), distance(c, DistanceMIPAlgorithm(solver=HiGHS))
(56, 28, 2)
Bivariate Bicycle codes
Bivariate Bicycle codes are a class of Abelian 2BGA codes formed by the direct product of two cyclic groups ℤₗ × ℤₘ
. The parameters l
and m
represent the orders of the first and second cyclic groups, respectively.
The ECC Zoo has an entry for this family.
A [[756, 16, ≤ 34]] code from Table 3 of (Bravyi et al., 2024):
julia> import Hecke: group_algebra, GF, abelian_group, gens; using QuantumClifford.ECC;
julia> l=21; m=18;
julia> GA = group_algebra(GF(2), abelian_group([l, m]));
julia> x, y = gens(GA);
julia> A = x^3 + y^10 + y^17;
julia> B = y^5 + x^3 + x^19;
julia> c = two_block_group_algebra_codes(A,B);
julia> code_n(c), code_k(c)
(756, 16)
Multivariate Bicycle code
The group algebra of the qubit multivariate bicycle (MB) code with r variables is 𝔽₂[𝐺ᵣ]
, where 𝐺ᵣ = ℤ/l₁ × ℤ/l₂ × ... × ℤ/lᵣ
.
A [[48, 4, 6]] Weight-6 TB-QLDPC code from Appendix A Table 2 of (Voss et al., 2024).
julia> import Hecke: group_algebra, GF, abelian_group, gens; using QuantumClifford.ECC;
julia> l=4; m=6;
julia> GA = group_algebra(GF(2), abelian_group([l, m]));
julia> x, y = gens(GA);
julia> z = x*y;
julia> A = x^3 + y^5;
julia> B = x + z^5 + y^5 + y^2;
julia> c = two_block_group_algebra_codes(A, B);
julia> import HiGHS
julia> code_n(c), code_k(c), distance(c, DistanceMIPAlgorithm(solver=HiGHS))
(48, 4, 6)
Coprime Bivariate Bicycle code
The coprime bivariate bicycle (BB) codes are defined by two polynomials 𝑎(𝑥,𝑦)
and 𝑏(𝑥,𝑦)
, where 𝑙
and 𝑚
are coprime, and can be expressed as univariate polynomials 𝑎(𝜋)
and 𝑏(𝜋)
, with generator 𝜋 = 𝑥𝑦
. They can be viewed as a special case of Lifted Product construction based on abelian group ℤₗ x ℤₘ
where ℤⱼ
cyclic group of order j
.
[[108, 12, 6]] coprime-bivariate bicycle (BB) code from Table 2 of (Wang and Mueller, 2024).
julia> import Hecke: group_algebra, GF, abelian_group, gens; using QuantumClifford.ECC;
julia> l=2; m=27;
julia> GA = group_algebra(GF(2), abelian_group([l*m]));
julia> 𝜋 = gens(GA)[1];
julia> A = 𝜋^2 + 𝜋^5 + 𝜋^44;
julia> B = 𝜋^8 + 𝜋^14 + 𝜋^47;
julia> c = two_block_group_algebra_codes(A, B);
julia> import HiGHS
julia> code_n(c), code_k(c), distance(c, DistanceMIPAlgorithm(solver=HiGHS))
(108, 12, 6)
Implemented in an extension requiring Oscar.jl
QuantumClifford.ECC.twobga_from_direct_product
— Methodtwobga_from_direct_product(
a_elts::Vector{Hecke.GroupAlgebraElem{Nemo.FqFieldElem, Hecke.GroupAlgebra{Nemo.FqFieldElem, Oscar.DirectProductGroup, Oscar.BasicGAPGroupElem{Oscar.DirectProductGroup}}}},
b_elts::Vector{Hecke.GroupAlgebraElem{Nemo.FqFieldElem, Hecke.GroupAlgebra{Nemo.FqFieldElem, Oscar.DirectProductGroup, Oscar.BasicGAPGroupElem{Oscar.DirectProductGroup}}}},
F2G::Hecke.GroupAlgebra{Nemo.FqFieldElem, Oscar.DirectProductGroup, Oscar.BasicGAPGroupElem{Oscar.DirectProductGroup}}
) -> Any
Constructing two block group algebra codes by specifying the direct product to be used. See also the more general two_block_group_algebra_codes
.
Two block group algebra codes are constructed by choosing a group (and specific generators), then choosing two polynomials made out of these generators, then piping these two polynomials as the elements of 1×1
matrices to the lifted product code constructors.
The Hecke library, for which we already have an extension, provides for a fairly easy way to construct such polynomials for many abelian and small groups. See two_block_group_algebra_codes
for those capabilities.
However, more esoteric groups can be specified as the direct product of other groups. To support arbitrary direct products we use Oscar, which builds upon Hecke. Oscar supports the direct product operation between two or more arbitrary general groups, including non-abelian groups such as alternating_group
, dihedral_group
, symmetric_group
, and even arbitrary finitely presented groups (e.g., free_group
). This capability is not available in Hecke.jl
. The 2BGA codes discovered in (Lin and Pryadko, 2024) rely on direct products of two or more general groups, which necessitate the use of Oscar.direct_product
.
This particular function is nothing more than a simple wrapper that takes care of argument conversions. Of note, the polynomials here are given as lists of monomials.
Of course, if you are comfortable with Oscar, you can use two_block_group_algebra_codes
directly.
See also: two_block_group_algebra_codes
, twobga_from_fp_group
Examples
The [[56, 28, 2]] abelian 2BGA code from Appendix C, Table II in (Lin and Pryadko, 2024) can be constructed using the direct product of two cyclic groups. Specifically, the group C₂₈
of order l = 28
can be represented as C₁₄ × C₂
, where the first group has order m = 14
and the second group has order n = 2
.
julia> import Oscar: cyclic_group, small_group_identification, describe, order
julia> import Hecke: gens, quo, group_algebra, GF, one, direct_product, sub
julia> using QuantumClifford, QuantumClifford.ECC
julia> m = 14; n = 2;
julia> C₁₄ = cyclic_group(m);
julia> C₂ = cyclic_group(n);
julia> G = direct_product(C₁₄, C₂);
julia> GA = group_algebra(GF(2), G);
julia> x, s = gens(GA)[1], gens(GA)[3];
julia> a = [one(GA), x^7];
julia> b = [one(GA), x^7, s, x^8, s * x^7, x];
julia> c = twobga_from_direct_product(a, b, GA);
julia> order(G)
28
julia> code_n(c), code_k(c)
(56, 28)
julia> describe(G), small_group_identification(G)
("C14 x C2", (28, 4))
When using the direct product, there isn't necessarily a unique set of generators. It is essential to verify that Oscar is providing you with the generators you expect, e.g. for a cycling group that you have the presentation Cₘ = ⟨x, s | xᵐ = s² = xsx⁻¹s⁻¹ = 1⟩
. For situations where the generators provided by Oscar are not the ones you want, you can also use twobga_from_fp_group
where you specify the group presentation directly.
As a verification that you have the correct generators, Oscar.sub
can be used to determine if H
is a subgroup of G
and to confirm that both C₁₄
and C₂
are subgroups of C₂₈
.
julia> order(gens(G)[1])
14
julia> order(gens(G)[3])
2
julia> x^14 == s^2 == x * s * x^-1 * s^-1
true
julia> H, _ = sub(G, [gens(G)[1], gens(G)[3]]);
julia> H == G
true
QuantumClifford.ECC.twobga_from_fp_group
— Methodtwobga_from_fp_group(
a_elts::Vector{Oscar.FPGroupElem},
b_elts::Vector{Oscar.FPGroupElem},
F2G::Hecke.GroupAlgebra{Nemo.FqFieldElem, Oscar.FPGroup, Oscar.FPGroupElem}
) -> Any
Constructing two block group algebra codes by specifying the group presentation.
Two block group algebra codes are constructed by choosing a group (and specific generators), then choosing two polynomials made out of these generators, then piping these two polynomials as the elements of 1×1
matrices to the lifted product code constructors.
The Hecke library, for which we already have an extension, provides for a fairly easy way to construct such polynomials for many abelian and small groups. See two_block_group_algebra_codes
for those capabilities.
However, more esoteric groups are usually specified by a group presentation ⟨S | R⟩
, where S
is a set of generators and R
is the relations those generators obey. To support arbitrary group presentations we use Oscar, which builds upon Hecke. We use Oscar.free_group
and quo
in order to first prepare the free group generated by S
, and then the group obeying also the relations R
, i.e. the ⟨S | R⟩
presentation.
After that point we proceed as usual, creating two polynomials of generators and piping them to two_block_group_algebra_codes
.
This particular function is nothing more than a simple wrapper that takes care of argument conversions. Of note, the polynomials here are given as lists of monomials.
Of course, if you are comfortable with Oscar, you can use two_block_group_algebra_codes
directly.
See also: two_block_group_algebra_codes
, twobga_from_direct_product
Examples
The [[96, 12, 10]] 2BGA code from Table I in (Lin and Pryadko, 2024) has the group presentation ⟨r, s | s⁶ = r⁸ = r⁻¹srs = 1⟩
(the group C₂ × (C₃ ⋉ C₈)
).
julia> import Oscar: free_group, small_group_identification, describe, order
julia> import Hecke: gens, quo, group_algebra, GF, one
julia> using QuantumClifford, QuantumClifford.ECC
julia> F = free_group(["r", "s"]);
julia> r, s = gens(F); # generators
julia> G, = quo(F, [s^6, r^8, r^(-1) * s * r * s]); # relations
julia> GA = group_algebra(GF(2), G);
julia> r, s = gens(G);
julia> a = [one(G), r, s^3 * r^2, s^2 * r^3];
julia> b = [one(G), r, s^4 * r^6, s^5 * r^3];
julia> c = twobga_from_fp_group(a, b, GA);
julia> order(G)
48
julia> code_n(c), code_k(c)
(96, 12)
julia> describe(G), small_group_identification(G)
("C2 x (C3 : C8)", (48, 9))
Cyclic Groups
Cyclic groups with specific group presentations, given by Cₘ = ⟨x, s | xᵐ = s² = xsx⁻¹s⁻¹ = 1⟩
, where the order is 2m
as seen in Table II of (Lin and Pryadko, 2024).
The [[56, 28, 2]] abelian 2BGA code from Appendix C, Table II in (Lin and Pryadko, 2024) is constructed using the group presentation ⟨x, s | xs = sx, xᵐ = s² = 1⟩
(the cyclic group C₂₈ = C₁₄ × C₂
).
julia> m = 14;
julia> F = free_group(["x", "s"]);
julia> x, s = gens(F); # generators
julia> G, = quo(F, [x^m, s^2, x * s * x^-1 * s^-1]); # relations
julia> GA = group_algebra(GF(2), G);
julia> x, s = gens(G);
julia> a = [one(G), x^7];
julia> b = [one(G), x^7, s, x^8, s * x^7, x];
julia> c = twobga_from_fp_group(a, b, GA);
julia> order(G)
28
julia> code_n(c), code_k(c)
(56, 28)
julia> describe(G), small_group_identification(G)
("C14 x C2", (28, 4))
Dihedral Groups
Dihedral (non-abelian) groups with group presentations given by Dₘ = ⟨r, s | rᵐ = s² = (rs)² = 1⟩
, where the order is 2m
.
The [[24, 8, 3]] 2BGA code from Appendix C, Table III in (Lin and Pryadko, 2024) is constructed by specifying a group presentation below (giving the group D₆ = C₆ ⋉ C₂
).
julia> m = 6;
julia> F = free_group(["r", "s"]);
julia> r, s = gens(F); # generators
julia> G, = quo(F, [r^m, s^2, (r*s)^2]); # relations
julia> GA = group_algebra(GF(2), G);
julia> r, s = gens(G);
julia> a = [one(G), r^4];
julia> b = [one(G), s*r^4, r^3, r^4, s*r^2, r];
julia> c = twobga_from_fp_group(a, b, GA);
julia> order(G)
12
julia> code_n(c), code_k(c)
(24, 8)
julia> describe(G), small_group_identification(G)
("D12", (12, 4))
Implemented in an extension requiring JuMP.jl
QECCore.distance
— Methoddistance(
code::QECCore.AbstractECC,
alg::QuantumClifford.ECC.DistanceMIPAlgorithm
) -> Any
Compute the distance of a code using mixed integer programming. See QuantumClifford.ECC.DistanceMIPAlgorithm
for configuration options.
Computes the minimum Hamming weight of a binary vector x
by solving an mixed integer program (MIP) that satisfies the following constraints:
- $\text{stab} \cdot x \equiv 0 \pmod{2}$: The binary vector
x
must have an
even overlap with each X
-check of the stabilizer binary representation stab
.
- $\text{logicOp} \cdot x \equiv 1 \pmod{2}$: The binary vector
x
must have
an odd overlap with logical-X
operator logicOp
on the i
-th logical qubit.
Specifically, it calculates the minimum Hamming weight $d_{Z}$ for the Z
-type logical operator. The minimum distance for X
-type logical operators is the same.
Background on Minimum Distance
For classical codes, the minimum distance, which measures a code's error-correcting capability, is equivalent to its minimum weight. This can be computed by generating all possible codewords from combinations of the generator matrix rows, calculating their weights, and finding the smallest. While accurate, this method takes exponential time. Vardy (Vardy, 1997) demonstrated that computing the minimum distance is NP-hard, and the corresponding decision problem is NP-complete, making polynomial-time algorithms unlikely.
For quantum codes, classical intuition does not always apply. The minimum distance is given by the minimum weight of a non-trivial logical operator. This is generally unrelated to the minimum distance of the corresponding stabilizer code when viewed as a classical, additive code. White and Grassl (White and Grassl, 2006) proposed mapping quantum codes to higher-dimensional classical linear codes. This mapping allows the minimum distance of the quantum additive code to be inferred from that of the classical linear code but increases parameters from n
to 3n
and d
to 2d
, adding complexity. Furthermore, once a minimal weight vector is identified, it is essential to verify whether it belongs to the Pauli group 𝒫ₙ
over n
qubits (Sabo, 2022).
Additionally, to illustrate how classical intuition can be misleading in this context, consider that the [[7, 1, 3]]
Steane code has a minimum distance of three, despite all its elements having a weight of four. This discrepancy occurs because stabilizer codes are defined by parity-check matrices, while their minimum distances are determined by the dual (Sabo, 2022).
julia> using QuantumClifford, QuantumClifford.ECC
julia> c = parity_checks(Steane7());
julia> stab_to_gf2(c)
6×14 Matrix{Bool}:
0 0 0 1 1 1 1 0 0 0 0 0 0 0
0 1 1 0 0 1 1 0 0 0 0 0 0 0
1 0 1 0 1 0 1 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 1 1 1 1
0 0 0 0 0 0 0 0 1 1 0 0 1 1
0 0 0 0 0 0 0 1 0 1 0 1 0 1
julia> minimum(sum(stab_to_gf2(c), dims=2))
4
julia> distance(Steane7())
3
The minimum distance problem for quantum codes is NP-hard, and this hardness extends to multiplicative and additive approximations, even when restricted to stabilizer or CSS codes, with the result established through a reduction from classical problems in the CWS framework using a 4-cycle free graph (Kapshikar and Kundu, 2023). Despite this, methods that improve on brute-force approaches are actively explored.
For a more in-depth background on minimum distance, see Chapter 3 of (Sabo, 2022).
Mixed Integer Programming
The MIP minimizes the Hamming weight w(x)
, defined as the number of nonzero elements in x
, while satisfying the constraints:
\[\begin{aligned} \text{Minimize} \quad & w(x) = \sum_{i=1}^n x_i, \\ \text{subject to} \quad & \text{stab} \cdot x \equiv 0 \pmod{2}, \\ & \text{logicOp} \cdot x \equiv 1 \pmod{2}, \\ & x_i \in \{0, 1\} \quad \text{for all } i. \end{aligned}\]
Here:
- $\text{stab}$ is the binary matrix representing the stabilizer group.
- $\text{logicOp}$ is the binary vector representing the logical-
X
operator. x
is the binary vector (decision variable) being optimized.
The optimal solution $w_{i}$ for each logical-X
operator corresponds to the minimum weight of a Pauli Z
-type operator satisfying the above conditions. The Z
-type distance is given by:
\[\begin{aligned} d_Z = \min(w_1, w_2, \dots, w_k), \end{aligned}\]
where k
is the number of logical qubits.
Example
A [[40, 8, 5]] 2BGA code with the minimum distance of 5 from Table 2 of (Lin and Pryadko, 2024).
julia> import Hecke: group_algebra, GF, abelian_group, gens; import HiGHS; import JuMP;
julia> using QuantumClifford, QuantumClifford.ECC
julia> l = 10; m = 2;
julia> GA = group_algebra(GF(2), abelian_group([l,m]));
julia> x, s = gens(GA);
julia> A = 1 + x^6;
julia> B = 1 + x^5 + s + x^6 + x + s*x^2;
julia> c = two_block_group_algebra_codes(A,B);
julia> code_n(c), code_k(c), distance(c, DistanceMIPAlgorithm(solver=HiGHS))
(40, 8, 5)
A [[48, 6, 8]] GB code with the minimum distance of 8 from (A3) in Appendix B of (Panteleev and Kalachev, 2021).
julia> l = 24;
julia> c1 = generalized_bicycle_codes([0, 2, 8, 15], [0, 2, 12, 17], l);
julia> code_n(c1), code_k(c1), distance(c1, DistanceMIPAlgorithm(solver=HiGHS))
(48, 6, 8)
Applications
Mixed-integer programming (MIP) is applied in quantum error correction, notably for decoding and minimum distance computation. Some applications are as follows:
The first usecase of the MIP approach was the code capacity Most Likely Error (MLE) decoder for color codes introduced in (Landahl et al., 2011).
For all quantum LDPC codes presented in (Panteleev and Kalachev, 2021), the lower and upper bounds on the minimum distance was obtained by reduction to a mixed integer linear program and using the GNU Linear Programming Kit ((Makhorin, 2008)).
For all the Bivariate Bicycle (BB) codes presented in (Bravyi et al., 2024), the code distance was calculated using the mixed integer programming approach.
(Lacroix et al., 2024) developed a MLE decoder that finds the most likely chain of Pauli errors given the observed error syndrome by solving a mixed-integer program using
HiGHS
package ((Huangfu and Hall, 2018)).(Cain et al., 2025) formulate maximum-likelihood decoding as a mixed-integer program maximizing $\prod_{j=1}^M p_j^{E_j}(1-p_j)^{1-E_j}$ (where binary variables $E_j \in {0,1}$ indicate error occurrence) subject to syndrome constraints, solved optimally via MIP solvers despite its NP-hard complexity.