Removed inputs from Pauli presimulation

CHANGELOG.md

@@ -15,6 +15,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
   - Introduced new wrapper methods in the `Pattern` class: `graphix.pattern.Pattern.extract_partial_order_layers`, `graphix.pattern.Pattern.extract_causal_flow` and `graphix.pattern.Pattern.extract_gflow`.
   - Introduced new module `graphix.flow._partial_order` with the function :func:`compute_topological_generations`.
 
+- #392: Added `graphix.pattern.Pattern.remove_input_nodes` method which removes the input nodes from the pattern and replaces them with N commands.
+
 - #385
   - Introduced `graphix.flow.core.XZCorrections.check_well_formed` which verifies the correctness of an XZ-corrections instance and raises an exception if incorrect.
   - Added XZ-correction exceptions to module `graphix.flow.core.exceptions`.
@@ -35,6 +37,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 ### Fixed
 
+- #392: `Pattern.remove_input_nodes` is required before the `Pattern.perform_pauli_measurements` method to ensure input nodes are removed and fixed in the |+> state.
+
 - #379: Removed unnecessary `meas_index` from API for rotation instructions `RZ`, `RY` and `RX`.
 
 - #347: Adapted existing method `graphix.opengraph.OpenGraph.isclose` to the new API introduced in #358.

docs/source/tutorial.rst

@@ -175,12 +175,13 @@ It is known that quantum circuit consisting of Pauli basis states, Clifford gate
 <https://en.wikipedia.org/wiki/Gottesman%E2%80%93Knill_theorem>`_; e.g. the graph state simulator runs in :math:`\mathcal{O}(n \log n)` time).
 The Pauli measurement part of the MBQC is exactly this, and they can be preprocessed by our graph state simulator :class:`~graphix.graphsim.GraphState` - see :doc:`lc-mbqc` for more detailed description.
 
-We can call this in a line by calling :meth:`~graphix.pattern.Pattern.perform_pauli_measurements()` of :class:`~graphix.pattern.Pattern` object, which acts as the optimization routine of the measurement pattern.
+We can call this in a line by calling :meth:`~graphix.pattern.Pattern.remove_input_nodes` followed by :meth:`~graphix.pattern.Pattern.perform_pauli_measurements()` (both methods of the :class:`~graphix.pattern.Pattern` object). The first method removes the input nodes, while the second method optimizes the measurement pattern.
 We get an updated measurement pattern without Pauli measurements as follows:
 
+>>> pattern.remove_input_nodes()
 >>> pattern.perform_pauli_measurements()
 >>> pattern
-Pattern(input_nodes=[0, 1], cmds=[N(3), N(7), E((0, 3)), E((1, 3)), E((1, 7)), M(0, Plane.YZ, 0.2907266109187514), M(1, Plane.YZ, 0.01258854060311348), C(3, Clifford.I), C(7, Clifford.I), Z(3, {0, 1, 5}), Z(7, {1, 5}), X(3, {2}), X(7, {2, 4, 6})], output_nodes=[3, 7])
+Pattern(input_nodes=[], cmds=[N(0), N(1), N(3), N(7), E((0, 3)), E((1, 3)), E((1, 7)), M(0, Plane.YZ, 0.2907266109187514), M(1, Plane.YZ, 0.01258854060311348), C(3, Clifford.I), C(7, Clifford.I), Z(3, {0, 1, 5}), Z(7, {1, 5}), X(3, {2}), X(7, {2, 4, 6})], output_nodes=[3, 7])
 
 
 Notice that all measurements with angle=0 (Pauli X measurements) disappeared - this means that a part of quantum computation was *classically* (and efficiently) preprocessed such that we only need much smaller quantum resource.
@@ -210,7 +211,7 @@ We can simply call :meth:`~graphix.pattern.Pattern.minimize_space()` to reduce t
 
 >>> pattern.minimize_space()
 >>> pattern
-Pattern(input_nodes=[0, 1], cmds=[N(3), E((0, 3)), M(0, Plane.YZ, 0.11120090987081546), E((1, 3)), N(7), E((1, 7)), M(1, Plane.YZ, 0.230565199664617), C(3, Clifford.I), C(7, Clifford.I), Z(3, {0, 1, 5}), Z(7, {1, 5}), X(3, {2}), X(7, {2, 4, 6})], output_nodes=[3, 7])
+Pattern(input_nodes=[], cmds=[N(0), N(1), N(3), E((0, 3)), M(0, Plane.YZ, 0.11120090987081546), E((1, 3)), N(7), E((1, 7)), M(1, Plane.YZ, 0.230565199664617), C(3, Clifford.I), C(7, Clifford.I), Z(3, {0, 1, 5}), Z(7, {1, 5}), X(3, {2}), X(7, {2, 4, 6})], output_nodes=[3, 7])
 
 
 With the original measurement pattern, the simulation should have proceeded as follows, with maximum of four qubits on the memory.

