From d2bbedc0db4a82db6eb6d97cc78e040548a832dc Mon Sep 17 00:00:00 2001 From: Noah Shutty Date: Sun, 17 Jul 2022 23:08:05 -0700 Subject: [PATCH] Added function to compute minimum distance of resulting codes. --- mindist.py | 152 +++++++++++++++++++++++++++++++++++++++++++++++ requirements.txt | 7 +++ solver.py | 30 +++++----- 3 files changed, 175 insertions(+), 14 deletions(-) create mode 100644 mindist.py diff --git a/mindist.py b/mindist.py new file mode 100644 index 0000000..f01d638 --- /dev/null +++ b/mindist.py @@ -0,0 +1,152 @@ +import numpy as np +import galois +import z3 + + +def xzdual(row): + assert len(row)%2==0 + num_qubits = len(row)//2 + return row[num_qubits:] + row[:num_qubits] + +def matrix_product(A, B, addfunc=lambda x,y: x+y): + m1 = len(A) + n1 = len(A[0]) + m2 = len(B) + n2 = len(B[0]) + assert n1 == m2 + C = [ + [# n2 columns + 0 for j in range(n2)] + for i in range(m1) # m1 rows + ] + for i in range(m1): + for j in range(n2): + for k in range(n1): + C[i][j] = addfunc(C[i][j], A[i][k] * B[k][j]) + return C + +def FastOr(a, b): + a_is_const = a is True or a is False + b_is_const = b is True or b is False + if a_is_const or b_is_const: + if a is True or b is True: + return True + if a is False: + return b + if b is False: + return a + return z3.Or(a, b) + +def FastAnd(a,b): + a_is_const = a is True or a is False + b_is_const = b is True or b is False + if a_is_const or b_is_const: + if a is True: + return b + if b is True: + return a + return False + else: + return z3.And(a, b) + +def FastNot(a): + a_is_const = a is True or a is False + if a_is_const: + return not a + else: + return z3.Not(a) + +def FastXor(a, b): + a_is_const = a is True or a is False + b_is_const = b is True or b is False + if a_is_const or b_is_const: + if a is False: + return b + if b is False: + return a + if a is True: + if b is False: + return True + if b is True: + return False + return FastNot(b) + if b is True: + if a is False: + return True + if a is True: + return False + return FastNot(a) + + return z3.Xor(a, b) + +def FastXorSum(li): + const = False + nonconsts = set(range(len(li))) + for i in range(len(li)): + val = li[i] + val_is_const = val is True or val is False + if val_is_const: + nonconsts.remove(i) + const = const != val + if not nonconsts: + return const + else: + nonconsts = list(nonconsts) + nonconst_sum = li[nonconsts[0]] + for i in nonconsts[1:]: + nonconst_sum = FastXor(nonconst_sum, li[i]) + if const: + return FastNot(nonconst_sum) + else: + return nonconst_sum + +def minimum_distance(num_qubits, x_stabs, z_stabs): + stabilizer_gens = [] + for stab in x_stabs: + gen = [0 for _ in range(2*num_qubits)] + for i in stab: + gen[i] = 1 + stabilizer_gens.append(gen) + for stab in z_stabs: + gen = [0 for _ in range(2*num_qubits)] + for i in stab: + gen[num_qubits + i] = 1 + stabilizer_gens.append(gen) + + num_gens = len(stabilizer_gens) + A = np.array(stabilizer_gens, dtype='int') + GF = galois.GF(2) + A = GF(A) + NS = A.null_space() + prod = np.array(matrix_product( + np.array(A, dtype=int), np.array(NS.transpose(), dtype=int), + addfunc=lambda x,y: (x+y)%2), dtype=int) + assert (prod == 0).all + dual_code_gens = np.array([ + xzdual(list(row)) for row in np.array(NS, dtype=int)], dtype='int') + x = {} + operator = [False for i in range(num_qubits*2)] + assert len(operator) == 2*num_qubits + for i, gen in enumerate(dual_code_gens): + assert len(operator) == 2*num_qubits + x[i] = z3.Bool(f'x_{i}') + for j in range(num_qubits*2): + operator[j] = FastXor(operator[j], FastAnd(x[i], bool(gen[j]))) + is_not_stabilizer = False + for i, gen in enumerate(dual_code_gens): + anticommutes = FastXorSum([ + FastAnd(bool(u), v) + for u,v in zip(xzdual(gen), operator)]) + is_not_stabilizer = FastOr(is_not_stabilizer, anticommutes) + total_weight = 0 + for j in range(num_qubits): + total_weight += z3.If(FastOr(operator[j], operator[num_qubits+j]), 1, 0) + constraints = [] + constraints.append(is_not_stabilizer) + opt = z3.Optimize() + handle = opt.minimize(total_weight) + for constraint in constraints: + opt.add(constraint) + assert opt.check() == z3.sat + assert opt.lower(handle) == opt.upper(handle) + return opt.upper(handle) diff --git a/requirements.txt b/requirements.txt index 04f0052..b069570 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1,2 +1,9 @@ +absl-py==1.1.0 +galois==0.0.30 +llvmlite==0.38.1 +numba==0.55.2 numpy==1.22.4 ortools==9.3.10497 +protobuf==4.21.2 +typing_extensions==4.3.0 +z3-solver==4.9.1.0 diff --git a/solver.py b/solver.py index 0a476df..03d40d6 100644 --- a/solver.py +++ b/solver.py @@ -3,9 +3,10 @@ from numpy.random import default_rng from ortools.sat.python import cp_model +from mindist import minimum_distance def random_graph(num_qubits, num_stabs, density, rng=default_rng()): - """ + """ Generate a random bipartite graph as a list of stabilizer adjacency. Each edge is added with probability `density`. """ @@ -70,17 +71,17 @@ def init_x_and_z_edges(self): ) def add_x_and_z_edge_constraints(self, activator, pauli, x, z): - """ - Decompose the constraints + """ + Decompose the constraints - x == activator AND pauli, + x == activator AND pauli, z == activator AND NOT(pauli), into the constraints - - Not(x) OR activator, - Not(z) OR activator, - Not(activator) OR x OR z, + + Not(x) OR activator, + Not(z) OR activator, + Not(activator) OR x OR z, Not(x) OR pauli, Not(z) OR NOT(pauli), """ @@ -115,11 +116,11 @@ def new_extra_var(self): def add_and_gate(self, left, right, target): """ Decompose the constraint - - left AND right == target - + + left AND right == target + into the contraints - + Not(target) OR left, Not(target) OR right, Not(left) OR Not(right) OR target, @@ -156,7 +157,7 @@ def add_commutation_constraint(self, stab0, stab1): def both_edges_are_active(self, qubit, stab0, stab1): """ - Add a variable representing that both edges (qubit, stab0) + Add a variable representing that both edges (qubit, stab0) and (qubit, stab1) are actives. """ acti0 = self.activator_var(qubit, stab0) @@ -238,7 +239,7 @@ def add_max_stab_deg_constraint(self, stab, qubits, max_deg): def with_balanced_stab_constraint(self): """ - Add the constraint that the numbers of X and Z stabilizers + Add the constraint that the numbers of X and Z stabilizers are the same (or differ by at most 1 if the number of stabilizers is odd). """ vars = [self.pauli_var(stab) for stab in range(self.num_stabs)] @@ -291,6 +292,7 @@ def build_stabilizers(generator, solver): print(f"It took {solver.WallTime()} seconds") if solver.StatusName(status) == "OPTIMAL": x_stabs, z_stabs = build_stabilizers(generator, solver) + print(f"Minimum distance is {minimum_distance(50, x_stabs, z_stabs)}") print() print(f"X stabilizers:") for s in x_stabs: