feat: Random qubit initialisation (#32)

* handle corrections on uncorrected qubits

* add random register manager tests

* Add random register manager

* Add documentation

* Repair long running tests

* Simplify docs source

* Combine corrections where possible

* Rename initialisation register list

* Add aditional documentation

pytket-mbqc-py/pytket_mbqc_py/__init__.py

@@ -5,3 +5,4 @@
 from .cnot_block import CNOTBlocksGraphCircuit
 from .graph_circuit import GraphCircuit
 from .qubit_manager import QubitManager
+from .random_register_manager import RandomRegisterManager

pytket-mbqc-py/pytket_mbqc_py/cnot_block.py

@@ -23,21 +23,25 @@ def __init__(
         n_physical_qubits: int,
         input_state: Tuple[int],
         n_layers: int,
+        n_registers: int,
     ) -> None:
         """Initialisation method.
 
         :param n_physical_qubits: The maximum number of physical qubits
             available. These qubits will be reused and so the total
             number of 'logical' qubits may be larger.
-        :type n_physical_qubits: int
         :param input_state: Integer tuple describing the input
             to the circuit. This will be a classical binary
             string and so the outcome is deterministic.
-        :type input_state: Tuple[int]
         :param n_layers: The number of layers of CNOT gates.
-        :type n_layers: int
+        :param n_registers: The number of classical registers to use
+            for state preparation information. Note that there should
+            be at least one register per logical qubit.
         """
 
+        # TODO: n_registers should be calculated from the given
+        # input variables.
+
         self.input_state = input_state
         self.n_layers = n_layers
 
@@ -56,7 +60,10 @@ def __init__(
             [[] for _ in range(n_rows)] for _ in range(n_layers)
         ]
 
-        super().__init__(n_physical_qubits=n_physical_qubits)
+        super().__init__(
+            n_physical_qubits=n_physical_qubits,
+            n_registers=n_registers,
+        )
 
         for layer in range(n_layers):
             # If this is the first layer then the control qubit of the first row needs

pytket-mbqc-py/pytket_mbqc_py/graph_circuit.py

@@ -5,29 +5,32 @@
 """
 
 from functools import reduce
-from typing import Dict, List, Tuple
+from typing import Dict, List, Tuple, Union
 
 import networkx as nx  # type:ignore
 from pytket import Qubit
 from pytket.circuit.logic_exp import BitLogicExp
 from pytket.unit_id import BitRegister
 
-from pytket_mbqc_py.qubit_manager import QubitManager
+from pytket_mbqc_py.random_register_manager import RandomRegisterManager
 
 
-class GraphCircuit(QubitManager):
+class GraphCircuit(RandomRegisterManager):
     """Class for the automated construction of MBQC computations.
     In particular only graphs with valid flow can be constructed.
     Graph state construction and measurement corrections are added
-    automatically. A child of
-    :py:class:`~pytket_mbqc_py.qubit_manager.QubitManager`.
+    automatically.
 
     :ivar entanglement_graph: Graph detailing the graph state entanglement.
     :ivar flow_graph: Graph describing the flow dependencies of the graph state.
     :ivar vertex_qubit: List mapping graph vertex to corresponding qubits.
     :ivar vertex_measured: List indicating if vertex has been measured.
     :ivar vertex_x_corr_reg: List mapping vertex index to
         the classical register where the required X correction is stored.
+    :ivar vertex_init_reg: List mapping vertex to the register describing
+        the state it was initialised in. In particular this is a 3 bit register
+        with the 0th entry giving the T rotation, the 1st giving the
+        S rotation, and the 2nd giving the Z rotation.
     """
 
     entanglement_graph: nx.Graph
@@ -37,12 +40,16 @@ class GraphCircuit(QubitManager):
     def __init__(
         self,
         n_physical_qubits: int,
+        n_registers: int,
     ) -> None:
         """Initialisation method. Creates tools to track
         the graph state structure and the measurement corrections.
 
         :param n_physical_qubits: The number of physical qubits available.
-        :type n_physical_qubits: int
+        :param n_registers: The number of state initialisation registers
+            to generate. Each register describes the state that the logical
+            qubit is initialised in. Note that the number of such registers
+            should be at least the number of logical qubits.
         """
         super().__init__(n_physical_qubits=n_physical_qubits)
 
@@ -60,6 +67,16 @@ def __init__(
         # vertex has been measured.
         self.vertex_x_corr_reg: List[BitRegister] = []
 
+        # Generate one random register per vertex.
+        # When qubits are added they will be initialised in this
+        # random register. This in except for the case of input qubits
+        # which are initialised in the 0 state, and in which case this
+        # register is overwritten.
+        self.vertex_init_reg = list(
+            self.generate_random_registers(n_registers=n_registers)
+        )
+        self.add_barrier(units=self.qubits)
+
     def get_outputs(self) -> Dict[int, Qubit]:
         """Return the output qubits. Output qubits are those that
         are unmeasured, and which do not have a flow. This should
@@ -101,7 +118,29 @@ def get_outputs(self) -> Dict[int, Qubit]:
 
         # We need to correct all unmeasured qubits.
         for vertex in output_qubits.keys():
-            self._apply_x_correction(vertex=vertex)
+            # Corrections are first applied to invert the initialisation
+            self.T(self.vertex_qubit[vertex], condition=self.vertex_init_reg[vertex][0])
+            self.S(self.vertex_qubit[vertex], condition=self.vertex_init_reg[vertex][0])
+
+            self.S(self.vertex_qubit[vertex], condition=self.vertex_init_reg[vertex][1])
+
+            self.Z(
+                self.vertex_qubit[vertex],
+                condition=(
+                    self.vertex_init_reg[vertex][2]
+                    ^ self.vertex_init_reg[vertex][1]
+                    ^ self.vertex_init_reg[vertex][0]
+                ),
+            )
+
+            # Apply X correction according to correction register.
+            # These are corrections resulting from the measurements
+            self.X(
+                self.vertex_qubit[vertex],
+                condition=self.vertex_x_corr_reg[vertex][0],
+            )
+
+            # Apply Z corrections resulting from measurements.
             self._apply_z_correction(vertex=vertex)
 
         return output_qubits
@@ -114,6 +153,9 @@ def _add_vertex(self, qubit: Qubit) -> int:
 
         :param qubit: Qubit to be added.
         :return: The vertex in the graphs corresponding to this qubit
+
+        :raises Exception: Raised if an insufficient number of initialisation
+            registers were initialised.
         """
         self.vertex_qubit.append(qubit)
         self.vertex_measured.append(False)
@@ -126,6 +168,12 @@ def _add_vertex(self, qubit: Qubit) -> int:
         self.vertex_x_corr_reg.append(x_corr_reg)
         self.add_c_register(register=x_corr_reg)
 
+        if index >= len(self.vertex_init_reg):
+            raise Exception(
+                "An insufficient number of initialisation registers "
+                + "were initialised."
+            )
+
         return index
 
     def add_input_vertex(self) -> Tuple[Qubit, int]:
@@ -140,6 +188,10 @@ def add_input_vertex(self) -> Tuple[Qubit, int]:
         qubit = self.get_qubit()
         index = self._add_vertex(qubit=qubit)
 
+        # In the case of input qubits, the initialisation is not random.
+        # As such the initialisation register should be set to 0.
+        self.add_c_setreg(0, self.vertex_init_reg[index])
+
         return (qubit, index)
 
     def add_graph_vertex(self) -> int:
@@ -151,6 +203,12 @@ def add_graph_vertex(self) -> int:
         self.H(qubit)
         index = self._add_vertex(qubit=qubit)
 
+        # The graph state is randomly initialised based on the
+        # initialisation register.
+        self.T(qubit, condition=self.vertex_init_reg[index][0])
+        self.S(qubit, condition=self.vertex_init_reg[index][1])
+        self.Z(qubit, condition=self.vertex_init_reg[index][2])
+
         return index
 
     @property
@@ -240,8 +298,8 @@ def add_edge(self, vertex_one: int, vertex_two: int) -> None:
         ):
             raise Exception(
                 "This does not define a valid flow. "
-                f"In particular {vertex_two} is the flow of {self.flow_graph.predecessors(vertex_two)}, "
-                f"some of which are measured before {vertex_one}."
+                f"In particular {vertex_two} is the flow of {list(self.flow_graph.predecessors(vertex_two))}, "
+                f"some of which are measured after {vertex_one}."
             )
 
         # If this is the first future of vertex_one then it is taken to be its flow.
@@ -260,34 +318,26 @@ def add_edge(self, vertex_one: int, vertex_two: int) -> None:
             v_of_edge=vertex_two,
         )
 
-    def _apply_x_correction(self, vertex: int) -> None:
-        """Apply X correction. This correction is drawn from
-        the x correction register which is altered
-        by the corrected measure method as appropriate.
-
-        :param vertex: The vertex to be corrected.
-        """
-        self.X(
-            self.vertex_qubit[vertex],
-            condition=self.vertex_x_corr_reg[vertex][0],
-        )
-
-    def _get_z_correction_expression(self, vertex: int) -> BitLogicExp:
+    def _get_z_correction_expression(self, vertex: int) -> Union[None, BitLogicExp]:
         """Create logical expression by taking the parity of
         the X corrections that have to be applied to the neighbouring
-        qubits.
+        qubits. If there are no neighbours then None will be returned.
 
         :param vertex: Vertex to be corrected.
         :return: Logical expression calculating the parity
             of the neighbouring x correction registers.
         """
-        return reduce(
-            lambda a, b: a ^ b,
-            [
-                self.vertex_x_corr_reg[neighbour][0]
-                for neighbour in self.entanglement_graph.neighbors(n=vertex)
-            ],
-        )
+
+        neighbour_reg_list = [
+            self.vertex_x_corr_reg[neighbour][0]
+            for neighbour in self.entanglement_graph.neighbors(n=vertex)
+        ]
+
+        # This happens of this vertex has no neighbours.
+        if len(neighbour_reg_list) == 0:
+            return None
+
+        return reduce(lambda a, b: a ^ b, neighbour_reg_list)
 
     def _apply_z_correction(self, vertex: int) -> None:
         """Apply Z correction on qubit. This correction is calculated
@@ -296,10 +346,12 @@ def _apply_z_correction(self, vertex: int) -> None:
 
         :param vertex: Vertex to be corrected.
         """
-        self.Z(
-            self.vertex_qubit[vertex],
-            condition=self._get_z_correction_expression(vertex=vertex),
-        )
+        condition = self._get_z_correction_expression(vertex=vertex)
+        if condition is not None:
+            self.Z(
+                self.vertex_qubit[vertex],
+                condition=condition,
+            )
 
     def _apply_classical_z_correction(self, vertex: int) -> None:
         """Apply Z correction on measurement result. This correction is calculated
@@ -308,11 +360,13 @@ def _apply_classical_z_correction(self, vertex: int) -> None:
 
         :param vertex: Vertex to be corrected.
         """
-        self.add_classicalexpbox_bit(
-            expression=self.qubit_meas_reg[self.vertex_qubit[vertex]][0]
-            ^ self._get_z_correction_expression(vertex=vertex),
-            target=[self.qubit_meas_reg[self.vertex_qubit[vertex]][0]],
-        )
+        condition = self._get_z_correction_expression(vertex=vertex)
+        if condition is not None:
+            self.add_classicalexpbox_bit(
+                expression=self.qubit_meas_reg[self.vertex_qubit[vertex]][0]
+                ^ condition,
+                target=[self.qubit_meas_reg[self.vertex_qubit[vertex]][0]],
+            )
 
     def corrected_measure(self, vertex: int, t_multiple: int = 0) -> None:
         """Perform a measurement, applying the appropriate corrections.
@@ -321,7 +375,6 @@ def corrected_measure(self, vertex: int, t_multiple: int = 0) -> None:
         :param vertex: Vertex to be measured.
         :param t_multiple: The angle in which to measure, defaults to 0.
             This defines the rotated hadamard basis to measure in.
-        :type t_multiple: int, optional
         :raises Exception: Raised if this vertex has already been measured.
         :raises Exception: Raised if there are vertex in the past of this
             one which have not been measured.
@@ -342,11 +395,68 @@ def corrected_measure(self, vertex: int, t_multiple: int = 0) -> None:
                 + f"are in the future of {vertex} and have already been measured."
             )
 
