diff --git a/abp/graphstate.py b/abp/graphstate.py index 6e57df2..3ff9747 100644 --- a/abp/graphstate.py +++ b/abp/graphstate.py @@ -118,20 +118,27 @@ class GraphState(object): def measure(self, node, basis, force=None): """ Measure in an arbitrary basis """ basis = clifford.by_name[basis] - old_basis = basis ha = clifford.conjugation_table[self.node[node]["vop"]] basis, phase = clifford.conjugate(basis, ha) - assert phase in (-1, 1) # TODO: remove - # TODO: wtf - force = force ^ 0x01 if force != -1 and phase == 0 else force + #TODO: wtf + #force = force ^ 0x01 if force != -1 and phase == 0 else force + if force != None and phase == 0+0j: + force = not force + + result = random.choice([0, 1]) + if which == clifford.by_name["px"]: + result = self.measure_x(node, result) + elif which == clifford.by_name["py"]: + result = self.measure_y(node, result) + elif which == clifford.by_name["pz"]: + result = self.measure_z(node, result) + else: + raise ValueError("You can only measure in PX, PY or PZ") - which = {1: self.measure_x, 2: - self.measure_y, 3: self.measure_z}[basis] - res = which(node, force) - res = res if phase == 1 else not res + if phase == 1: + result = not result - # TODO: put the asserts from graphsim.cpp into tests return res def toggle_edges(a, b): @@ -142,14 +149,11 @@ class GraphState(object): done.add((i, j), (j, i)) self.toggle_edge(i, j) - def measure_x(self, node, force=None): + def measure_x(self, node, result): """ Measure the graph in the X-basis """ if len(self.adj[node]) == 0: return 0 - # Flip a coin - result = force if force != None else random.choice([0, 1]) - # Pick a vertex friend = next(self.adj[node].iterkeys()) @@ -179,11 +183,8 @@ class GraphState(object): return result - def measure_y(self, node, force=None): + def measure_y(self, node, result): """ Measure the graph in the Y-basis """ - # Flip a coin - result = force if force != None else random.choice([0, 1]) - # Do some rotations for neighbour in self.adj[node]: # NB: should these be hermitian_conjugated? @@ -197,11 +198,8 @@ class GraphState(object): self.act_local_rotation(node, "msqz" if result else "msqz_h") return result - def measure_z(self, node, force=None): + def measure_z(self, node, result): """ Measure the graph in the Z-basis """ - # Flip a coin - result = force if force != None else random.choice([0, 1]) - # Disconnect for neighbour in self.adj[node]: self.del_edge(node, neighbour) @@ -268,8 +266,10 @@ class GraphState(object): output[a][b] = "Z" else: output[a][b] = "I" + # TODO: signs return output def __eq__(self, other): """ Check equality between graphs """ return self.adj == other.adj and self.node == other.node + diff --git a/abp/qi.py b/abp/qi.py index e86f544..9c0fa13 100644 --- a/abp/qi.py +++ b/abp/qi.py @@ -9,12 +9,13 @@ And a circuit-model simulator import numpy as np import itertools as it + def hermitian_conjugate(u): """ Shortcut to the Hermitian conjugate """ return np.conjugate(np.transpose(u)) -# Constants -ir2 = 1/np.sqrt(2) +# Constants +ir2 = 1 / np.sqrt(2) # Operators id = np.array(np.eye(2, dtype=complex)) px = np.array([[0, 1], [1, 0]], dtype=complex) @@ -22,23 +23,29 @@ py = np.array([[0, -1j], [1j, 0]], dtype=complex) pz = np.array([[1, 0], [0, -1]], dtype=complex) ha = hadamard = np.array([[1, 1], [1, -1]], dtype=complex) * ir2 ph = np.array([[1, 0], [0, 1j]], dtype=complex) -t = np.array([[1, 0], [0, np.exp(1j*np.pi/4)]], dtype=complex) - -sqx = np.array([[ 1.+0.j, -0.+1.j], [-0.+1.j, 1.-0.j]], dtype=complex)*ir2 -msqx = np.array([[ 1.+0.j, 0.-1.j], [ 0.-1.j, 1.-0.j]], dtype=complex)*ir2 -sqy = np.array([[ 1.+0.j, 1.+0.j], [-1.-0.j, 1.-0.j]], dtype=complex)*ir2 -msqy = np.array([[ 1.+0.j, -1.-0.j], [ 1.+0.j, 1.-0.j]], dtype=complex)*ir2 -sqz = np.array([[ 1.+1.j, 0.+0.j], [ 0.+0.j, 1.-1.j]], dtype=complex)*ir2 -msqz = np.array([[ 1.-1.j, 0.+0.j], [ 0.+0.j, 1.+1.j]], dtype=complex)*ir2 +t = np.array([[1, 0], [0, np.exp(1j * np.pi / 4)]], dtype=complex) + +sqx = np.array( + [[1. + 0.j, -0. + 1.j], [-0. + 1.j, 1. - 0.j]], dtype=complex) * ir2 +msqx = np.array( + [[1. + 0.j, 0. - 1.j], [0. - 1.j, 1. - 0.j]], dtype=complex) * ir2 +sqy = np.array( + [[1. + 0.j, 1. + 0.j], [-1. - 0.j, 1. - 0.j]], dtype=complex) * ir2 +msqy = np.array( + [[1. + 0.j, -1. - 0.j], [1. + 0.j, 1. - 0.j]], dtype=complex) * ir2 +sqz = np.array( + [[1. + 1.j, 0. + 0.j], [0. + 0.j, 1. - 1.j]], dtype=complex) * ir2 +msqz = np.array( + [[1. - 1.j, 0. + 0.j], [0. + 0.j, 1. + 1.j]], dtype=complex) * ir2 # CZ gate cz = np.array(np.eye(4), dtype=complex) -cz[3,3]=-1 +cz[3, 3] = -1 # States -zero = np.array([[1],[0]], dtype=complex) -one = np.array([[0],[1]], dtype=complex) -plus = np.array([[1],[1]], dtype=complex) / np.sqrt(2) +zero = np.array([[1], [0]], dtype=complex) +one = np.array([[0], [1]], dtype=complex) +plus = np.array([[1], [1]], dtype=complex) / np.sqrt(2) bond = cz.dot(np.kron(plus, plus)) nobond = np.kron(plus, plus) @@ -60,11 +67,12 @@ def normalize_global_phase(m): class CircuitModel(object): + def __init__(self, nqubits): self.nqubits = nqubits - self.d = 2**nqubits + self.d = 2 ** nqubits self.state = np.zeros((self.d, 1), dtype=complex) - self.state[0, 0]=1 + self.state[0, 0] = 1 def act_cz(self, control, target): """ Act a CU somewhere """ @@ -86,19 +94,18 @@ class CircuitModel(object): output = np.zeros((self.d, 1), dtype=complex) for i, v in enumerate(self.state): q = i & where > 0 - output[i] += v*ha[q, q] - output[i ^ where] += v*ha[not q, q] + output[i] += v * ha[q, q] + output[i ^ where] += v * ha[not q, q] self.state = output - def act_local_rotation(self, qubit, u): """ Act a local unitary somwhere """ where = 1 << qubit output = np.zeros((self.d, 1), dtype=complex) for i, v in enumerate(self.state): q = i & where > 0 - output[i] += v*u[q, q] # TODO this is probably wrong - output[i ^ where] += v*u[not q, q] + output[i] += v * u[q, q] # TODO this is probably wrong + output[i ^ where] += v * u[not q, q] self.state = output def __eq__(self, other): @@ -108,12 +115,10 @@ class CircuitModel(object): b = normalize_global_phase(other.state) return np.allclose(a, b) - def __str__(self): s = "" for i in range(self.d): label = bin(i)[2:].rjust(self.nqubits, "0") - if abs(self.state[i, 0])>0.00001: + if abs(self.state[i, 0]) > 0.00001: s += "|{}>: {}\n".format(label, self.state[i, 0].round(3)) return s - diff --git a/tests/test_measurement.py b/tests/test_measurement.py index 686eaf9..28a7b0d 100644 --- a/tests/test_measurement.py +++ b/tests/test_measurement.py @@ -1,7 +1,10 @@ +import numpy as np from abp import GraphState +from abp import qi from anders_briegel import graphsim -def test_measurements(): +def test_single_qubit_measurements(): + """ Various simple tests of measurements """ # Test that measuring |0> in Z gives 0 g = GraphState([0]) @@ -20,7 +23,9 @@ def test_measurements(): g.act_local_rotation(0, "pz") assert g.measure(0, "px") == 1, "Measuring |-> in X gives 1" - # Test random outcomes + +def test_random_outcomes(): + """ Testing random behaviour """ ones = 0 for i in range(1000): g = GraphState([0]) @@ -28,6 +33,22 @@ def test_measurements(): ones += g.measure(0, "pz") assert 400 < ones < 600, "This is a probabilistic test!" +def test_projection(): + """ Test that projection works correctly """ + g = GraphState([0]) + g.act_local_rotation(0, "hadamard") + g.measure(0, "pz", 0) + print g.to_state_vector() + assert np.allclose(g.to_state_vector().state, qi.zero) + + g = GraphState([0]) + g.act_local_rotation(0, "hadamard") + print g.to_state_vector() + g.measure(0, "pz", 1) + print g.to_state_vector() + assert np.allclose(g.to_state_vector().state, qi.one) + + def test_z_measurement_against_ab():