Rework error correction examples

README.md

@@ -5,7 +5,7 @@ This repository contains my various exercises and research projects in the field
 List of projects:
 * Exercises with Cirq:
   * [Phase estimation](Phase%20estimation.ipynb) (2019)
-  * [Error correction](Quantum%20error%20correction%20with%20Cirq.ipynb) (2019)
+  * [Error correction](error_correction/) (2019, upd. in 2026)
 * Linear algebra:
   * [Schmidt decomposition of a vector](Schmidt%20decomposition%20of%20a%20vector.ipynb) (2019)
   * [Schmidt decomposition of a 4x4 matrix](Schmidt%20decomposition%20of%204x4%20matrix.ipynb) (2019)

cirq_sparse_sim/sparse_sim.py

@@ -694,6 +694,35 @@ class SparseSimulatorTrialResult(
 ):
     """Final sparse simulation state, materialized as a vector only on request."""
 
+    def bloch_vector_of(self, qubit: cirq.Qid) -> np.ndarray:
+        """Return the Bloch vector of a qubit in the final state."""
+        state = self._get_merged_sim_state()
+        qubit_index = self.qubit_map[qubit]
+        sparse_state = state.sparse_state
+        axis = sparse_state.axis_by_qubit[state.qubits[qubit_index]]
+        bit = 1 << axis
+        amplitudes = dict(
+            zip(
+                sparse_state.basis_states,
+                sparse_state.amplitudes,
+                strict=True,
+            )
+        )
+
+        z = 0.0
+        coherence = 0.0j
+        for basis, amplitude in amplitudes.items():
+            probability = abs(amplitude) ** 2
+            if basis & bit:
+                z -= probability
+            else:
+                z += probability
+                coherence += amplitude * np.conj(amplitudes.get(basis | bit, 0.0j))
+
+        return np.array(
+            [2 * coherence.real, -2 * coherence.imag, z], dtype=np.float32
+        )
+
     @property
     def final_state_vector(self) -> np.ndarray:
         return self._get_merged_sim_state().state_vector()

cirq_sparse_sim/sparse_sim_test.py

@@ -1463,6 +1463,49 @@ def test_simulate_random_circuits_matches_dense_simulator(seed: int) -> None:
     )
 
 
+def test_simulation_result_bloch_vectors_match_dense_simulator() -> None:
+    simulator = SparseSimulator()
+    qubits = simulator.qubit_manager.qalloc(3)
+    circuit = cirq.Circuit(
+        cirq.ry(0.7)(qubits[0]),
+        cirq.rz(-1.2)(qubits[0]),
+        cirq.H(qubits[1]),
+        cirq.CNOT(qubits[1], qubits[2]),
+    )
+    order = [qubits[2], qubits[0], qubits[1]]
+
+    actual = simulator.simulate(circuit, qubit_order=order)
+    expected = cirq.Simulator(dtype=np.complex128).simulate(
+        circuit, qubit_order=order
+    )
+
+    for qubit in qubits:
+        assert actual.bloch_vector_of(qubit).dtype == np.float32
+        np.testing.assert_allclose(
+            actual.bloch_vector_of(qubit),
+            expected.bloch_vector_of(qubit),
+            rtol=0,
+            atol=1e-7,
+        )
+    with pytest.raises(KeyError):
+        actual.bloch_vector_of(cirq.NamedQubit("missing"))
+
+
+def test_simulation_result_bloch_vector_stays_sparse() -> None:
+    simulator = SparseSimulator()
+    qubits = simulator.qubit_manager.qalloc(128)
+    result = simulator.simulate(
+        cirq.Circuit(cirq.X(qubits[-1])), qubit_order=qubits
+    )
+
+    np.testing.assert_array_equal(
+        result.bloch_vector_of(qubits[0]), np.array([0, 0, 1], dtype=np.float32)
+    )
+    np.testing.assert_array_equal(
+        result.bloch_vector_of(qubits[-1]), np.array([0, 0, -1], dtype=np.float32)
+    )
+
+
 def test_simulate_measurements_and_final_state_are_independent_snapshots() -> None:
     simulator = SparseSimulator(seed=1)
     qubits = simulator.qubit_manager.qalloc(2)

error_correction/README.md

@@ -0,0 +1,6 @@
+### Error correction
+
+This folder contains some very basic exercises with error correction codes.
+
+Notebooks must be run with path set to repository root, e.g.
+`PYTHONPATH="$PWD" jupyter notebook`

error_correction/channels.py

@@ -0,0 +1,34 @@
+import numpy as np
+import cirq
+
+
+class BitFlipChannel:
+    def __init__(self, flip_prob):
+        self.flip_prob = flip_prob
+
+    def transmit(self, ct, q):
+        if np.random.rand() < self.flip_prob:
+            ct.append(cirq.X(q))
+        return q
+
+
+class PhaseFlipChannel:
+    def __init__(self, flip_prob):
+        self.flip_prob = flip_prob
+
+    def transmit(self, ct, q):
+        if np.random.rand() < self.flip_prob:
+            ct.append(cirq.Z(q))
+        return q
+
+
+class ArbitraryErrorChannel:
+    def __init__(self, flip_prob):
+        self.flip_prob = flip_prob
+
+    def transmit(self, ct, q):
+        if np.random.rand() < self.flip_prob:
+            ct.append(cirq.Rx(rads=np.random.rand() * 2 * np.pi).on(q))
+            ct.append(cirq.Ry(rads=np.random.rand() * 2 * np.pi).on(q))
+            ct.append(cirq.Rz(rads=np.random.rand() * 2 * np.pi).on(q))
+        return q