-        # This is actually optional for graph vertices
-        # as the X correction commutes with the hadamard
-        # basis measurement.
-        self._apply_x_correction(vertex=vertex)
+        # Apply X correction according to correction register.
+        # This is to correct for measurement outcomes.
+        self.X(
+            self.vertex_qubit[vertex],
+            condition=self.vertex_x_corr_reg[vertex][0],
+        )
+
+        self.T(
+            self.vertex_qubit[vertex],
+            # Required to invert random T from initialisation.
+            condition=self.vertex_init_reg[vertex][0],
+        )
+        self.S(
+            self.vertex_qubit[vertex],
+            # Required to invert random T from initialisation.
+            # This additional term is required to account for the case where
+            # the correcting T is commuted through an X correction.
+            condition=(
+                self.vertex_init_reg[vertex][0] & self.vertex_x_corr_reg[vertex][0]
+            ),
+        )
+
+        self.S(
+            self.vertex_qubit[vertex],
+            # Required to invert random T from initialisation.
+            condition=self.vertex_init_reg[vertex][0],
+        )
+
+        self.S(
+            self.vertex_qubit[vertex],
+            # Required to invert random S from initialisation.
+            condition=self.vertex_init_reg[vertex][1],
+        )
+
+        self.Z(
+            self.vertex_qubit[vertex],
+            condition=(
+                # Required to invert random T from initialisation.
+                self.vertex_init_reg[vertex][0]
+                # Required to invert random S from initialisation.
+                ^ self.vertex_init_reg[vertex][1]
+                # Required to invert random Z from initialisation.
+                ^ self.vertex_init_reg[vertex][2]
+                # Required to invert random S from initialisation.
+                # This additional term is required to account for the case where
+                # the correcting S is commuted through an X correction.
+                ^ (self.vertex_init_reg[vertex][1] & self.vertex_x_corr_reg[vertex][0])
+                # Required to invert random T from initialisation.
+                # This additional term is required to account for the case where
+                # the correcting S is commuted through an X correction.
+                ^ (self.vertex_init_reg[vertex][0] & self.vertex_x_corr_reg[vertex][0])
+                # Required to invert random T from initialisation.
+                # This additional term is required to account for the case where
+                # the correcting T is commuted through an X correction.
+                ^ (self.vertex_init_reg[vertex][0] & self.vertex_x_corr_reg[vertex][0])
+            ),
+        )
 
+        # Rotate measurement basis.
+        # TODO: These measurements should be combined with the above
+        # so that the measurement angles are hidden by the initialisation
+        # angles.
         inverse_t_multiple = 8 - t_multiple
         inverse_t_multiple = inverse_t_multiple % 8
         if inverse_t_multiple // 4: