# We'll implement Pauli propagation and enumerate error patterns up to two errors for the given circuit.
from itertools import product, combinations
# Represent Pauli on n qubits via (x,z) bit vectors and a phase integer m in {0,1,2,3} corresponding to i^m factor
# We'll ignore phase for commutation tests except when needed; but to keep fullness we can track phase.
n = 5 # qubits 0..4, with 4 as ancilla
class Pauli:
def __init__(self, x=None, z=None, m=0):
if x is None:
x = [0]*n
if z is None:
z = [0]*n
self.x = x[:] # list of 0/1
self.z = z[:] # list of 0/1
self.m = m % 4
def copy(self):
return Pauli(self.x, self.z, self.m)
def __mul__(self, other):
# Multiply two Paulis, returning new Pauli; track phase using symplectic product
assert isinstance(other, Pauli)
x = [(a ^ b) for a,b in zip(self.x, other.x)]
z = [(a ^ b) for a,b in zip(self.z, other.z)]
# phase: i^{m1+m2 + 2*symp(self,other)}? Actually multiplication rule:
# P(a) P(b) = i^{[a,b]} P(a+b), where [a,b] in Z_4 computed from symplectic product.
# More precisely, the phase exponent increment is sum over j of (a_z[j]*b_x[j] - a_x[j]*b_z[j]) mod 4.
inc = 0
for j in range(n):
inc += self.z[j]*other.x[j]
inc -= self.x[j]*other.z[j]
m = (self.m + other.m + 2*inc) % 4 # factor of i^{2*inc} = (-1)^{inc}; but to keep i exponents consistent use 2*inc mod 4
# Actually correct formula is: P1*P2 = i^{omega(P1,P2)} P3 with omega in Z, but using 2*inc mod 4 yields correct +/- sign; ignoring global i-phase suffices.
return Pauli(x,z,m)
def commutes_with(self, other):
# Compute commutation: P1 P2 = (-1)^{symp} P2 P1
s = 0
for j in range(n):
s ^= (self.x[j] & other.z[j]) ^ (self.z[j] & other.x[j])
return (s == 0)
def conj_H(self, q):
# Conjugate by H on qubit q: X<->Z, Y -> -Y (phase)
xq, zq = self.x[q], self.z[q]
self.x[q], self.z[q] = zq, xq
if xq and zq: # Y
self.m = (-self.m) % 4 # a sign flip; but phase tracking isn't critical
return self
def conj_S(self, q):
# Phase gate if needed; not used
xq, zq = self.x[q], self.z[q]
if xq and zq:
self.m = (self.m + 2) % 4 # Y -> -Y
self.z[q] ^= self.x[q]
return self
def conj_CNOT(self, c, t):
# Conjugation rules for Pauli by CNOT(c,t):
# X_c -> X_c X_t; Z_c -> Z_c; X_t -> X_t; Z_t -> Z_c Z_t
# In symplectic representation for conjugation: x_t ^= x_c; z_c ^= z_t
# Wait careful: Conjugation mapping of Paulis P under U: U P U^
# Standard mapping on binary vectors: x' = x (with x_t ^= x_c), z' = z (with z_c ^= z_t).
xc, zc = self.x[c], self.z[c]
xt, zt = self.x[t], self.z[t]
# Update Z on control: z_c ^= z_t
self.z[c] ^= zt
# Update X on target: x_t ^= x_c
self.x[t] ^= xc
# Phase handling for Y terms may introduce sign; ignore phase
return self
def has_X_or_Y(self, q):
return self.x[q] == 1
def has_Z_or_Y(self, q):
return self.z[q] == 1
def __repr__(self):
# Return string like 'XIZYI'
s = ''
for j in range(n):
if self.x[j]==0 and self.z[j]==0:
s += 'I'
elif self.x[j]==1 and self.z[j]==0:
s += 'X'
elif self.x[j]==0 and self.z[j]==1:
s += 'Z'
else:
s += 'Y'
return s
# Helper: create single-qubit Pauli on q with label in {'I','X','Y','Z'}
def single_pauli(q,label):
x = [0]*n; z=[0]*n
if label=='I':
pass
elif label=='X':
x[q]=1
elif label=='Z':
z[q]=1
elif label=='Y':
x[q]=1; z[q]=1
else:
raise ValueError
return Pauli(x,z)
# Build two-qubit Paulis on given pair (a,b), excluding II
labels = ['I','X','Y','Z']
def twoqubit_paulis(a,b):
out=[]
for la in labels:
for lb in labels:
if la=='I' and lb=='I':
continue
P = single_pauli(a,la)
P = P * single_pauli(b,lb)
out.append((la+lb,P))
return out
# Circuit description: H1, CNOT12, CNOT10, CNOT23, CNOT34, CNOT04
# We'll need to propagate an error inserted after each CNOT to the measurement end (i.e., through later CNOTs)
# We'll implement a function that, given an insertion location index l in 1..5 and a two-qubit Pauli error on the two qubits of that gate,
# returns the propagated Pauli at the end of the circuit.
# Prepare list of gates as functions to conjugate; We'll store as tuples ('H',q) or ('CNOT',c,t)
GATES = [('H',1), ('CNOT',1,2), ('CNOT',1,0), ('CNOT',2,3), ('CNOT',3,4), ('CNOT',0,4)]
# Map insertion locations to gate indices: after gate index i (1..5) - corresponds to after GATES[i]
# Actually GATES length is 6; but only after the CNOTs at indices 1..5 (1-based). We'll compute properly with 0-based indexing.
CNOT_indices = [1,2,3,4,5] # indices in GATES list where CNOTs occur
# For each insertion location l, return the list of subsequent gates (with higher index)
def propagate_error_from_location(l_index_zero_based, error_P):
# Conjugate error_P by subsequent gates: for gate k = l_index+1 .. len(GATES)-1
P = error_P.copy()
for k in range(l_index_zero_based+1, len(GATES)):
g = GATES[k]
if g[0]=='H':
P.conj_H(g[1])
elif g[0]=='CNOT':
P.conj_CNOT(g[1], g[2])
return P
# Build convenience: mapping from CNOT location index 0..4 to (c,t) qubit indices, and to the gate index in GATES
cnots = [(1,2),(1,0),(2,3),(3,4),(0,4)]
cnots_gate_index = [i for i,g in enumerate(GATES) if g[0]=='CNOT']
assert cnots == [(GATES[i][1], GATES[i][2]) for i in cnots_gate_index]
# Define stabilizers on data qubits 0..3: S1=XXXX, S2=ZZZZ
S1 = Pauli([1,1,1,1,0],[0,0,0,0,0]) # XXXXII on 5 qubits
S2 = Pauli([0,0,0,0,0],[1,1,1,1,0]) # ZZZZ I
# Define logicals: X_A = X I X I; X_B = X X I I; Z_A = Z Z I I; Z_B = Z I Z I, on qubits 0..3
X_A = Pauli([1,0,1,0,0],[0,0,0,0,0])
X_B = Pauli([1,1,0,0,0],[0,0,0,0,0])
Z_A = Pauli([0,0,0,0,0],[1,1,0,0,0])
Z_B = Pauli([0,0,0,0,0],[1,0,1,0,0])
# Helper: Check if Pauli on data qubits is in normalizer of S (commutes with both S1,S2)
def commutes_with_stabilizers(P):
return P.commutes_with(S1) and P.commutes_with(S2)
# Helper: reduce Pauli on data qubits modulo stabilizer generators S1,S2 to a canonical representative in normalizer quotient.
# We'll decide if it includes logical X factors. We'll do by solving for decomposition in basis of {S1,S2,X_A,X_B,Z_A,Z_B} ignoring overall phase.
# We'll find booleans xA,xB,zA,zB,s1,s2 such that P ~ S1^s1 S2^s2 X_A^xA X_B^xB Z_A^zA Z_B^zB
# We'll do Gaussian elimination on binary symplectic representation restricted to data qubits 0..3.
import itertools
# Build basis list for data part: using 4-qubit Paulis represented as Pauli on 5 qubits but ancilla factor must be I on 4.
BASIS = [S1, S2, X_A, X_B, Z_A, Z_B]
# We can try brute force: test all 64 combinations of s1,s2,xA,xB,zA,zB to see which equals P up to a phase on data qubits.
def classify_logical_effect(P):
# P is a Pauli on 5 qubits; consider its data part P_data with ancilla I.
# Construct a P0 that is I on ancilla and equals P on data; We'll ignore ancilla content for classification.
# To compare equality up to phase, we can multiply P by candidate^-1 and check if the result acts trivially on data and possibly Z on ancilla only; but we only care about data equality.
# We'll compare x,z on qubits 0..3.
px = P.x[:4]
pz = P.z[:4]
# Brute force search
for s1 in [0,1]:
for s2 in [0,1]:
for xA_ in [0,1]:
for xB_ in [0,1]:
for zA_ in [0,1]:
for zB_ in [0,1]:
Q = Pauli([0]*5,[0]*5)
for gen,exp in zip(BASIS, [s1,s2,xA_,xB_,zA_,zB_]):
if exp:
Q = Q * gen
# Compare data part
if Q.x[:4]==px and Q.z[:4]==pz:
return dict(s1=s1,s2=s2,xA=xA_,xB=xB_,zA=zA_,zB=zB_)
# If not found, P is not in the subgroup generated by basis; could be not in normalizer
return None
# Acceptance tests
def ancilla_accepted(P):
# Accept if no X or Y on ancilla (i.e., x4==0)
return P.x[4]==0
def stabilizers_accepted(P):n # Accept if commutes with S1,S2
return commutes_with_stabilizers(P)
# Enumerate single errors: verify none are accepted-harmful
single_detected=0
single_accepted=0
single_accepted_harmful=0
for loc_index, (c,t) in enumerate(cnots):
# Gate index in GATES for this CNOT
gate_idx = cnots_gate_index[loc_index]
# Two-qubit paulis on (c,t)
for (lab, Err) in twoqubit_paulis(c,t):
P = propagate_error_from_location(gate_idx, Err)
anc_ok = ancilla_accepted(P)
stab_ok = commutes_with_stabilizers(P)
accepted = anc_ok and stab_ok
if accepted:
single_accepted += 1
cls = classify_logical_effect(P)
# Harmful: includes xA or xB = 1 (Y counted as X times Z)
harmful = (cls is not None) and ((cls['xA']%2==1) or (cls['xB']%2==1))
if harmful:
single_accepted_harmful += 1
else:
single_detected += 1
single_detected, single_accepted, single_accepted_harmful