error_correction/protocols.py

@@ -0,0 +1,141 @@
+import cirq
+
+from cirq import Circuit, Qid
+
+from abc import ABC, abstractmethod
+
+
+class QecProtocol(ABC):
+    @abstractmethod
+    def encode(self, circuit: Circuit, qubit: list[Qid]) -> list[Qid]: ...
+
+    @abstractmethod
+    def decode(self, circuit: Circuit, qubits: list[Qid]) -> list[Qid]: ...
+
+
+class NoEncodingProtocol(QecProtocol):
+    def __init__(self):
+        self.name = "No encoding"
+
+    def encode(self, circuit: Circuit, qubits: list[Qid]) -> list[Qid]:
+        return qubits
+
+    def decode(self, circuit: Circuit, qubits: list[Qid]) -> list[Qid]:
+        return qubits
+
+
+class ThreeQubitBitFlipProtocol(QecProtocol):
+    def __init__(self):
+        self.name = "3 qubit bit flip protocol"
+
+    def encode(self, circuit, qubits: list[Qid]) -> list[Qid]:
+        assert len(qubits) == 1
+        q0 = cirq.NamedQubit("aux_%d" % len(circuit.all_qubits()))
+        circuit.append(cirq.CNOT(qubits[0], q0))
+        q1 = cirq.NamedQubit("aux_%d" % len(circuit.all_qubits()))
+        circuit.append(cirq.CNOT(qubits[0], q1))
+        return [qubits[0], q0, q1]
+
+    def decode(self, circuit, qubits: list[Qid]) -> list[Qid]:
+        def X(target):
+            circuit.append(cirq.X(qubits[target]))
+
+        def CCNOT(target):
+            i1, i2 = 1, 2
+            if target == 1:
+                i1, i2 = 0, 2
+            if target == 2:
+                i1, i2 = 0, 1
+            circuit.append(cirq.CCNOT(qubits[i1], qubits[i2], qubits[target]))
+
+        X(1)
+        CCNOT(0)
+        X(1)
+        CCNOT(1)
+        CCNOT(0)
+        X(0)
+        CCNOT(2)
+        X(0)
+        X(2)
+        CCNOT(0)
+        X(0)
+        X(2)
+        CCNOT(2)
+        X(0)
+        X(1)
+        CCNOT(0)
+        X(1)
+        CCNOT(1)
+        CCNOT(0)
+
+        # Measurement is not needed for Bit-Flip, but is needed so we can use it in Shor Code.
+        circuit.append(cirq.measure(qubits[0]))
+        circuit.append(cirq.measure(qubits[1]))
+
+        return [qubits[2]]
+
+
+class ThreeQubitPhaseFlipProtocol(QecProtocol):
+    def __init__(self):
+        self.name = "3 qubit phase flip protocol"
+        self.bf_protocol = ThreeQubitBitFlipProtocol()
+
+    def encode(self, circuit, qubits: list[Qid]) -> list[Qid]:
+        qubits = self.bf_protocol.encode(circuit, qubits)
+        for q in qubits:
+            circuit.append(cirq.H(q))
+        return qubits
+
+    def decode(self, circuit, qubits: list[Qid]) -> list[Qid]:
+        for q in qubits:
+            circuit.append(cirq.H(q))
+        return self.bf_protocol.decode(circuit, qubits)
+
+
+class ShorProtocol(QecProtocol):
+    def __init__(self):
+        self.name = "Shor Code"
+        self.bf_protocol = ThreeQubitBitFlipProtocol()
+        self.pf_protocol = ThreeQubitPhaseFlipProtocol()
+
+    def encode(self, circuit, qubits: list[Qid]) -> list[Qid]:
+        assert len(qubits) == 1
+        result = []
+        qubits1 = self.pf_protocol.encode(circuit, qubits)
+        for q in qubits1:
+            result += self.bf_protocol.encode(circuit, [q])
+        return result
+
+    def decode(self, circuit, qubits: list[Qid]) -> list[Qid]:
+        return self.pf_protocol.decode(
+            circuit,
+            [
+                self.bf_protocol.decode(circuit, qubits[0:3])[0],
+                self.bf_protocol.decode(circuit, qubits[3:6])[0],
+                self.bf_protocol.decode(circuit, qubits[6:9])[0],
+            ],
+        )
+
+
+class NineQubitBitFlipProtocol(QecProtocol):
+    def __init__(self):
+        self.name = "9 qubit bit flip protocol"
+        self.bf_protocol_1 = ThreeQubitBitFlipProtocol()
+        self.bf_protocol_2 = ThreeQubitBitFlipProtocol()
+
+    def encode(self, circuit, qubit):
+        result = []
+        qubits1 = self.bf_protocol_1.encode(circuit, qubit)
+        for q in qubits1:
+            result += self.bf_protocol_2.encode(circuit, q)
+        return result
+
+    def decode(self, circuit, qubits):
+        return self.bf_protocol_1.decode(
+            circuit,
+            [
+                self.bf_protocol_2.decode(circuit, qubits[0:3]),
+                self.bf_protocol_2.decode(circuit, qubits[3:6]),
+                self.bf_protocol_2.decode(circuit, qubits[6:9]),
+            ],
+        )

error_correction/protocols_test.py

@@ -0,0 +1,9 @@
+from error_correction.channels import ArbitraryErrorChannel
+from error_correction.protocols import ShorProtocol
+from error_correction.utils import test_protocol
+
+
+def test_shor_code():
+    protocol = ShorProtocol()
+    channel = ArbitraryErrorChannel(0.01)
+    assert test_protocol(protocol, channel, num_experiments=10) <= 0.1