examples/deutsch_jozsa.py

@@ -71,8 +71,9 @@
 print(pattern.to_ascii(left_to_right=True, limit=15))
 
 # %%
-# Now we preprocess all Pauli measurements
+# Now we preprocess all Pauli measurements, which requires that we move inputs to N commands
 
+pattern.remove_input_nodes()
 pattern.perform_pauli_measurements()
 print(
     pattern.to_ascii(

examples/mbqc_vqe.py

@@ -89,6 +89,7 @@ def build_mbqc_pattern(self, params: Iterable[Angle]) -> Pattern:
         pattern = circuit.transpile().pattern
         pattern.standardize()
         pattern.shift_signals()
+        pattern.remove_input_nodes()
         pattern.perform_pauli_measurements()  # Perform Pauli measurements
         return pattern
 

examples/qaoa.py

@@ -40,6 +40,7 @@
 # %%
 # perform Pauli measurements and plot the new (minimal) graph to perform the same quantum computation
 
+pattern.remove_input_nodes()
 pattern.perform_pauli_measurements()
 pattern.draw_graph(flow_from_pattern=False)
 

examples/qft_with_tn.py

@@ -62,6 +62,7 @@ def qft(circuit, n):
 # %%
 # Using efficient graph state simulator `graphix.graphsim`, we can classically preprocess Pauli measurements.
 # We are currently improving the speed of this process by using rust-based graph manipulation backend.
+pattern.remove_input_nodes()
 pattern.perform_pauli_measurements()
 
 

examples/tn_simulation.py

@@ -70,7 +70,7 @@ def ansatz(circuit, n, gamma, beta, iterations):
 
 # %%
 # Optimizing by performing Pauli measurements in the pattern using efficient stabilizer simulator.
-
+pattern.remove_input_nodes()
 pattern.perform_pauli_measurements()
 
 # %%
@@ -181,6 +181,7 @@ def cost(params, n, ham, quantum_iter, slice_index, opt=None):
     pattern = circuit.transpile().pattern
     pattern.standardize()
     pattern.shift_signals()
+    pattern.remove_input_nodes()
     pattern.perform_pauli_measurements()
     mbqc_tn = pattern.simulate_pattern(backend="tensornetwork", graph_prep="parallel")
     exp_val = 0

examples/visualization.py

@@ -37,8 +37,8 @@
 
 # %%
 # next, show the gflow:
-
-pattern.perform_pauli_measurements(leave_input=True)
+pattern.remove_input_nodes()
+pattern.perform_pauli_measurements()
 pattern.draw_graph(flow_from_pattern=False, show_measurement_planes=True, node_distance=(1, 0.6))
 
 

graphix/pattern.py

@@ -124,6 +124,7 @@ def add(self, cmd: Command) -> None:
         cmd : :class:`graphix.command.Command`
             MBQC command.
         """
+        self._pauli_preprocessed = False
         if cmd.kind == CommandKind.N:
             self.__n_node += 1
             self.__output_nodes.append(cmd.node)
@@ -137,6 +138,7 @@ def extend(self, *cmds: Command | Iterable[Command]) -> None:
 
         :param cmds: sequences of commands
         """
+        self._pauli_preprocessed = False
         for item in cmds:
             if isinstance(item, Iterable):
                 for cmd in item:
@@ -146,6 +148,7 @@ def extend(self, *cmds: Command | Iterable[Command]) -> None:
 
     def clear(self) -> None:
         """Clear the sequence of pattern commands."""
+        self._pauli_preprocessed = False
         self.__n_node = len(self.__input_nodes)
         self.__seq = []
         self.__output_nodes = list(self.__input_nodes)
@@ -157,6 +160,7 @@ def replace(self, cmds: list[Command], input_nodes: list[int] | None = None) ->
 
         :param input_nodes: optional, list of input qubits (by default, keep the same input nodes as before)
         """
+        self._pauli_preprocessed = False
         if input_nodes is not None:
             self.__input_nodes = list(input_nodes)
         self.clear()
@@ -202,6 +206,7 @@ def compose(
         - Input (and, respectively, output) nodes in the returned pattern have the order of the pattern `self` followed by those of the pattern `other`. Merged nodes are removed.
         - If `preserve_mapping = True` and :math:`|M_1| = |I_2| = |O_2|`, then the outputs of the returned pattern are the outputs of pattern `self`, where the nth merged output is replaced by the output of pattern `other` corresponding to its nth input instead.
         """
+        self._pauli_preprocessed = False
         nodes_p1 = self.extract_nodes() | self.results.keys()  # Results contain preprocessed Pauli nodes
         nodes_p2 = other.extract_nodes() | other.results.keys()
 
@@ -1402,14 +1407,21 @@ def simulate_pattern(
         sim.run(input_state, rng=rng)
         return sim.backend.state
 
-    def perform_pauli_measurements(self, leave_input: bool = False, ignore_pauli_with_deps: bool = False) -> None:
+    def remove_input_nodes(self) -> None:
+        """Remove the input nodes from the pattern and replace them with N commands.
+
+        This removes the possibility of choosing the input state, fixing the input state to the plus state.
+        .. seealso:: :class:`graphix.command.N`
+        """
+        self.__seq[0:0] = [command.N(node=node) for node in self.input_nodes]
+        empty_nodes: list[int] = []
+        self.__input_nodes = empty_nodes
+
+    def perform_pauli_measurements(self, ignore_pauli_with_deps: bool = False) -> None:
         """Perform Pauli measurements in the pattern using efficient stabilizer simulator.
 
         Parameters
         ----------
-        leave_input : bool
-            Optional (*False* by default).
-            If *True*, measurements on input nodes are preserved as-is in the pattern.
         ignore_pauli_with_deps : bool
             Optional (*False* by default).
             If *True*, Pauli measurements with domains depending on other measures are preserved as-is in the pattern.
@@ -1418,9 +1430,9 @@ def perform_pauli_measurements(self, leave_input: bool = False, ignore_pauli_wit
         .. seealso:: :func:`measure_pauli`
 
         """
-        # if not ignore_pauli_with_deps:
-        #    self.move_pauli_measurements_to_the_front()
-        measure_pauli(self, leave_input, copy=False, ignore_pauli_with_deps=ignore_pauli_with_deps)
+        if self.input_nodes:
+            raise ValueError("Remove inputs with `self.remove_input_nodes()` before performing Pauli presimulation.")
+        measure_pauli(self, copy=False, ignore_pauli_with_deps=ignore_pauli_with_deps)
 
     def draw_graph(
         self,
@@ -1652,23 +1664,17 @@ def __str__(self) -> str:
         assert_never(self.reason)
 
 
-def measure_pauli(
-    pattern: Pattern, leave_input: bool, *, copy: bool = False, ignore_pauli_with_deps: bool = False
-) -> Pattern:
+def measure_pauli(pattern: Pattern, *, copy: bool = False, ignore_pauli_with_deps: bool = False) -> Pattern:
     """Perform Pauli measurement of a pattern by fast graph state simulator.
 
-    Uses the decorated-graph method implemented in graphix.graphsim to perform
-    the measurements in Pauli bases, and then sort remaining nodes back into
-    pattern together with Clifford commands.
+    Uses the decorated-graph method implemented in graphix.graphsim to perform the measurements in Pauli bases, and then sort remaining nodes back into
+    pattern together with Clifford commands. Users are required to ensure there are no input nodes with :func:`graphix.pattern.Pattern.remove_input_nodes` before using this function.
 
     TODO: non-XY plane measurements in original pattern
 
     Parameters
     ----------
     pattern : graphix.pattern.Pattern object
-    leave_input : bool
-        True: input nodes will not be removed
-        False: all the nodes measured in Pauli bases will be removed
     copy : bool
         True: changes will be applied to new copied object and will be returned
         False: changes will be applied to the supplied Pattern object
@@ -1684,20 +1690,20 @@ def measure_pauli(
         only returned if copy argument is True.
 
 
+    .. seealso:: :class:`graphix.pattern.Pattern.remove_input_nodes`
     .. seealso:: :class:`graphix.graphsim.GraphState`
     """
+    pat = Pattern() if copy else pattern
     standardized_pattern = optimization.StandardizedPattern.from_pattern(pattern)
     if not ignore_pauli_with_deps:
         standardized_pattern = standardized_pattern.perform_pauli_pushing()
     output_nodes = set(pattern.output_nodes)
     graph = standardized_pattern.extract_graph()
     graph_state = GraphState(nodes=graph.nodes, edges=graph.edges, vops=standardized_pattern.c_dict)
-    results: dict[int, Outcome] = {}
-    to_measure, non_pauli_meas = pauli_nodes(standardized_pattern, leave_input)
-    if not leave_input and len(list(set(pattern.input_nodes) & {i[0].node for i in to_measure})) > 0:
-        new_inputs = []
-    else:
-        new_inputs = pattern.input_nodes
+    results: dict[int, Outcome] = pat.results
+    to_measure, non_pauli_meas = pauli_nodes(standardized_pattern)
+    if not to_measure:
+        return pattern
     for cmd in to_measure:
         pattern_cmd = cmd[0]
         measurement_basis = cmd[1]
@@ -1746,7 +1752,7 @@ def measure_pauli(
     # update command sequence
     vops = graph_state.get_vops()
     new_seq: list[Command] = []
-    new_seq.extend(command.N(node=index) for index in set(graph_state.nodes) - set(new_inputs))
+    new_seq.extend(command.N(node=index) for index in set(graph_state.nodes))
     new_seq.extend(command.E(nodes=edge) for edge in graph_state.edges)
     new_seq.extend(
         cmd.clifford(Clifford(vops[cmd.node])) for cmd in standardized_pattern.m_list if cmd.node in graph_state.nodes
@@ -1758,26 +1764,20 @@ def measure_pauli(
     )
     new_seq.extend(command.Z(node=node, domain=set(domain)) for node, domain in standardized_pattern.z_dict.items())
     new_seq.extend(command.X(node=node, domain=set(domain)) for node, domain in standardized_pattern.x_dict.items())
-
-    pat = Pattern() if copy else pattern
-
-    pat.replace(new_seq, input_nodes=new_inputs)
+    pat.replace(new_seq, input_nodes=[])
     pat.reorder_output_nodes(standardized_pattern.output_nodes)
     assert pat.n_node == len(graph_state.nodes)
     pat.results = results
     pat._pauli_preprocessed = True
     return pat
 
 
-def pauli_nodes(
-    pattern: optimization.StandardizedPattern, leave_input: bool
-) -> tuple[list[tuple[command.M, PauliMeasurement]], set[int]]:
+def pauli_nodes(pattern: optimization.StandardizedPattern) -> tuple[list[tuple[command.M, PauliMeasurement]], set[int]]:
     """Return the list of measurement commands that are in Pauli bases and that are not dependent on any non-Pauli measurements.
 
     Parameters
     ----------
     pattern : optimization.StandardizedPattern
-    leave_input : bool
 
     Returns
     -------
@@ -1790,7 +1790,7 @@ def pauli_nodes(
     non_pauli_node: set[int] = set()
     for cmd in pattern.m_list:
         pm = PauliMeasurement.try_from(cmd.plane, cmd.angle)  # None returned if the measurement is not in Pauli basis
-        if pm is not None and (cmd.node not in pattern.input_nodes or not leave_input):
+        if pm is not None:
             # Pauli measurement to be removed
             if pm.axis == Axis.X:
                 if cmd.t_domain & non_pauli_node:  # cmd depend on non-Pauli measurement

tests/test_gflow.py

@@ -538,6 +538,7 @@ def test_rand_circ_gflow(self, fx_rng: Generator) -> None:
         pattern = circ.transpile().pattern
         pattern.standardize()
         pattern.shift_signals()
+        pattern.remove_input_nodes()
         pattern.perform_pauli_measurements()
         graph = pattern.extract_graph()
         input_ = set()

tests/test_optimization.py

@@ -61,6 +61,7 @@ def test_incorporate_pauli_results(fx_bg: PCG64, jumps: int) -> None:
     pattern = circuit.transpile().pattern
     pattern.standardize()
     pattern.shift_signals()
+    pattern.remove_input_nodes()
     pattern.perform_pauli_measurements()
     pattern2 = incorporate_pauli_results(pattern)
     state = pattern.simulate_pattern(rng=rng)
@@ -78,6 +79,7 @@ def test_flow_after_pauli_preprocessing(fx_bg: PCG64, jumps: int) -> None:
     pattern.standardize()
     pattern.shift_signals()
     # pattern.move_pauli_measurements_to_the_front()
+    pattern.remove_input_nodes()
     pattern.perform_pauli_measurements()
     pattern2 = incorporate_pauli_results(pattern)
     gflow = pattern2.extract_gflow()
@@ -93,6 +95,7 @@ def test_remove_useless_domains(fx_bg: PCG64, jumps: int) -> None:
     pattern = circuit.transpile().pattern
     pattern.standardize()
     pattern.shift_signals()
+    pattern.remove_input_nodes()
     pattern.perform_pauli_measurements()
     pattern2 = remove_useless_domains(pattern)
     state = pattern.simulate_pattern(rng=rng)

tests/test_parameter.py

@@ -158,6 +158,7 @@ def test_random_circuit_with_parameters(fx_bg: PCG64, jumps: int, use_xreplace:
     pattern = circuit.transpile().pattern
     pattern.standardize()
     pattern.shift_signals()
+    pattern.remove_input_nodes()
     pattern.perform_pauli_measurements()
     pattern.minimize_space()
     assignment: dict[Parameter, float] = {alpha: rng.uniform(high=2), beta: rng.uniform(high=2)}