diff --git a/CHANGELOG.md b/CHANGELOG.md index fe318797..1de6d275 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -8,10 +8,21 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ## Unreleased ### Added - +- #358: Refactor of flow tools - Part I + - New module `graphix.flow.core` which introduces classes `PauliFlow`, `GFlow`, `CausalFlow` and `XZCorrections` allowing a finer analysis of MBQC flows. This module subsumes `graphix.generator` which has been removed and part of `graphix.gflow` which will be removed in the future. + - New module `graphix.flow._find_cflow` with the existing causal-flow finding algorithm. + - New module `graphix.flow._find_gpflow` with the existing g- and Pauli-flow finding algorithm introduced in #337. + - New abstract types `graphix.fundamentals.AbstractMeasurement` and `graphix.fundamentals.AbstractPlanarMeasurement` which serve as an umbrella of the existing types `graphix.measurements.Measurement`, `graphix.fundamentals.Plane` and `graphix.fundamentals.Axis`. + - New method `graphix.pattern.Pattern.extract_opengraph` which subsumes the static method `graphix.opengraph.OpenGraph.from_pattern`. + - New methods of `graphix.opengraph.OpenGraph` which allow to extract a causal, g- or Pauli flow. ### Fixed ### Changed + - #358: Refactor of flow tools - Part I + - API for the `graphix.opengraph.OpenGraph` class: + - `OpenGraphs` are parametrically typed so that they can be defined on planes and axes mappings in addition to measurements mappings. + - Attribute names are now `graph`, `input_nodes`, `output_nodes` and `measurements`. + ## [0.3.3] - 2025-10-23 diff --git a/docs/source/data.rst b/docs/source/data.rst index c8545d86..0e84328d 100644 --- a/docs/source/data.rst +++ b/docs/source/data.rst @@ -62,6 +62,18 @@ This module defines standard data structure for Pauli operators. .. autoclass:: Pauli +:mod:`graphix.measurements` module +++++++++++++++++++++++++++++++++++ + +This module defines data structures for single-qubit measurements in MBQC. + +.. automodule:: graphix.measurements + +.. currentmodule:: graphix.measurements + +.. autoclass:: Measurement + :members: + :mod:`graphix.instruction` module +++++++++++++++++++++++++++++++++ diff --git a/docs/source/generator.rst b/docs/source/generator.rst index 083e4d8e..1a991595 100644 --- a/docs/source/generator.rst +++ b/docs/source/generator.rst @@ -41,12 +41,3 @@ Pattern Generation .. autoclass:: TranspileResult .. autoclass:: SimulateResult - -:mod:`graphix.generator` module -+++++++++++++++++++++++++++++++ - -.. automodule:: graphix.generator - -.. currentmodule:: graphix.generator - -.. autofunction:: graphix.generator.generate_from_graph diff --git a/docs/source/open_graph.rst b/docs/source/open_graph.rst index 41f02da6..ab1dbdbf 100644 --- a/docs/source/open_graph.rst +++ b/docs/source/open_graph.rst @@ -9,5 +9,3 @@ This module defines classes for defining MBQC patterns as Open Graphs. .. currentmodule:: graphix.opengraph .. autoclass:: OpenGraph - -.. autoclass:: Measurement diff --git a/graphix/__init__.py b/graphix/__init__.py index 782a7ef7..cb8aedb6 100644 --- a/graphix/__init__.py +++ b/graphix/__init__.py @@ -2,10 +2,9 @@ from __future__ import annotations -from graphix.generator import generate_from_graph from graphix.graphsim import GraphState from graphix.pattern import Pattern from graphix.sim.statevec import Statevec from graphix.transpiler import Circuit -__all__ = ["Circuit", "GraphState", "Pattern", "Statevec", "generate_from_graph"] +__all__ = ["Circuit", "GraphState", "Pattern", "Statevec"] diff --git a/graphix/find_pflow.py b/graphix/find_pflow.py index 6a0b8e94..609c00cd 100644 --- a/graphix/find_pflow.py +++ b/graphix/find_pflow.py @@ -24,6 +24,7 @@ if TYPE_CHECKING: from collections.abc import Set as AbstractSet + from graphix.measurements import Measurement from graphix.opengraph import OpenGraph @@ -44,14 +45,14 @@ class OpenGraphIndex: At initialization, `non_outputs_optim` is a copy of `non_outputs`. The nodes corresponding to zero-rows of the order-demand matrix are removed for calculating the P matrix more efficiently in the `:func: _find_pflow_general` routine. """ - def __init__(self, og: OpenGraph) -> None: + def __init__(self, og: OpenGraph[Measurement]) -> None: self.og = og - nodes = set(og.inside.nodes) + nodes = set(og.graph.nodes) # Nodes don't need to be sorted. We do it for debugging purposes, so we can check the matrices in intermediate steps of the algorithm. - nodes_non_input = sorted(nodes - set(og.inputs)) - nodes_non_output = sorted(nodes - set(og.outputs)) + nodes_non_input = sorted(nodes - set(og.input_nodes)) + nodes_non_output = sorted(nodes - set(og.output_nodes)) self.non_inputs = NodeIndex() self.non_inputs.extend(nodes_non_input) @@ -84,7 +85,7 @@ def _compute_reduced_adj(ogi: OpenGraphIndex) -> MatGF2: See Definition 3.3 in Mitosek and Backens, 2024 (arXiv:2410.23439). """ - graph = ogi.og.inside + graph = ogi.og.graph row_tags = ogi.non_outputs col_tags = ogi.non_inputs @@ -119,7 +120,7 @@ def _compute_pflow_matrices(ogi: OpenGraphIndex) -> tuple[MatGF2, MatGF2]: flow_demand_matrix = _compute_reduced_adj(ogi) order_demand_matrix = flow_demand_matrix.copy() - inputs_set = set(ogi.og.inputs) + inputs_set = set(ogi.og.input_nodes) meas = ogi.og.measurements row_tags = ogi.non_outputs @@ -212,7 +213,7 @@ def _compute_p_matrix(ogi: OpenGraphIndex, nb_matrix: MatGF2) -> MatGF2 | None: See Theorem 4.4, steps 8 - 12 in Mitosek and Backens, 2024 (arXiv:2410.23439). """ n_no = len(ogi.non_outputs) # number of columns of P matrix. - n_oi_diff = len(ogi.og.outputs) - len(ogi.og.inputs) # number of rows of P matrix. + n_oi_diff = len(ogi.og.output_nodes) - len(ogi.og.input_nodes) # number of rows of P matrix. n_no_optim = len(ogi.non_outputs_optim) # number of rows and columns of the third block of the K_{LS} matrix. # Steps 8, 9 and 10 @@ -414,7 +415,7 @@ def _find_pflow_general(ogi: OpenGraphIndex) -> tuple[MatGF2, MatGF2] | None: See Theorem 4.4 and Algorithm 3 in Mitosek and Backens, 2024 (arXiv:2410.23439). """ n_no = len(ogi.non_outputs) - n_oi_diff = len(ogi.og.outputs) - len(ogi.og.inputs) + n_oi_diff = len(ogi.og.output_nodes) - len(ogi.og.input_nodes) # Steps 1 and 2 flow_demand_matrix, order_demand_matrix = _compute_pflow_matrices(ogi) @@ -552,7 +553,7 @@ def _cnc_matrices2pflow( if (topo_gen := _compute_topological_generations(ordering_matrix)) is None: return None # The NC matrix is not a DAG, therefore there's no flow. - l_k = dict.fromkeys(ogi.og.outputs, 0) # Output nodes are always in layer 0. + l_k = dict.fromkeys(ogi.og.output_nodes, 0) # Output nodes are always in layer 0. # If m >_c n, with >_c the flow order for two nodes m, n, then layer(n) > layer(m). # Therefore, we iterate the topological sort of the graph in _reverse_ order to obtain the order of measurements. @@ -570,7 +571,7 @@ def _cnc_matrices2pflow( return pf, l_k -def find_pflow(og: OpenGraph) -> tuple[dict[int, set[int]], dict[int, int]] | None: +def find_pflow(og: OpenGraph[Measurement]) -> tuple[dict[int, set[int]], dict[int, int]] | None: """Return a Pauli flow of the input open graph if it exists. Parameters @@ -592,8 +593,8 @@ def find_pflow(og: OpenGraph) -> tuple[dict[int, set[int]], dict[int, int]] | No ----- See Theorems 3.1, 4.2 and 4.4, and Algorithms 2 and 3 in Mitosek and Backens, 2024 (arXiv:2410.23439). """ - ni = len(og.inputs) - no = len(og.outputs) + ni = len(og.input_nodes) + no = len(og.output_nodes) if ni > no: return None diff --git a/graphix/flow/__init__.py b/graphix/flow/__init__.py new file mode 100644 index 00000000..cf9a3529 --- /dev/null +++ b/graphix/flow/__init__.py @@ -0,0 +1 @@ +"""Flow classes and flow-finding algorithms.""" diff --git a/graphix/flow/_find_cflow.py b/graphix/flow/_find_cflow.py new file mode 100644 index 00000000..56cb375c --- /dev/null +++ b/graphix/flow/_find_cflow.py @@ -0,0 +1,120 @@ +"""Causal flow finding algorithm. + +This module implements Algorithm 1 from Ref. [1]. For a given labelled open graph (G, I, O, meas_plane), this algorithm finds a causal flow [2] in polynomial time with the number of nodes, :math:`O(N^2)`. + +References +---------- +[1] Mhalla and Perdrix, (2008), Finding Optimal Flows Efficiently, doi.org/10.1007/978-3-540-70575-8_70 +[2] Browne et al., 2007 New J. Phys. 9 250 (arXiv:quant-ph/0702212) +""" + +from __future__ import annotations + +from typing import TYPE_CHECKING + +from graphix.flow.core import CausalFlow +from graphix.fundamentals import Plane + +if TYPE_CHECKING: + from collections.abc import Set as AbstractSet + + from graphix.opengraph import OpenGraph, _PM_co + + +def find_cflow(og: OpenGraph[_PM_co]) -> CausalFlow[_PM_co] | None: + """Return the causal flow of the input open graph if it exists. + + Parameters + ---------- + og : OpenGraph[_PM_co] + Open graph whose causal flow is calculated. + + Returns + ------- + CausalFlow[_PM_co] | None + A causal flow object if the open graph has causal flow, `None` otherwise. + + Notes + ----- + - See Definition 2, Theorem 1 and Algorithm 1 in Ref. [1]. + - The open graph instance must be of parametric type `Measurement` or `Plane` since the causal flow is only defined on open graphs with :math:`XY` measurements. + + References + ---------- + [1] Mhalla and Perdrix, (2008), Finding Optimal Flows Efficiently, doi.org/10.1007/978-3-540-70575-8_70 + """ + for measurement in og.measurements.values(): + if measurement.to_plane() in {Plane.XZ, Plane.YZ}: + return None + + corrected_nodes = set(og.output_nodes) + corrector_candidates = corrected_nodes - set(og.input_nodes) + non_input_nodes = og.graph.nodes - set(og.input_nodes) + + cf: dict[int, frozenset[int]] = {} + # Output nodes are always in layer 0. If the open graph has flow, it must have outputs, so we never end up with an empty set at `layers[0]`. + layers: list[frozenset[int]] = [ + frozenset(corrected_nodes) + ] # A copy is necessary because `corrected_nodes` is mutable and changes during the algorithm. + + return _flow_aux(og, non_input_nodes, corrected_nodes, corrector_candidates, cf, layers) + + +def _flow_aux( + og: OpenGraph[_PM_co], + non_input_nodes: AbstractSet[int], + corrected_nodes: AbstractSet[int], + corrector_candidates: AbstractSet[int], + cf: dict[int, frozenset[int]], + layers: list[frozenset[int]], +) -> CausalFlow[_PM_co] | None: + """Find one layer of the causal flow. + + Parameters + ---------- + og : OpenGraph[_PM_co] + Open graph whose causal flow is calculated. + non_input_nodes : AbstractSet[int] + Non-input nodes of the input open graph. This parameter remains constant throughout the execution of the algorithm and can be derived from `og` at any time. It is passed as an argument to avoid unnecessary recalculations. + corrected_nodes : AbstractSet[int] + Nodes which have already been corrected. + corrector_candidates : AbstractSet[int] + Nodes which could correct a node at the time of calling the function. This set can never contain input nodes, uncorrected nodes or nodes which already correct another node. + cf : dict[int, frozenset[int]] + Causal flow correction function. `cf[i]` is the one-qubit set correcting the measurement of qubit `i`. + layers : list[frozenset[int]] + Partial order between corrected qubits in a layer form. The set `layers[i]` comprises the nodes in layer `i`. Nodes in layer `i` are "larger" in the partial order than nodes in layer `i+1`. + + + Returns + ------- + CausalFlow[_PM_co] | None + A causal flow object if the open graph has causal flow, `None` otherwise. + """ + corrected_nodes_new: set[int] = set() + corrector_nodes_new: set[int] = set() + curr_layer: set[int] = set() + + non_corrected_nodes = og.graph.nodes - corrected_nodes + + if corrected_nodes == set(og.graph.nodes): + return CausalFlow(og, cf, tuple(layers)) + + for p in corrector_candidates: + non_corrected_neighbors = og.neighbors({p}) & non_corrected_nodes + if len(non_corrected_neighbors) == 1: + (q,) = non_corrected_neighbors + cf[q] = frozenset({p}) + curr_layer.add(q) + corrected_nodes_new |= {q} + corrector_nodes_new |= {p} + + layers.append(frozenset(curr_layer)) + + if len(corrected_nodes_new) == 0: + return None + + corrected_nodes |= corrected_nodes_new + corrector_candidates = (corrector_candidates - corrector_nodes_new) | (corrected_nodes_new & non_input_nodes) + + return _flow_aux(og, non_input_nodes, corrected_nodes, corrector_candidates, cf, layers) diff --git a/graphix/flow/_find_gpflow.py b/graphix/flow/_find_gpflow.py new file mode 100644 index 00000000..65e3ac91 --- /dev/null +++ b/graphix/flow/_find_gpflow.py @@ -0,0 +1,670 @@ +"""Pauli flow finding algorithm. + +This module implements the algorithm presented in [1]. For a given labelled open graph (G, I, O, meas_plane), this algorithm finds a maximally delayed Pauli flow [2] in polynomial time with the number of nodes, :math:`O(N^3)`. +If the input graph does not have Pauli measurements, the algorithm returns a generalised flow (gflow) if it exists by definition. + +References +---------- +[1] Mitosek and Backens, 2024 (arXiv:2410.23439). +[2] Browne et al., 2007 New J. Phys. 9 250 (arXiv:quant-ph/0702212) +""" + +from __future__ import annotations + +from copy import deepcopy +from dataclasses import dataclass +from functools import cached_property +from typing import TYPE_CHECKING, Generic, TypeVar + +import numpy as np +from typing_extensions import override + +from graphix._linalg import MatGF2, solve_f2_linear_system +from graphix.fundamentals import AbstractMeasurement, AbstractPlanarMeasurement, Axis, Plane +from graphix.sim.base_backend import NodeIndex + +if TYPE_CHECKING: + from collections.abc import Set as AbstractSet + + from graphix.opengraph import OpenGraph + + +_M_co = TypeVar("_M_co", bound=AbstractMeasurement, covariant=True) +_PM_co = TypeVar("_PM_co", bound=AbstractPlanarMeasurement, covariant=True) + + +class AlgebraicOpenGraph(Generic[_M_co]): + """A class for providing an algebraic representation of open graphs as introduced in [1]. In particular, it allows managing the mapping between node labels of the graph and the relevant matrix indices. The flow-demand and order-demand matrices are cached properties. + + It reuses the class `:class: graphix.sim.base_backend.NodeIndex` introduced for managing the mapping between node numbers and qubit indices in the internal state of the backend. + + Attributes + ---------- + og (OpenGraph) + non_inputs (NodeIndex) : Mapping between matrix indices and non-input nodes (labelled with integers). + non_outputs (NodeIndex) : Mapping between matrix indices and non-output nodes (labelled with integers). + non_outputs_optim (NodeIndex) : Mapping between matrix indices and a subset of non-output nodes (labelled with integers). + + Notes + ----- + At initialization, `non_outputs_optim` is a copy of `non_outputs`. The nodes corresponding to zero-rows of the order-demand matrix are removed for calculating the :math:`P` matrix more efficiently in the `:func: _compute_correction_matrix_general` routine. + + References + ---------- + [1] Mitosek and Backens, 2024 (arXiv:2410.23439). + """ + + def __init__(self, og: OpenGraph[_M_co]) -> None: + """Initialize AlgebraicOpenGraph objects. + + Parameters + ---------- + og : OpenGraph[_M_co] + The open graph in its standard representation. + """ + self.og = og + nodes = set(og.graph.nodes) + + # Nodes don't need to be sorted. We do it for debugging purposes, so we can check the matrices in intermediate steps of the algorithm. + + nodes_non_input = sorted(nodes - set(og.input_nodes)) + nodes_non_output = sorted(nodes - set(og.output_nodes)) + + self.non_inputs = NodeIndex() + self.non_inputs.extend(nodes_non_input) + + self.non_outputs = NodeIndex() + self.non_outputs.extend(nodes_non_output) + + # Needs to be a deep copy because it may be modified during runtime. + self.non_outputs_optim = deepcopy(self.non_outputs) + + @property + def flow_demand_matrix(self) -> MatGF2: + """Return the flow-demand matrix. + + Returns + ------- + MatGF2 + Flow-demand matrix + + Notes + ----- + See Definition 3.4 and Algorithm 1 in Mitosek and Backens, 2024 (arXiv:2410.23439). + """ + return self._og_matrices[0] + + @property + def order_demand_matrix(self) -> MatGF2: + """Return the flow-demand matrix. + + Returns + ------- + MatGF2 + Order-demand matrix + + Notes + ----- + See Definition 3.5 and Algorithm 1 in Mitosek and Backens, 2024 (arXiv:2410.23439). + """ + return self._og_matrices[1] + + def _compute_reduced_adj(self) -> MatGF2: + r"""Return the reduced adjacency matrix (RAdj) of the open graph. + + Returns + ------- + adj_red : MatGF2 + Reduced adjacency matrix. + + Notes + ----- + The adjacency matrix of a graph :math:`Adj_G` is an :math:`n \times n` matrix. + + The RAdj matrix of an open graph OG is an :math:`(n - n_O) \times (n - n_I)` submatrix of :math:`Adj_G` constructed by removing the output rows and input columns of :math:`Adj_G`. + + See Definition 3.3 in Mitosek and Backens, 2024 (arXiv:2410.23439). + """ + graph = self.og.graph + row_tags = self.non_outputs + col_tags = self.non_inputs + + adj_red = np.zeros((len(row_tags), len(col_tags)), dtype=np.uint8).view(MatGF2) + + for n1, n2 in graph.edges: + for u, v in ((n1, n2), (n2, n1)): + if u in row_tags and v in col_tags: + i, j = row_tags.index(u), col_tags.index(v) + adj_red[i, j] = 1 + + return adj_red + + @cached_property + def _og_matrices(self) -> tuple[MatGF2, MatGF2]: + r"""Construct the flow-demand and order-demand matrices. + + Returns + ------- + flow_demand_matrix : MatGF2 + order_demand_matrix : MatGF2 + + Notes + ----- + - See Definitions 3.4 and 3.5, and Algorithm 1 in Mitosek and Backens, 2024 (arXiv:2410.23439). + """ + flow_demand_matrix = self._compute_reduced_adj() + order_demand_matrix = flow_demand_matrix.copy() + + inputs_set = set(self.og.input_nodes) + + row_tags = self.non_outputs + col_tags = self.non_inputs + + for v in row_tags: # v is a node tag + i = row_tags.index(v) + plane_axis_v = self._get_measurement_label(v) + + if plane_axis_v in {Plane.YZ, Plane.XZ, Axis.Z}: + flow_demand_matrix[i, :] = 0 # Set row corresponding to node v to 0 + if plane_axis_v in {Plane.YZ, Plane.XZ, Axis.Y, Axis.Z} and v not in inputs_set: + j = col_tags.index(v) + flow_demand_matrix[i, j] = 1 # Set element (v, v) = 0 + if plane_axis_v in {Plane.XY, Axis.X, Axis.Y, Axis.Z}: + order_demand_matrix[i, :] = 0 # Set row corresponding to node v to 0 + if plane_axis_v in {Plane.XY, Plane.XZ} and v not in inputs_set: + j = col_tags.index(v) + order_demand_matrix[i, j] = 1 # Set element (v, v) = 1 + + return flow_demand_matrix, order_demand_matrix + + def _get_measurement_label(self, node: int) -> Plane | Axis: + """Return the measurement label (plane or axis) of a node in the open graph. + + Parameters + ---------- + node : int + Measured node. + + Returns + ------- + Plane | Axis + Measurement label. + + Notes + ----- + Measurements with a Pauli angle are intepreted as `Axis` instances. + """ + return self.og.measurements[node].to_plane_or_axis() + + +class PlanarAlgebraicOpenGraph(AlgebraicOpenGraph[_PM_co]): + """A subclass of `AlgebraicOpenGraph`. + + This class differs from its parent class only in that Pauli measurements are interpreted as `Plane` instances (instead of `Axis`) when constructing the flow-demand and order-demand matrices. This allows to verify if open graphs with measurements along Pauli angles interpreted as planes have generalised flow. + + """ + + @override + def _get_measurement_label(self, node: int) -> Plane: + """Return the measurement label (plane) of a node in the open graph. + + Parameters + ---------- + node : int + Measured node. + + Returns + ------- + Plane + Measurement label. + + Notes + ----- + Measurements with a Pauli angle are intepreted as `Plane` instances. + """ + return self.og.measurements[node].to_plane() + + +@dataclass(frozen=True) # `NamedTuple` does not support multiple inheritance in Python 3.9 and 3.10 +class CorrectionMatrix(Generic[_M_co]): + r"""A dataclass to bundle the correction matrix and its associated open graph. + + Attributes + ---------- + aog (AgebraicOpenGraph) : Open graph in an algebraic representation. + c_matrix (MatGF2) : Matrix encoding the correction function of a Pauli (or generalised) flow, :math:`C`. + + Notes + ----- + The correction matrix :math:`C` is an :math:`(n - n_I) \times (n - n_O)` matrix related to the correction function :math:`c(v) = \{u \in I^c|C_{u,v} = 1\}`, where :math:`I^c` are the non-input nodes of `aog`. In other words, the column :math:`v` of :math:`C` encodes the correction set of :math:`v`, :math:`c(v)`. + + See Definition 3.6 in Mitosek and Backens, 2024 (arXiv:2410.23439). + """ + + aog: AlgebraicOpenGraph[_M_co] + c_matrix: MatGF2 + + def to_correction_function(self) -> dict[int, frozenset[int]]: + r"""Transform the correction matrix into a correction function. + + Returns + ------- + correction_function : dict[int, frozenset[int]] + Pauli (or generalised) flow correction function. `correction_function[i]` is the set of qubits correcting the measurement of qubit `i`. + """ + row_tags = self.aog.non_inputs + col_tags = self.aog.non_outputs + correction_function: dict[int, frozenset[int]] = {} + for node in col_tags: + i = col_tags.index(node) + correction_set = {row_tags[j] for j in np.flatnonzero(self.c_matrix[:, i])} + correction_function[node] = frozenset(correction_set) + return correction_function + + +def _compute_p_matrix(aog: AlgebraicOpenGraph[_M_co], nb_matrix: MatGF2) -> MatGF2 | None: + r"""Perform the steps 8 - 12 of the general case (larger number of outputs than inputs) algorithm. + + Parameters + ---------- + aog : AlgebraicOpenGraph + Open graph for which the matrix :math:`P` is computed. + nb_matrix : MatGF2 + Matrix :math:`N_B` + + Returns + ------- + p_matrix : MatGF2 + Matrix encoding the correction function. + + or `None` + if the input open graph does not have Pauli flow. + + Notes + ----- + See Theorem 4.4, steps 8 - 12 in Mitosek and Backens, 2024 (arXiv:2410.23439). + """ + n_no = len(aog.non_outputs) # number of columns of P matrix. + n_oi_diff = len(aog.og.output_nodes) - len(aog.og.input_nodes) # number of rows of P matrix. + n_no_optim = len(aog.non_outputs_optim) # number of rows and columns of the third block of the K_{LS} matrix. + + # Steps 8, 9 and 10 + kils_matrix = np.concatenate( + (nb_matrix[:, n_no:], nb_matrix[:, :n_no], np.eye(n_no_optim, dtype=np.uint8)), axis=1 + ).view(MatGF2) # N_R | N_L | 1 matrix. + kls_matrix = kils_matrix.gauss_elimination(ncols=n_oi_diff, copy=True) # RREF form is not needed, only REF. + + # Step 11 + p_matrix = np.zeros((n_oi_diff, n_no), dtype=np.uint8).view(MatGF2) + solved_nodes: set[int] = set() + non_outputs_set = set(aog.non_outputs) + + # Step 12 + while solved_nodes != non_outputs_set: + solvable_nodes = _find_solvable_nodes(aog, kls_matrix, non_outputs_set, solved_nodes, n_oi_diff) # Step 12.a + if not solvable_nodes: + return None + + _update_p_matrix(aog, kls_matrix, p_matrix, solvable_nodes, n_oi_diff) # Steps 12.b, 12.c + _update_kls_matrix(aog, kls_matrix, kils_matrix, solvable_nodes, n_oi_diff, n_no, n_no_optim) # Step 12.d + solved_nodes.update(solvable_nodes) + + return p_matrix + + +def _find_solvable_nodes( + aog: AlgebraicOpenGraph[_M_co], + kls_matrix: MatGF2, + non_outputs_set: AbstractSet[int], + solved_nodes: AbstractSet[int], + n_oi_diff: int, +) -> set[int]: + """Return the set nodes whose associated linear system is solvable. + + A node is solvable if: + - It has not been solved yet. + - Its column in the second block of :math:`K_{LS}` (which determines the constants in each equation) has only zeros where it intersects rows for which all the coefficients in the first block are 0s. + + See Theorem 4.4, step 12.a in Mitosek and Backens, 2024 (arXiv:2410.23439). + """ + solvable_nodes: set[int] = set() + + row_idxs = np.flatnonzero( + ~kls_matrix[:, :n_oi_diff].any(axis=1) + ) # Row indices of the 0-rows in the first block of K_{LS}. + if row_idxs.size: + for v in non_outputs_set - solved_nodes: + j = n_oi_diff + aog.non_outputs.index(v) # `n_oi_diff` is the column offset from the first block of K_{LS}. + if not kls_matrix[row_idxs, j].any(): + solvable_nodes.add(v) + else: + # If the first block of K_{LS} does not have 0-rows, all non-solved nodes are solvable. + solvable_nodes = set(non_outputs_set - solved_nodes) + + return solvable_nodes + + +def _update_p_matrix( + aog: AlgebraicOpenGraph[_M_co], + kls_matrix: MatGF2, + p_matrix: MatGF2, + solvable_nodes: AbstractSet[int], + n_oi_diff: int, +) -> None: + """Update `p_matrix`. + + The solution of the linear system associated with node :math:`v` in `solvable_nodes` corresponds to the column of `p_matrix` associated with node :math:`v`. + + See Theorem 4.4, steps 12.b and 12.c in Mitosek and Backens, 2024 (arXiv:2410.23439). + """ + for v in solvable_nodes: + j = aog.non_outputs.index(v) + j_shift = n_oi_diff + j # `n_oi_diff` is the column offset from the first block of K_{LS}. + mat = MatGF2(kls_matrix[:, :n_oi_diff]) # First block of K_{LS}, in row echelon form. + b = MatGF2(kls_matrix[:, j_shift]) + x = solve_f2_linear_system(mat, b) + p_matrix[:, j] = x + + +def _update_kls_matrix( + aog: AlgebraicOpenGraph[_M_co], + kls_matrix: MatGF2, + kils_matrix: MatGF2, + solvable_nodes: AbstractSet[int], + n_oi_diff: int, + n_no: int, + n_no_optim: int, +) -> None: + """Update `kls_matrix`. + + Bring the linear system encoded in :math:`K_{LS}` to the row-echelon form (REF) that would be achieved by Gaussian elimination if the row and column vectors corresponding to vertices in `solvable_nodes` where not included in the starting matrix. + + See Theorem 4.4, step 12.d in Mitosek and Backens, 2024 (arXiv:2410.23439). + """ + shift = n_oi_diff + n_no # `n_oi_diff` + `n_no` is the column offset from the first two blocks of K_{LS}. + row_permutation: list[int] + + def reorder(old_pos: int, new_pos: int) -> None: # Used in step 12.d.vi + """Reorder the elements of `row_permutation`. + + The element at `old_pos` is placed on the right of the element at `new_pos`. + Example: + ``` + row_permutation = [0, 1, 2, 3, 4] + reorder(1, 3) -> [0, 2, 3, 1, 4] + reorder(2, -1) -> [2, 0, 1, 3, 4] + ``` + """ + val = row_permutation.pop(old_pos) + row_permutation.insert(new_pos + (new_pos < old_pos), val) + + for v in solvable_nodes: + if ( + v in aog.non_outputs_optim + ): # if `v` corresponded to a zero row in N_B, it was not present in `kls_matrix` because we removed it in the optimization process, so there's no need to do Gaussian elimination for that vertex. + # Step 12.d.ii + j = aog.non_outputs_optim.index(v) + j_shift = shift + j + row_idxs = np.flatnonzero( + kls_matrix[:, j_shift] + ).tolist() # Row indices with 1s in column of node `v` in third block. + + # `row_idxs` can't be empty: + # The third block of K_{LS} is initially the identity matrix, so all columns have initially a 1. Row permutations and row additions in the Gaussian elimination routine can't remove all 1s from a given column. + k = row_idxs.pop() + + # Step 12.d.iii + kls_matrix[row_idxs] ^= kls_matrix[k] # Adding a row to previous rows preserves REF. + + # Step 12.d.iv + kls_matrix[k] ^= kils_matrix[j] # Row `k` may now break REF. + + # Step 12.d.v + pivots: list[np.int_] = [] # Store pivots for next step. + for i, row in enumerate(kls_matrix): + if i != k: + col_idxs = np.flatnonzero(row[:n_oi_diff]) # Column indices with 1s in first block. + if col_idxs.size == 0: + # Row `i` has all zeros in the first block. Only row `k` can break REF, so rows below have all zeros in the first block too. + break + pivots.append(p := col_idxs[0]) + if kls_matrix[k, p]: # Row `k` has a 1 in the column corresponding to the leading 1 of row `i`. + kls_matrix[k] ^= row + + row_permutation = list(range(n_no_optim)) # Row indices of `kls_matrix`. + n_pivots = len(pivots) + + col_idxs = np.flatnonzero(kls_matrix[k, :n_oi_diff]) + pk = col_idxs[0] if col_idxs.size else None # Pivot of row `k`. + + if pk and k >= n_pivots: # Row `k` is non-zero in the FB (first block) and it's among zero rows. + # Find row `new_pos` s.t. `pivots[new_pos] <= pk < pivots[new_pos+1]`. + new_pos = ( + int(np.argmax(np.array(pivots) > pk) - 1) if pivots else -1 + ) # `pivots` can be empty. If so, we bring row `k` to the top since it's non-zero. + elif pk: # Row `k` is non-zero in the FB and it's among non-zero rows. + # Find row `new_pos` s.t. `pivots[new_pos] <= pk < pivots[new_pos+1]` + new_pos = int(np.argmax(np.array(pivots) > pk) - 1) + # We skipped row `k` in loop of step 12.d.v, so `pivots[j]` can be the pivot of row `j` or `j+1`. + if new_pos >= k: + new_pos += 1 + elif k < n_pivots: # Row `k` is zero in the first block and it's among non-zero rows. + new_pos = ( + n_pivots # Move row `k` to the top of the zeros block (i.e., below the row of the last pivot). + ) + else: # Row `k` is zero in the first block and it's among zero rows. + new_pos = k # Do nothing. + + if new_pos != k: + reorder(k, new_pos) # Modify `row_permutation` in-place. + kls_matrix[:] = kls_matrix[ + row_permutation + ] # `[:]` is crucial to modify the data pointed by `kls_matrix`. + + +def _compute_correction_matrix_general_case( + aog: AlgebraicOpenGraph[_M_co], flow_demand_matrix: MatGF2, order_demand_matrix: MatGF2 +) -> MatGF2 | None: + r"""Construct the generalized correction matrix :math:`C'C^B` for an open graph with larger number of outputs than inputs. + + Parameters + ---------- + aog : AlgebraicOpenGraph + Open graph for which :math:`C'C^B` and :math:`NC'C^B` are computed. + flow_demand_matrix: MatGF2 + Flow demand matrix :math:`M` (a property of the open graph). + order_demand_matrix: MatGF2 + Order demand matrix :math:`N` (a property of the open graph). + + Returns + ------- + correction_matrix : MatGF2 + Matrix encoding the correction function. + + or `None` + if the input open graph does not have Pauli flow. + + Notes + ----- + - The function returns `None` if + a) The flow-demand matrix is not invertible, or + b) Not all linear systems of equations associated to the non-output nodes are solvable, + meaning that `aog` does not have Pauli flow. + Condition (b) is satisfied when the flow-demand matrix :math:`M` does not have a right inverse :math:`C` such that :math:`NC` represents a directed acyclical graph (DAG). + + See Theorem 4.4 and Algorithm 3 in Mitosek and Backens, 2024 (arXiv:2410.23439). + """ + n_no = len(aog.non_outputs) + n_oi_diff = len(aog.og.output_nodes) - len(aog.og.input_nodes) + + # Steps 3 and 4 + correction_matrix_0 = flow_demand_matrix.right_inverse() # C0 matrix. + if correction_matrix_0 is None: + return None # The flow-demand matrix is not invertible, therefore there's no flow. + + # Steps 5, 6 and 7 + ker_flow_demand_matrix = flow_demand_matrix.null_space().transpose() # F matrix. + c_prime_matrix = np.concatenate((correction_matrix_0, ker_flow_demand_matrix), axis=1).view(MatGF2) + + row_idxs = np.flatnonzero(order_demand_matrix.any(axis=1)) # Row indices of the non-zero rows. + + if row_idxs.size: + # The p-matrix finding algorithm runs on the `order_demand_matrix` without the zero rows. + # This optimization is allowed because: + # - The zero rows remain zero after the change of basis (multiplication by `c_prime_matrix`). + # - The zero rows remain zero after gaussian elimination. + # - Removing the zero rows does not change the solvability condition of the open graph nodes. + nb_matrix_optim = ( + order_demand_matrix[row_idxs].view(MatGF2).mat_mul(c_prime_matrix) + ) # `view` is used to keep mypy happy without copying data. + for i in set(range(order_demand_matrix.shape[0])).difference(row_idxs): + aog.non_outputs_optim.remove(aog.non_outputs[i]) # Update the node-index mapping. + + # Steps 8 - 12 + if (p_matrix := _compute_p_matrix(aog, nb_matrix_optim)) is None: + return None + else: + # If all rows of `order_demand_matrix` are zero, any matrix will solve the associated linear system of equations. + p_matrix = np.zeros((n_oi_diff, n_no), dtype=np.uint8).view(MatGF2) + + # Step 13 + cb_matrix = np.concatenate((np.eye(n_no, dtype=np.uint8), p_matrix), axis=0).view(MatGF2) + + return c_prime_matrix.mat_mul(cb_matrix) + + +def _compute_topological_generations(ordering_matrix: MatGF2) -> list[list[int]] | None: + """Stratify the directed acyclic graph (DAG) represented by the ordering matrix into generations. + + Parameters + ---------- + ordering_matrix : MatGF2 + Matrix encoding the partial ordering between nodes interpreted as the adjacency matrix of a directed graph. + + Returns + ------- + list[list[int]] + topological generations. Integers represent the indices of the matrix `ordering_matrix`, not the labelling of the nodes. + + or `None` + if `ordering_matrix` is not a DAG. + + Notes + ----- + This function is adapted from `:func: networkx.algorithms.dag.topological_generations` so that it works directly on the adjacency matrix (which is the output of the Pauli-flow finding algorithm) instead of a `:class: nx.DiGraph` object. This avoids calling the function `nx.from_numpy_array` which can be expensive for certain graph instances. + + Here we use the convention that the element `ordering_matrix[i,j]` represents a link `j -> i`. NetworkX uses the opposite convention. + """ + adj_mat = ordering_matrix + + indegree_map: dict[int, int] = {} + zero_indegree: list[int] = [] + neighbors = {node: set(np.flatnonzero(row).astype(int)) for node, row in enumerate(adj_mat.T)} + for node, col in enumerate(adj_mat): + parents = np.flatnonzero(col) + if parents.size: + indegree_map[node] = parents.size + else: + zero_indegree.append(node) + + generations: list[list[int]] = [] + + while zero_indegree: + this_generation = zero_indegree + zero_indegree = [] + for node in this_generation: + for child in neighbors[node]: + indegree_map[child] -= 1 + if indegree_map[child] == 0: + zero_indegree.append(child) + del indegree_map[child] + generations.append(this_generation) + + if indegree_map: + return None + return generations + + +def compute_partial_order_layers(correction_matrix: CorrectionMatrix[_M_co]) -> tuple[frozenset[int], ...] | None: + r"""Compute the partial order compatible with the correction matrix if it exists. + + Parameters + ---------- + correction_matrix : CorrectionMatrix[_M_co] + Algebraic representation of the correction function. + + Returns + ------- + layers : tuple[frozenset[int], ...] + Partial order between corrected qubits in a layer form. The frozenset `layers[i]` comprises the nodes in layer `i`. Nodes in layer `i` are "larger" in the partial order than nodes in layer `i+1`. Output nodes are always in layer 0. + + or `None` + If the correction matrix is not compatible with a partial order on the the open graph, in which case the associated ordering matrix is not a DAG. In the context of the flow-finding algorithm, this means that the input open graph does not have Pauli (or generalised) flow. + + Notes + ----- + - The partial order of the Pauli (or generalised) flow :math:`<_c` is the transitive closure of :math:`\lhd_c`, where the latter is related to the ordering matrix :math:`NC` as :math:`v \lhd_c w \Leftrightarrow (NC)_{w,v} = 1`, for :math:`v, w, \in O^c` two non-output nodes of `aog`. The ordering matrix is the product of the order-demand and the correction matrices and it is the adjacency matrix of the directed acyclical graph encoding the partial order. + + - If the open graph has flow, it must have outputs, so `layers[0]` always contains a finite set of nodes. + + See Lemma 3.12, and Theorem 3.1 in Mitosek and Backens, 2024 (arXiv:2410.23439). + """ + aog, c_matrix = correction_matrix.aog, correction_matrix.c_matrix + ordering_matrix = aog.order_demand_matrix.mat_mul(c_matrix) + + if (topo_gen := _compute_topological_generations(ordering_matrix)) is None: + return None # The NC matrix is not a DAG, therefore there's no flow. + + layers = [ + frozenset(aog.og.output_nodes) + ] # Output nodes are always in layer 0. If the open graph has flow, it must have outputs, so we never end up with an empty set at `layers[0]`. + + # If m >_c n, with >_c the flow partial order for two nodes m, n, then layer(n) > layer(m). + # Therefore, we iterate the topological sort of the graph in _reverse_ order to obtain the order of measurements. + col_tags = aog.non_outputs + layers.extend(frozenset({col_tags[i] for i in idx_layer}) for idx_layer in reversed(topo_gen)) + + return tuple(layers) + + +def compute_correction_matrix(aog: AlgebraicOpenGraph[_M_co]) -> CorrectionMatrix[_M_co] | None: + """Return the correction matrix of the input open graph if it exists. + + Parameters + ---------- + aog : AlgebraicOpenGraph[_M_co] + Algberaic representation of the open graph whose correction matrix is calculated. + + Returns + ------- + correction_matrix : CorrectionMatrix[_M_co] + Algebraic representation of the correction function. + + or `None` + if the input open graph does not have Pauli (or generalised) flow. + + Notes + ----- + - In the case of open graphs with equal number of inputs and outputs, the function only returns `None` when the flow-demand matrix is not invertible (meaning that `aog` does not have Pauli flow). The additional condition for the existence of Pauli flow that the ordering matrix :math:`NC` must encode a directed acyclic graph (DAG) is verified by :func:`compute_partial_order`, which is called from the `graphix.flow.core.PauliFlow` constructor. + + - See Definitions 3.4, 3.5 and 3.6, Theorems 3.1, 4.2 and 4.4, and Algorithms 2 and 3 in Mitosek and Backens, 2024 (arXiv:2410.23439). + """ + ni = len(aog.og.input_nodes) + no = len(aog.og.output_nodes) + + if ni > no: + return None + + # Steps 1 and 2 + # Flow-demand and order-demand matrices are cached properties of `aog`. + flow_demand_matrix, order_demand_matrix = aog._og_matrices + + if ni == no: + correction_matrix = flow_demand_matrix.right_inverse() + else: + correction_matrix = _compute_correction_matrix_general_case(aog, flow_demand_matrix, order_demand_matrix) + + if correction_matrix is None: + return None + + return CorrectionMatrix(aog, correction_matrix) diff --git a/graphix/flow/core.py b/graphix/flow/core.py new file mode 100644 index 00000000..94c88008 --- /dev/null +++ b/graphix/flow/core.py @@ -0,0 +1,496 @@ +"""Class for flow objects and XZ-corrections.""" + +from __future__ import annotations + +from collections.abc import Sequence +from copy import copy +from dataclasses import dataclass +from typing import TYPE_CHECKING, Generic + +import networkx as nx + +# override introduced in Python 3.12 +from typing_extensions import override + +import graphix.pattern +from graphix.command import E, M, N, X, Z +from graphix.flow._find_gpflow import CorrectionMatrix, _M_co, _PM_co, compute_partial_order_layers + +if TYPE_CHECKING: + from collections.abc import Mapping + from collections.abc import Set as AbstractSet + from typing import Self + + from graphix.measurements import Measurement + from graphix.opengraph import OpenGraph + from graphix.pattern import Pattern + +TotalOrder = Sequence[int] + + +@dataclass(frozen=True) +class XZCorrections(Generic[_M_co]): + """An unmutable dataclass providing a representation of XZ-corrections. + + Attributes + ---------- + og : OpenGraph[_M_co] + The open graph with respect to which the XZ-corrections are defined. + x_corrections : Mapping[int, AbstractSet[int]] + Mapping of X-corrections: in each (`key`, `value`) pair, `key` is a measured node, and `value` is the set of nodes on which an X-correction must be applied depending on the measurement result of `key`. + z_corrections : Mapping[int, AbstractSet[int]] + Mapping of Z-corrections: in each (`key`, `value`) pair, `key` is a measured node, and `value` is the set of nodes on which an Z-correction must be applied depending on the measurement result of `key`. + partial_order_layers : Sequence[AbstractSet[int]] + Partial order between the open graph's nodes in a layer form determined by the corrections. The set `layers[i]` comprises the nodes in layer `i`. Nodes in layer `i` are "larger" in the partial order than nodes in layer `i+1`. If the open graph has output nodes, they are always in layer 0. Non-corrected, measured nodes are always in the last layer. + + Notes + ----- + The XZ-corrections mappings define a partial order, therefore, only `og`, `x_corrections` and `z_corrections` are necessary to initialize an `XZCorrections` instance (see :func:`XZCorrections.from_measured_nodes_mapping`). However, XZ-corrections are often extracted from a flow whose partial order is known and can be used to construct a pattern, so it can also be passed as an argument to the `dataclass` constructor. The correctness of the input parameters is not verified automatically. + + """ + + og: OpenGraph[_M_co] + x_corrections: Mapping[int, AbstractSet[int]] # {domain: nodes} + z_corrections: Mapping[int, AbstractSet[int]] # {domain: nodes} + partial_order_layers: Sequence[AbstractSet[int]] + + @staticmethod + def from_measured_nodes_mapping( + og: OpenGraph[_M_co], + x_corrections: Mapping[int, AbstractSet[int]] | None = None, + z_corrections: Mapping[int, AbstractSet[int]] | None = None, + ) -> XZCorrections[_M_co]: + """Create an `XZCorrections` instance from the XZ-corrections mappings. + + Parameters + ---------- + og : OpenGraph[_M_co] + Open graph with respect to which the corrections are defined. + x_corrections : Mapping[int, AbstractSet[int]] | None + Mapping of X-corrections: in each (`key`, `value`) pair, `key` is a measured node, and `value` is the set of nodes on which an X-correction must be applied depending on the measurement result of `key`. + z_corrections : Mapping[int, AbstractSet[int]] | None + Mapping of X-corrections: in each (`key`, `value`) pair, `key` is a measured node, and `value` is the set of nodes on which an X-correction must be applied depending on the measurement result of `key`. + + Returns + ------- + XZCorrections[_M_co] + + Notes + ----- + This method computes the partial order induced by the XZ-corrections. + """ + x_corrections = x_corrections or {} + z_corrections = z_corrections or {} + + nodes_set = set(og.graph.nodes) + outputs_set = frozenset(og.output_nodes) + non_outputs_set = nodes_set - outputs_set + + if not non_outputs_set.issuperset(x_corrections): + raise ValueError("Keys of input X-corrections contain non-measured nodes.") + if not set(z_corrections).issubset(non_outputs_set): + raise ValueError("Keys of input Z-corrections contain non-measured nodes.") + + dag = _corrections_to_dag(x_corrections, z_corrections) + partial_order_layers = _dag_to_partial_order_layers(dag) + + if partial_order_layers is None: + raise ValueError( + "Input XZ-corrections are not runnable since the induced directed graph contains closed loops." + ) + + # If there're no corrections, the partial order has 2 layers only: outputs and measured nodes. + if len(partial_order_layers) == 0: + partial_order_layers = [outputs_set] if outputs_set else [] + if non_outputs_set: + partial_order_layers.append(frozenset(non_outputs_set)) + return XZCorrections(og, x_corrections, z_corrections, tuple(partial_order_layers)) + + # If the open graph has outputs, the first element in the output of `_dag_to_partial_order_layers(dag)` may or may not contain all or some output nodes. + if outputs_set: + if measured_layer_0 := partial_order_layers[0] - outputs_set: + # `partial_order_layers[0]` contains (some or all) outputs and measured nodes + partial_order_layers = [ + outputs_set, + frozenset(measured_layer_0), + *partial_order_layers[1:], + ] + else: + # `partial_order_layers[0]` contains only (some or all) outputs + partial_order_layers = [ + outputs_set, + *partial_order_layers[1:], + ] + + ordered_nodes = frozenset.union(*partial_order_layers) + + if not ordered_nodes.issubset(nodes_set): + raise ValueError("Values of input mapping contain labels which are not nodes of the input open graph.") + + # We include all the non-output nodes not involved in the corrections in the last layer (first measured nodes). + if unordered_nodes := frozenset(nodes_set - ordered_nodes): + partial_order_layers[-1] |= unordered_nodes + + return XZCorrections(og, x_corrections, z_corrections, tuple(partial_order_layers)) + + def to_pattern( + self: XZCorrections[Measurement], + total_measurement_order: TotalOrder | None = None, + ) -> Pattern: + """Generate a unique pattern from an instance of `XZCorrections[Measurement]`. + + Parameters + ---------- + total_measurement_order : TotalOrder | None + Ordered sequence of all the non-output nodes in the open graph indicating the measurement order. This parameter must be compatible with the partial order induced by the XZ-corrections. + Optional, defaults to `None`. If `None` an arbitrary total order compatible with `self.partial_order_layers` is generated. + + Returns + ------- + Pattern + + Notes + ----- + - The `XZCorrections` instance must be of parametric type `Measurement` to allow for a pattern extraction, otherwise the underlying open graph does not contain information about the measurement angles. + + - The resulting pattern is guaranteed to be runnable if the `XZCorrections` object is well formed, but does not need to be deterministic. It will be deterministic if the XZ-corrections were inferred from a flow. In this case, this routine follows the recipe in Theorems 1, 2 and 4 in Ref. [1]. + + References + ---------- + [1] Browne et al., NJP 9, 250 (2007). + """ + if total_measurement_order is None: + total_measurement_order = self.generate_total_measurement_order() + elif not self.is_compatible(total_measurement_order): + raise ValueError( + "The input total measurement order is not compatible with the partial order induced by the XZ-corrections." + ) + + pattern = graphix.pattern.Pattern(input_nodes=self.og.input_nodes) + non_input_nodes = set(self.og.graph.nodes) - set(self.og.input_nodes) + + for i in non_input_nodes: + pattern.add(N(node=i)) + for e in self.og.graph.edges: + pattern.add(E(nodes=e)) + + for measured_node in total_measurement_order: + measurement = self.og.measurements[measured_node] + pattern.add(M(node=measured_node, plane=measurement.plane, angle=measurement.angle)) + + for corrected_node in self.z_corrections.get(measured_node, []): + pattern.add(Z(node=corrected_node, domain={measured_node})) + + for corrected_node in self.x_corrections.get(measured_node, []): + pattern.add(X(node=corrected_node, domain={measured_node})) + + pattern.reorder_output_nodes(self.og.output_nodes) + return pattern + + def generate_total_measurement_order(self) -> TotalOrder: + """Generate a sequence of all the non-output nodes in the open graph in an arbitrary order compatible with the intrinsic partial order of the XZ-corrections. + + Returns + ------- + TotalOrder + """ + shift = 1 if self.og.output_nodes else 0 + total_order = [node for layer in reversed(self.partial_order_layers[shift:]) for node in layer] + + assert set(total_order) == set(self.og.graph.nodes) - set(self.og.output_nodes) + return total_order + + def extract_dag(self) -> nx.DiGraph[int]: + """Extract the directed graph induced by the XZ-corrections. + + Returns + ------- + nx.DiGraph[int] + Directed graph in which an edge `i -> j` represents a correction applied to qubit `j`, conditioned on the measurement outcome of qubit `i`. + + Notes + ----- + - Not all nodes of the underlying open graph are nodes of the returned directed graph, but only those involved in a correction, either as corrected qubits or belonging to a correction domain. + - The output of this method is not guaranteed to be a directed acyclical graph (i.e., a directed graph without any loops). This is only the case if the `XZCorrections` object is well formed, which is verified by the method :func:`XZCorrections.is_wellformed`. + """ + return _corrections_to_dag(self.x_corrections, self.z_corrections) + + def is_compatible(self, total_measurement_order: TotalOrder) -> bool: + """Verify if a given total measurement order is compatible with the intrisic partial order of the XZ-corrections. + + Parameters + ---------- + total_measurement_order: TotalOrder + An ordered sequence of all the non-output nodes in the open graph. + + Returns + ------- + bool + `True` if `total_measurement_order` is compatible with `self.partial_order_layers`, `False` otherwise. + """ + non_outputs_set = set(self.og.graph.nodes) - set(self.og.output_nodes) + + if set(total_measurement_order) != non_outputs_set: + print("The input total measurement order does not contain all non-output nodes.") + return False + + if len(total_measurement_order) != len(non_outputs_set): + print("The input total measurement order contains duplicates.") + return False + + shift = 1 if self.og.output_nodes else 0 + measured_layers = list(reversed(self.partial_order_layers[shift:])) + + i = 0 + n_measured_layers = len(measured_layers) + layer = measured_layers[0] + + for node in total_measurement_order: + while node not in layer: + i += 1 + if i == n_measured_layers: + return False + layer = measured_layers[i] + + return True + + +@dataclass(frozen=True) +class PauliFlow(Generic[_M_co]): + """An unmutable dataclass providing a representation of a Pauli flow. + + Attributes + ---------- + og : OpenGraph[_M_co] + The open graph with respect to which the Pauli flow is defined. + correction_function : Mapping[int, AbstractSet[int] + Pauli flow correction function. `correction_function[i]` is the set of qubits correcting the measurement of qubit `i`. + partial_order_layers : Sequence[AbstractSet[int]] + Partial order between the open graph's nodes in a layer form. The set `layers[i]` comprises the nodes in layer `i`. Nodes in layer `i` are "larger" in the partial order than nodes in layer `i+1`. Output nodes are always in layer 0. + + Notes + ----- + - See Definition 5 in Ref. [1] for a definition of Pauli flow. + + - The flow's correction function defines a partial order (see Def. 2.8 and 2.9, Lemma 2.11 and Theorem 2.12 in Ref. [2]), therefore, only `og` and `correction_function` are necessary to initialize an `PauliFlow` instance (see :func:`PauliFlow.from_correction_matrix`). However, flow-finding algorithms generate a partial order in a layer form, which is necessary to extract the flow's XZ-corrections, so it is stored as an attribute. + + - A correct flow can only exist on an open graph with output nodes, so `layers[0]` always contains a finite set of nodes. + + References + ---------- + [1] Browne et al., 2007 New J. Phys. 9 250 (arXiv:quant-ph/0702212). + [2] Mitosek and Backens, 2024 (arXiv:2410.23439). + + """ + + og: OpenGraph[_M_co] + correction_function: Mapping[int, AbstractSet[int]] + partial_order_layers: Sequence[AbstractSet[int]] + + @classmethod + def try_from_correction_matrix(cls, correction_matrix: CorrectionMatrix[_M_co]) -> Self | None: + """Initialize a Pauli flow object from a matrix encoding a correction function. + + Attributes + ---------- + correction_matrix : CorrectionMatrix[_M_co] + Algebraic representation of the correction function. + + Returns + ------- + Self | None + A Pauli flow if it exists, `None` otherwise. + + Notes + ----- + This method verifies if there exists a partial measurement order on the input open graph compatible with the input correction matrix. See Lemma 3.12, and Theorem 3.1 in Ref. [1]. Failure to find a partial order implies the non-existence of a Pauli flow if the correction matrix was calculated by means of Algorithms 2 and 3 in [1]. + + References + ---------- + [1] Mitosek and Backens, 2024 (arXiv:2410.23439). + """ + correction_function = correction_matrix.to_correction_function() + partial_order_layers = compute_partial_order_layers(correction_matrix) + if partial_order_layers is None: + return None + + return cls(correction_matrix.aog.og, correction_function, partial_order_layers) + + def to_corrections(self) -> XZCorrections[_M_co]: + """Compute the X and Z corrections induced by the Pauli flow encoded in `self`. + + Returns + ------- + XZCorrections[_M_co] + + Notes + ----- + This method partially implements Theorem 4 in [1]. The generated X and Z corrections can be used to obtain a robustly deterministic pattern on the underlying open graph. + + References + ---------- + [1] Browne et al., 2007 New J. Phys. 9 250 (arXiv:quant-ph/0702212). + """ + future = copy(self.partial_order_layers[0]) # Sets are mutable + x_corrections: dict[int, AbstractSet[int]] = {} # {domain: nodes} + z_corrections: dict[int, AbstractSet[int]] = {} # {domain: nodes} + + for layer in self.partial_order_layers[1:]: + for measured_node in layer: + correcting_set = self.correction_function[measured_node] + # Conditionals avoid storing empty correction sets + if x_corrected_nodes := correcting_set & future: + x_corrections[measured_node] = x_corrected_nodes + if z_corrected_nodes := self.og.odd_neighbors(correcting_set) & future: + z_corrections[measured_node] = z_corrected_nodes + + future |= layer + + return XZCorrections(self.og, x_corrections, z_corrections, self.partial_order_layers) + + +@dataclass(frozen=True) +class GFlow(PauliFlow[_PM_co], Generic[_PM_co]): + """An unmutable subclass of `PauliFlow` providing a representation of a generalised flow (gflow). + + This class differs from its parent class in the following: + - It cannot be constructed from `OpenGraph[Axis]` instances, since the gflow is only defined for planar measurements. + - The extraction of XZ-corrections from the gflow does not require knowledge on the partial order. + - The method :func:`GFlow.is_well_formed` verifies the definition of gflow (Definition 2.36 in Ref. [1]). + + References + ---------- + [1] Backens et al., Quantum 5, 421 (2021), doi.org/10.22331/q-2021-03-25-421 + + """ + + @override + def to_corrections(self) -> XZCorrections[_PM_co]: + r"""Compute the XZ-corrections induced by the generalised flow encoded in `self`. + + Returns + ------- + XZCorrections[_PM_co] + + Notes + ----- + - This function partially implements Theorem 2 in Ref. [1]. The generated XZ-corrections can be used to obtain a robustly deterministic pattern on the underlying open graph. + + - Contrary to the overridden method in the parent class, here we do not need any information on the partial order to build the corrections since a valid correction function :math:`g` guarantees that both :math:`g(i)\setminus \{i\}` and :math:`Odd(g(i))` are in the future of :math:`i`. + + References + ---------- + [1] Browne et al., 2007 New J. Phys. 9 250 (arXiv:quant-ph/0702212). + """ + x_corrections: dict[int, AbstractSet[int]] = {} # {domain: nodes} + z_corrections: dict[int, AbstractSet[int]] = {} # {domain: nodes} + + for measured_node, correcting_set in self.correction_function.items(): + # Conditionals avoid storing empty correction sets + if x_corrected_nodes := correcting_set - {measured_node}: + x_corrections[measured_node] = x_corrected_nodes + if z_corrected_nodes := self.og.odd_neighbors(correcting_set) - {measured_node}: + z_corrections[measured_node] = z_corrected_nodes + + return XZCorrections(self.og, x_corrections, z_corrections, self.partial_order_layers) + + +@dataclass(frozen=True) +class CausalFlow(GFlow[_PM_co], Generic[_PM_co]): + """An unmutable subclass of `GFlow` providing a representation of a causal flow. + + This class differs from its parent class in the following: + - The extraction of XZ-corrections from the causal flow does assumes that correction sets have one element only. + - The method :func:`CausalFlow.is_well_formed` verifies the definition of causal flow (Definition 2 in Ref. [1]). + + References + ---------- + [1] Browne et al., 2007 New J. Phys. 9 250 (arXiv:quant-ph/0702212). + + """ + + @override + @classmethod + def try_from_correction_matrix(cls, correction_matrix: CorrectionMatrix[_PM_co]) -> None: + raise NotImplementedError("Initialization of a causal flow from a correction matrix is not supported.") + + @override + def to_corrections(self) -> XZCorrections[_PM_co]: + r"""Compute the XZ-corrections induced by the causal flow encoded in `self`. + + Returns + ------- + XZCorrections[_PM_co] + + Notes + ----- + This function partially implements Theorem 1 in Ref. [1]. The generated XZ-corrections can be used to obtain a robustly deterministic pattern on the underlying open graph. + + References + ---------- + [1] Browne et al., 2007 New J. Phys. 9 250 (arXiv:quant-ph/0702212). + """ + x_corrections: dict[int, AbstractSet[int]] = {} # {domain: nodes} + z_corrections: dict[int, AbstractSet[int]] = {} # {domain: nodes} + + for measured_node, correcting_set in self.correction_function.items(): + # Conditionals avoid storing empty correction sets + if x_corrected_nodes := correcting_set: + x_corrections[measured_node] = x_corrected_nodes + if z_corrected_nodes := self.og.neighbors(correcting_set) - {measured_node}: + z_corrections[measured_node] = z_corrected_nodes + + return XZCorrections(self.og, x_corrections, z_corrections, self.partial_order_layers) + + +def _corrections_to_dag( + x_corrections: Mapping[int, AbstractSet[int]], z_corrections: Mapping[int, AbstractSet[int]] +) -> nx.DiGraph[int]: + """Convert an XZ-corrections mapping into a directed graph. + + Parameters + ---------- + x_corrections : Mapping[int, AbstractSet[int]] + Mapping of X-corrections: in each (`key`, `value`) pair, `key` is a measured node, and `value` is the set of nodes on which an X-correction must be applied depending on the measurement result of `key`. + z_corrections : Mapping[int, AbstractSet[int]] + Mapping of Z-corrections: in each (`key`, `value`) pair, `key` is a measured node, and `value` is the set of nodes on which an Z-correction must be applied depending on the measurement result of `key`. + + Returns + ------- + nx.DiGraph[int] + Directed graph in which an edge `i -> j` represents a correction applied to qubit `j`, conditioned on the measurement outcome of qubit `i`. + + Notes + ----- + See :func:`XZCorrections.extract_dag`. + """ + relations = ( + (measured_node, corrected_node) + for corrections in (x_corrections, z_corrections) + for measured_node, corrected_nodes in corrections.items() + for corrected_node in corrected_nodes + ) + + return nx.DiGraph(relations) + + +def _dag_to_partial_order_layers(dag: nx.DiGraph[int]) -> list[frozenset[int]] | None: + """Return the partial order encoded in a directed graph in a layer form if it exists. + + Parameters + ---------- + dag : nx.DiGraph[int] + A directed graph. + + Returns + ------- + list[set[int]] | None + Partial order between corrected qubits in a layer form or `None` if the input directed graph is not acyclical. + The set `layers[i]` comprises the nodes in layer `i`. Nodes in layer `i` are "larger" in the partial order than nodes in layer `i+1`. + """ + try: + topo_gen = reversed(list(nx.topological_generations(dag))) + except nx.NetworkXUnfeasible: + return None + + return [frozenset(layer) for layer in topo_gen] diff --git a/graphix/fundamentals.py b/graphix/fundamentals.py index 3f5de010..2b0733d1 100644 --- a/graphix/fundamentals.py +++ b/graphix/fundamentals.py @@ -5,11 +5,15 @@ import enum import sys import typing -from enum import Enum +from abc import ABC, ABCMeta, abstractmethod +from enum import Enum, EnumMeta from typing import TYPE_CHECKING, SupportsComplex, SupportsFloat, SupportsIndex, overload import typing_extensions +# override introduced in Python 3.12 +from typing_extensions import override + from graphix.ops import Ops from graphix.parameter import cos_sin from graphix.repr_mixins import EnumReprMixin @@ -214,7 +218,50 @@ def matrix(self) -> npt.NDArray[np.complex128]: typing_extensions.assert_never(self) -class Axis(EnumReprMixin, Enum): +class CustomMeta(ABCMeta, EnumMeta): + """Custom metaclass to allow multiple inheritance from `Enum` and `ABC`.""" + + +class AbstractMeasurement(ABC): + """Abstract base class for measurement objects. + + Measurement objects are: + - :class:`graphix.measurements.Measurement`. + - :class:`graphix.fundamentals.Plane`. + - :class:`graphix.fundamentals.Axis`. + + """ + + @abstractmethod + def to_plane_or_axis(self) -> Plane | Axis: + """Return the plane or axis of a measurement object. + + Returns + ------- + Plane | Axis + """ + + +class AbstractPlanarMeasurement(AbstractMeasurement): + """Abstract base class for planar measurement objects. + + Planar measurement objects are: + - :class:`graphix.measurements.Measurement`. + - :class:`graphix.fundamentals.Plane`. + + """ + + @abstractmethod + def to_plane(self) -> Plane: + """Return the plane of a measurement object. + + Returns + ------- + Plane + """ + + +class Axis(AbstractMeasurement, EnumReprMixin, Enum, metaclass=CustomMeta): """Axis: *X*, *Y* or *Z*.""" X = enum.auto() @@ -232,8 +279,12 @@ def matrix(self) -> npt.NDArray[np.complex128]: return Ops.Z typing_extensions.assert_never(self) + @override + def to_plane_or_axis(self) -> Axis: + return self -class Plane(EnumReprMixin, Enum): + +class Plane(AbstractPlanarMeasurement, EnumReprMixin, Enum, metaclass=CustomMeta): # TODO: Refactor using match """Plane: *XY*, *YZ* or *XZ*.""" @@ -317,3 +368,23 @@ def from_axes(a: Axis, b: Axis) -> Plane: return Plane.XZ assert a == b raise ValueError(f"Cannot make a plane giving the same axis {a} twice.") + + @override + def to_plane_or_axis(self) -> Plane: + """Return the plane. + + Returns + ------- + Plane + """ + return self + + @override + def to_plane(self) -> Plane: + """Return the plane. + + Returns + ------- + Plane + """ + return self diff --git a/graphix/generator.py b/graphix/generator.py deleted file mode 100644 index 39cf4b50..00000000 --- a/graphix/generator.py +++ /dev/null @@ -1,189 +0,0 @@ -"""MBQC pattern generator.""" - -from __future__ import annotations - -from typing import TYPE_CHECKING - -import graphix.gflow -from graphix.command import E, M, N, X, Z -from graphix.fundamentals import Plane -from graphix.pattern import Pattern - -if TYPE_CHECKING: - from collections.abc import Iterable, Mapping - from collections.abc import Set as AbstractSet - - import networkx as nx - - from graphix.parameter import ExpressionOrFloat - - -def generate_from_graph( - graph: nx.Graph[int], - angles: Mapping[int, ExpressionOrFloat], - inputs: Iterable[int], - outputs: Iterable[int], - meas_planes: Mapping[int, Plane] | None = None, -) -> Pattern: - r"""Generate the measurement pattern from open graph and measurement angles. - - This function takes an open graph ``G = (nodes, edges, input, outputs)``, - specified by :class:`networkx.Graph` and two lists specifying input and output nodes. - Currently we support XY-plane measurements. - - Searches for the flow in the open graph using :func:`graphix.gflow.find_flow` and if found, - construct the measurement pattern according to the theorem 1 of [NJP 9, 250 (2007)]. - - Then, if no flow was found, searches for gflow using :func:`graphix.gflow.find_gflow`, - from which measurement pattern can be constructed from theorem 2 of [NJP 9, 250 (2007)]. - - Then, if no gflow was found, searches for Pauli flow using :func:`graphix.gflow.find_pauliflow`, - from which measurement pattern can be constructed from theorem 4 of [NJP 9, 250 (2007)]. - - The constructed measurement pattern deterministically realize the unitary embedding - - .. math:: - - U = \left( \prod_i \langle +_{\alpha_i} |_i \right) E_G N_{I^C}, - - where the measurements (bras) with always :math:`\langle+|` bases determined by the measurement - angles :math:`\alpha_i` are applied to the measuring nodes, - i.e. the randomness of the measurement is eliminated by the added byproduct commands. - - .. seealso:: :func:`graphix.gflow.find_flow` :func:`graphix.gflow.find_gflow` :func:`graphix.gflow.find_pauliflow` :class:`graphix.pattern.Pattern` - - Parameters - ---------- - graph : :class:`networkx.Graph` - Graph on which MBQC should be performed - angles : dict - measurement angles for each nodes on the graph (unit of pi), except output nodes - inputs : list - list of node indices for input nodes - outputs : list - list of node indices for output nodes - meas_planes : dict - optional: measurement planes for each nodes on the graph, except output nodes - - Returns - ------- - pattern : graphix.pattern.Pattern - constructed pattern. - """ - inputs_set = set(inputs) - outputs_set = set(outputs) - - measuring_nodes = set(graph.nodes) - outputs_set - - meas_planes = dict.fromkeys(measuring_nodes, Plane.XY) if not meas_planes else dict(meas_planes) - - # search for flow first - f, l_k = graphix.gflow.find_flow(graph, inputs_set, outputs_set, meas_planes=meas_planes) - if f is not None and l_k is not None: - # flow found - pattern = _flow2pattern(graph, angles, inputs, f, l_k) - pattern.reorder_output_nodes(outputs) - return pattern - - # no flow found - we try gflow - g, l_k = graphix.gflow.find_gflow(graph, inputs_set, outputs_set, meas_planes=meas_planes) - if g is not None and l_k is not None: - # gflow found - pattern = _gflow2pattern(graph, angles, inputs, meas_planes, g, l_k) - pattern.reorder_output_nodes(outputs) - return pattern - - # no flow or gflow found - we try pflow - p, l_k = graphix.gflow.find_pauliflow(graph, inputs_set, outputs_set, meas_planes=meas_planes, meas_angles=angles) - if p is not None and l_k is not None: - # pflow found - pattern = _pflow2pattern(graph, angles, inputs, meas_planes, p, l_k) - pattern.reorder_output_nodes(outputs) - return pattern - - raise ValueError("no flow or gflow or pflow found") - - -def _flow2pattern( - graph: nx.Graph[int], - angles: Mapping[int, ExpressionOrFloat], - inputs: Iterable[int], - f: Mapping[int, AbstractSet[int]], - l_k: Mapping[int, int], -) -> Pattern: - """Construct a measurement pattern from a causal flow according to the theorem 1 of [NJP 9, 250 (2007)].""" - depth, layers = graphix.gflow.get_layers(l_k) - pattern = Pattern(input_nodes=inputs) - for i in set(graph.nodes) - set(inputs): - pattern.add(N(node=i)) - for e in graph.edges: - pattern.add(E(nodes=e)) - measured: list[int] = [] - for i in range(depth, 0, -1): # i from depth, depth-1, ... 1 - for j in layers[i]: - measured.append(j) - pattern.add(M(node=j, angle=angles[j])) - neighbors: set[int] = set() - for k in f[j]: - neighbors |= set(graph.neighbors(k)) - for k in neighbors - {j}: - # if k not in measured: - pattern.add(Z(node=k, domain={j})) - (fj,) = f[j] - pattern.add(X(node=fj, domain={j})) - return pattern - - -def _gflow2pattern( - graph: nx.Graph[int], - angles: Mapping[int, ExpressionOrFloat], - inputs: Iterable[int], - meas_planes: Mapping[int, Plane], - g: Mapping[int, AbstractSet[int]], - l_k: Mapping[int, int], -) -> Pattern: - """Construct a measurement pattern from a generalized flow according to the theorem 2 of [NJP 9, 250 (2007)].""" - depth, layers = graphix.gflow.get_layers(l_k) - pattern = Pattern(input_nodes=inputs) - for i in set(graph.nodes) - set(inputs): - pattern.add(N(node=i)) - for e in graph.edges: - pattern.add(E(nodes=e)) - for i in range(depth, 0, -1): # i from depth, depth-1, ... 1 - for j in layers[i]: - pattern.add(M(node=j, plane=meas_planes[j], angle=angles[j])) - odd_neighbors = graphix.gflow.find_odd_neighbor(graph, g[j]) - for k in odd_neighbors - {j}: - pattern.add(Z(node=k, domain={j})) - for k in g[j] - {j}: - pattern.add(X(node=k, domain={j})) - return pattern - - -def _pflow2pattern( - graph: nx.Graph[int], - angles: Mapping[int, ExpressionOrFloat], - inputs: Iterable[int], - meas_planes: Mapping[int, Plane], - p: Mapping[int, AbstractSet[int]], - l_k: Mapping[int, int], -) -> Pattern: - """Construct a measurement pattern from a Pauli flow according to the theorem 4 of [NJP 9, 250 (2007)].""" - depth, layers = graphix.gflow.get_layers(l_k) - pattern = Pattern(input_nodes=inputs) - for i in set(graph.nodes) - set(inputs): - pattern.add(N(node=i)) - for e in graph.edges: - pattern.add(E(nodes=e)) - for i in range(depth, 0, -1): # i from depth, depth-1, ... 1 - for j in layers[i]: - pattern.add(M(node=j, plane=meas_planes[j], angle=angles[j])) - odd_neighbors = graphix.gflow.find_odd_neighbor(graph, p[j]) - future_nodes: set[int] = set.union( - *(nodes for (layer, nodes) in layers.items() if layer < i) - ) # {k | k > j}, with "j" last corrected node and ">" the Pauli flow ordering - for k in odd_neighbors & future_nodes: - pattern.add(Z(node=k, domain={j})) - for k in p[j] & future_nodes: - pattern.add(X(node=k, domain={j})) - return pattern diff --git a/graphix/gflow.py b/graphix/gflow.py index 439c6367..62a7c849 100644 --- a/graphix/gflow.py +++ b/graphix/gflow.py @@ -16,9 +16,9 @@ from typing_extensions import assert_never +import graphix.find_pflow import graphix.opengraph from graphix.command import CommandKind -from graphix.find_pflow import find_pflow as _find_pflow from graphix.fundamentals import Axis, Plane from graphix.measurements import Measurement, PauliMeasurement from graphix.parameter import Placeholder @@ -86,12 +86,12 @@ def find_gflow( """ meas = {node: Measurement(Placeholder("Angle"), plane) for node, plane in meas_planes.items()} og = graphix.opengraph.OpenGraph( - inside=graph, - inputs=list(iset), - outputs=list(oset), + graph=graph, + input_nodes=list(iset), + output_nodes=list(oset), measurements=meas, ) - gf = _find_pflow(og) + gf = graphix.find_pflow.find_pflow(og) if gf is None: return None, None # This is to comply with old API. It will be change in the future to `None`` return gf[0], gf[1] @@ -271,12 +271,12 @@ def find_pauliflow( """ meas = {node: Measurement(angle, meas_planes[node]) for node, angle in meas_angles.items()} og = graphix.opengraph.OpenGraph( - inside=graph, - inputs=list(iset), - outputs=list(oset), + graph=graph, + input_nodes=list(iset), + output_nodes=list(oset), measurements=meas, ) - pf = _find_pflow(og) + pf = graphix.find_pflow.find_pflow(og) if pf is None: return None, None # This is to comply with old API. It will be change in the future to `None`` return pf[0], pf[1] diff --git a/graphix/measurements.py b/graphix/measurements.py index fa382705..179c26f4 100644 --- a/graphix/measurements.py +++ b/graphix/measurements.py @@ -2,14 +2,14 @@ from __future__ import annotations -import dataclasses import math +from dataclasses import dataclass from typing import Literal, NamedTuple, SupportsInt from typing_extensions import TypeAlias # TypeAlias introduced in Python 3.10 from graphix import utils -from graphix.fundamentals import Axis, Plane, Sign +from graphix.fundamentals import AbstractPlanarMeasurement, Axis, Plane, Sign # Ruff suggests to move this import to a type-checking block, but dataclass requires it here from graphix.parameter import ExpressionOrFloat # noqa: TC001 @@ -27,7 +27,7 @@ def toggle_outcome(outcome: Outcome) -> Outcome: return 1 if outcome == 0 else 0 -@dataclasses.dataclass +@dataclass class Domains: """Represent `X^sZ^t` where s and t are XOR of results from given sets of indices.""" @@ -35,11 +35,16 @@ class Domains: t_domain: set[int] -class Measurement(NamedTuple): - """An MBQC measurement. +@dataclass +class Measurement(AbstractPlanarMeasurement): + r"""An MBQC measurement. - :param angle: the angle of the measurement. Should be between [0, 2) - :param plane: the measurement plane + Attributes + ---------- + angle : Expressionor Float + The angle of the measurement in units of :math:`\pi`. Should be between [0, 2). + plane : graphix.fundamentals.Plane + The measurement plane. """ angle: ExpressionOrFloat @@ -50,7 +55,7 @@ def isclose(self, other: Measurement, rel_tol: float = 1e-09, abs_tol: float = 0 Example ------- - >>> from graphix.opengraph import Measurement + >>> from graphix.measurements import Measurement >>> from graphix.fundamentals import Plane >>> Measurement(0.0, Plane.XY).isclose(Measurement(0.0, Plane.XY)) True @@ -65,6 +70,30 @@ def isclose(self, other: Measurement, rel_tol: float = 1e-09, abs_tol: float = 0 else self.angle == other.angle ) and self.plane == other.plane + def to_plane_or_axis(self) -> Plane | Axis: + """Return the measurements's plane or axis. + + Returns + ------- + Plane | Axis + + Notes + ----- + Measurements with Pauli angles (i.e., ``self.angle == n/2`` with ``n`` an integer) are interpreted as `Axis` instances. + """ + if pm := PauliMeasurement.try_from(self.plane, self.angle): + return pm.axis + return self.plane + + def to_plane(self) -> Plane: + """Return the measurement's plane. + + Returns + ------- + Plane + """ + return self.plane + class PauliMeasurement(NamedTuple): """Pauli measurement.""" diff --git a/graphix/opengraph.py b/graphix/opengraph.py index e6bbadd4..540aab3d 100644 --- a/graphix/opengraph.py +++ b/graphix/opengraph.py @@ -1,68 +1,92 @@ -"""Provides a class for open graphs.""" +"""Class for open graph states.""" from __future__ import annotations from dataclasses import dataclass -from typing import TYPE_CHECKING +from typing import TYPE_CHECKING, Generic, TypeVar import networkx as nx -import graphix.generator -from graphix.measurements import Measurement +from graphix.flow._find_cflow import find_cflow +from graphix.flow._find_gpflow import AlgebraicOpenGraph, PlanarAlgebraicOpenGraph, compute_correction_matrix +from graphix.flow.core import GFlow, PauliFlow +from graphix.fundamentals import AbstractMeasurement, AbstractPlanarMeasurement if TYPE_CHECKING: - from collections.abc import Iterable, Mapping + from collections.abc import Collection, Iterable, Mapping, Sequence + from graphix.flow.core import CausalFlow + from graphix.measurements import Measurement from graphix.pattern import Pattern +# TODO: Maybe move these definitions to graphix.fundamentals and graphix.measurements ? Now they are redefined in graphix.flow._find_gpflow, not very elegant. +_M_co = TypeVar("_M_co", bound=AbstractMeasurement, covariant=True) +_PM_co = TypeVar("_PM_co", bound=AbstractPlanarMeasurement, covariant=True) -@dataclass(frozen=True) -class OpenGraph: - """Open graph contains the graph, measurement, and input and output nodes. - - This is the graph we wish to implement deterministically. - :param inside: the underlying :class:`networkx.Graph` state - :param measurements: a dictionary whose key is the ID of a node and the - value is the measurement at that node - :param inputs: an ordered list of node IDs that are inputs to the graph - :param outputs: an ordered list of node IDs that are outputs of the graph +@dataclass(frozen=True) +class OpenGraph(Generic[_M_co]): + """An unmutable dataclass providing a representation of open graph states. + + Attributes + ---------- + graph : networkx.Graph[int] + The underlying resource-state graph. Nodes represent qubits and edges represent the application of :math:`CZ` gate on the linked nodes. + input_nodes : Sequence[int] + An ordered sequence of node labels corresponding to the open graph inputs. + output_nodes : Sequence[int] + An ordered sequence of node labels corresponding to the open graph outputs. + measurements : Mapping[int, _M_co] + A mapping between the non-output nodes of the open graph (``key``) and their corresponding measurement label (``value``). Measurement labels can be specified as `Measurement`, `Plane` or `Axis` instances. + + Notes + ----- + The inputs and outputs of `OpenGraph` instances in Graphix are defined as ordered sequences of node labels. This contrasts the usual definition of open graphs in the literature, where inputs and outputs are unordered sets of nodes labels. This restriction facilitates the interplay with `Pattern` objects, where the order of input and output nodes represents a choice of Hilbert space basis. Example ------- >>> import networkx as nx >>> from graphix.fundamentals import Plane - >>> from graphix.opengraph import OpenGraph, Measurement - >>> - >>> inside_graph = nx.Graph([(0, 1), (1, 2), (2, 0)]) + >>> from graphix.opengraph import OpenGraph + >>> from graphix.measurements import Measurement >>> + >>> graph = nx.Graph([(0, 1), (1, 2)]) >>> measurements = {i: Measurement(0.5 * i, Plane.XY) for i in range(2)} - >>> inputs = [0] - >>> outputs = [2] - >>> og = OpenGraph(inside_graph, measurements, inputs, outputs) + >>> input_nodes = [0] + >>> output_nodes = [2] + >>> og = OpenGraph(graph, input_nodes, output_nodes, measurements) """ - inside: nx.Graph[int] - measurements: dict[int, Measurement] - inputs: list[int] # Inputs are ordered - outputs: list[int] # Outputs are ordered + graph: nx.Graph[int] + input_nodes: Sequence[int] + output_nodes: Sequence[int] + measurements: Mapping[int, _M_co] def __post_init__(self) -> None: - """Validate the open graph.""" - if not all(node in self.inside.nodes for node in self.measurements): + """Validate the correctness of the open graph.""" + all_nodes = set(self.graph.nodes) + inputs = set(self.input_nodes) + outputs = set(self.output_nodes) + + if not set(self.measurements).issubset(all_nodes): raise ValueError("All measured nodes must be part of the graph's nodes.") - if not all(node in self.inside.nodes for node in self.inputs): + if not inputs.issubset(all_nodes): raise ValueError("All input nodes must be part of the graph's nodes.") - if not all(node in self.inside.nodes for node in self.outputs): + if not outputs.issubset(all_nodes): raise ValueError("All output nodes must be part of the graph's nodes.") - if any(node in self.outputs for node in self.measurements): - raise ValueError("Output node cannot be measured.") - if len(set(self.inputs)) != len(self.inputs): + if outputs & self.measurements.keys(): + raise ValueError("Output nodes cannot be measured.") + if all_nodes - outputs != self.measurements.keys(): + raise ValueError("All non-output nodes must be measured.") + if len(inputs) != len(self.input_nodes): raise ValueError("Input nodes contain duplicates.") - if len(set(self.outputs)) != len(self.outputs): + if len(outputs) != len(self.output_nodes): raise ValueError("Output nodes contain duplicates.") - def isclose(self, other: OpenGraph, rel_tol: float = 1e-09, abs_tol: float = 0.0) -> bool: + # TODO: Up docstrings and generalise to any type + def isclose( + self: OpenGraph[Measurement], other: OpenGraph[Measurement], rel_tol: float = 1e-09, abs_tol: float = 0.0 + ) -> bool: """Return `True` if two open graphs implement approximately the same unitary operator. Ensures the structure of the graphs are the same and all @@ -71,10 +95,10 @@ def isclose(self, other: OpenGraph, rel_tol: float = 1e-09, abs_tol: float = 0.0 This doesn't check they are equal up to an isomorphism. """ - if not nx.utils.graphs_equal(self.inside, other.inside): + if not nx.utils.graphs_equal(self.graph, other.graph): return False - if self.inputs != other.inputs or self.outputs != other.outputs: + if self.input_nodes != other.input_nodes or self.output_nodes != other.output_nodes: return False if set(self.measurements.keys()) != set(other.measurements.keys()): @@ -85,38 +109,230 @@ def isclose(self, other: OpenGraph, rel_tol: float = 1e-09, abs_tol: float = 0.0 for node, m in self.measurements.items() ) - @staticmethod - def from_pattern(pattern: Pattern) -> OpenGraph: - """Initialise an `OpenGraph` object based on the resource-state graph associated with the measurement pattern.""" - graph = pattern.extract_graph() + def to_pattern(self: OpenGraph[Measurement]) -> Pattern: + """Extract a deterministic pattern from an `OpenGraph[Measurement]` if it exists. + + Returns + ------- + Pattern + A deterministic pattern on the open graph. + + Raises + ------ + OpenGraphError + If the open graph does not have flow. + + Notes + ----- + - The open graph instance must be of parametric type `Measurement` to allow for a pattern extraction, otherwise it does not contain information about the measurement angles. + + - This method proceeds by searching a flow on the open graph and converting it into a pattern as prescripted in Ref. [1]. + It first attempts to find a causal flow because the corresponding flow-finding algorithm has lower complexity. If it fails, it attemps to find a Pauli flow because this property is more general than a generalised flow, and the corresponding flow-finding algorithms have the same complexity in the current implementation. + + References + ---------- + [1] Browne et al., NJP 9, 250 (2007) + """ + for extractor in (self.find_causal_flow, self.find_pauli_flow): + flow = extractor() + if flow is not None: + return flow.to_corrections().to_pattern() + + raise OpenGraphError("The open graph does not have flow. It does not support a deterministic pattern.") + + def neighbors(self, nodes: Collection[int]) -> set[int]: + """Return the set containing the neighborhood of a set of nodes in the open graph. + + Parameters + ---------- + nodes : Collection[int] + Set of nodes whose neighborhood is to be found + + Returns + ------- + neighbors_set : set[int] + Neighborhood of set `nodes`. + """ + neighbors_set: set[int] = set() + for node in nodes: + neighbors_set |= set(self.graph.neighbors(node)) + return neighbors_set - inputs = pattern.input_nodes - outputs = pattern.output_nodes + def odd_neighbors(self, nodes: Collection[int]) -> set[int]: + """Return the set containing the odd neighborhood of a set of nodes in the open graph. - meas_planes = pattern.get_meas_plane() - meas_angles = pattern.get_angles() - meas = {node: Measurement(meas_angles[node], meas_planes[node]) for node in meas_angles} + Parameters + ---------- + nodes : Collection[int] + Set of nodes whose odd neighborhood is to be found - return OpenGraph(graph, meas, inputs, outputs) + Returns + ------- + odd_neighbors_set : set[int] + Odd neighborhood of set `nodes`. + """ + odd_neighbors_set: set[int] = set() + for node in nodes: + odd_neighbors_set ^= self.neighbors([node]) + return odd_neighbors_set - def to_pattern(self) -> Pattern: - """Convert the `OpenGraph` into a `Pattern`. + def extract_causal_flow(self: OpenGraph[_PM_co]) -> CausalFlow[_PM_co]: + """Try to extract a causal flow on the open graph. - Will raise an exception if the open graph does not have flow, gflow, or - Pauli flow. - The pattern will be generated using maximally-delayed flow. + This method is a wrapper over :func:`OpenGraph.find_causal_flow` with a single return type. + + Returns + ------- + CausalFlow[_PM_co] + A causal flow object if the open graph has causal flow. + + Raises + ------ + OpenGraphError + If the open graph does not have a causal flow. + + See Also + -------- + :func:`OpenGraph.find_causal_flow` """ - g = self.inside.copy() - inputs = self.inputs - outputs = self.outputs - meas = self.measurements + cf = self.find_causal_flow() + if cf is None: + raise OpenGraphError("The open graph does not have a causal flow.") + return cf + + def extract_gflow(self: OpenGraph[_PM_co]) -> GFlow[_PM_co]: + r"""Try to extract a maximally delayed generalised flow (gflow) on the open graph. + + This method is a wrapper over :func:`OpenGraph.find_gflow` with a single return type. - angles = {node: m.angle for node, m in meas.items()} - planes = {node: m.plane for node, m in meas.items()} + Returns + ------- + GFlow[_PM_co] + A gflow object if the open graph has gflow. + + Raises + ------ + OpenGraphError + If the open graph does not have a gflow. + + See Also + -------- + :func:`OpenGraph.find_gflow` + """ + gf = self.find_gflow() + if gf is None: + raise OpenGraphError("The open graph does not have a gflow.") + return gf - return graphix.generator.generate_from_graph(g, angles, inputs, outputs, planes) + def extract_pauli_flow(self: OpenGraph[_M_co]) -> PauliFlow[_M_co]: + r"""Try to extract a maximally delayed Pauli on the open graph. - def compose(self, other: OpenGraph, mapping: Mapping[int, int]) -> tuple[OpenGraph, dict[int, int]]: + This method is a wrapper over :func:`OpenGraph.find_pauli_flow` with a single return type. + + Returns + ------- + PauliFlow[_M_co] + A Pauli flow object if the open graph has Pauli flow. + + Raises + ------ + OpenGraphError + If the open graph does not have a Pauli flow. + + See Also + -------- + :func:`OpenGraph.find_pauli_flow` + """ + pf = self.find_pauli_flow() + if pf is None: + raise OpenGraphError("The open graph does not have a Pauli flow.") + return pf + + def find_causal_flow(self: OpenGraph[_PM_co]) -> CausalFlow[_PM_co] | None: + """Return a causal flow on the open graph if it exists. + + Returns + ------- + CausalFlow[_PM_co] | None + A causal flow object if the open graph has causal flow or ``None`` otherwise. + + See Also + -------- + :func:`OpenGraph.extract_causal_flow` + + Notes + ----- + - The open graph instance must be of parametric type `Measurement` or `Plane` since the causal flow is only defined on open graphs with :math:`XY` measurements. + - This function implements the algorithm presented in Ref. [1] with polynomial complexity on the number of nodes, :math:`O(N^2)`. + + References + ---------- + [1] Mhalla and Perdrix, (2008), Finding Optimal Flows Efficiently, doi.org/10.1007/978-3-540-70575-8_70 + """ + return find_cflow(self) + + def find_gflow(self: OpenGraph[_PM_co]) -> GFlow[_PM_co] | None: + r"""Return a maximally delayed Pauli on the open graph if it exists. + + Returns + ------- + GFlow[_PM_co] | None + A gflow object if the open graph has gflow or ``None`` otherwise. + + See Also + -------- + :func:`OpenGraph.extract_gflow` + + Notes + ----- + - The open graph instance must be of parametric type `Measurement` or `Plane` since the gflow is only defined on open graphs with planar measurements. Measurement instances with a Pauli angle (integer multiple of :math:`\pi/2`) are interpreted as `Plane` instances, in contrast with :func:`OpenGraph.find_pauli_flow`. + - This function implements the algorithm presented in Ref. [1] with polynomial complexity on the number of nodes, :math:`O(N^3)`. + + References + ---------- + [1] Mitosek and Backens, 2024 (arXiv:2410.23439). + """ + aog = PlanarAlgebraicOpenGraph(self) + correction_matrix = compute_correction_matrix(aog) + if correction_matrix is None: + return None + return GFlow.try_from_correction_matrix( + correction_matrix + ) # The constructor returns `None` if the correction matrix is not compatible with any partial order on the open graph. + + def find_pauli_flow(self: OpenGraph[_M_co]) -> PauliFlow[_M_co] | None: + r"""Return a maximally delayed Pauli on the open graph if it exists. + + Returns + ------- + PauliFlow[_M_co] | None + A Pauli flow object if the open graph has Pauli flow or ``None`` otherwise. + + See Also + -------- + :func:`OpenGraph.extract_pauli_flow` + + Notes + ----- + - Measurement instances with a Pauli angle (integer multiple of :math:`\pi/2`) are interpreted as `Axis` instances, in contrast with :func:`OpenGraph.find_gflow`. + - This function implements the algorithm presented in Ref. [1] with polynomial complexity on the number of nodes, :math:`O(N^3)`. + + References + ---------- + [1] Mitosek and Backens, 2024 (arXiv:2410.23439). + """ + aog = AlgebraicOpenGraph(self) + correction_matrix = compute_correction_matrix(aog) + if correction_matrix is None: + return None + return PauliFlow.try_from_correction_matrix( + correction_matrix + ) # The constructor returns `None` if the correction matrix is not compatible with any partial order on the open graph. + + # TODO: Generalise `compose` to any type of OpenGraph + def compose( + self: OpenGraph[Measurement], other: OpenGraph[Measurement], mapping: Mapping[int, int] + ) -> tuple[OpenGraph[Measurement], dict[int, int]]: r"""Compose two open graphs by merging subsets of nodes from `self` and `other`, and relabeling the nodes of `other` that were not merged. Parameters @@ -148,7 +364,7 @@ def compose(self, other: OpenGraph, mapping: Mapping[int, int]) -> tuple[OpenGra - If only one node of the pair `{v:u}` is measured, this measure is assigned to :math:`u \in V` in the resulting open graph. - Input (and, respectively, output) nodes in the returned open graph have the order of the open graph `self` followed by those of the open graph `other`. Merged nodes are removed, except when they are input (or output) nodes in both open graphs, in which case, they appear in the order they originally had in the graph `self`. """ - if not (mapping.keys() <= other.inside.nodes): + if not (mapping.keys() <= other.graph.nodes): raise ValueError("Keys of mapping must be correspond to nodes of other.") if len(mapping) != len(set(mapping.values())): raise ValueError("Values in mapping contain duplicates.") @@ -160,18 +376,18 @@ def compose(self, other: OpenGraph, mapping: Mapping[int, int]) -> tuple[OpenGra ): raise ValueError(f"Attempted to merge nodes {v}:{u} but have different measurements") - shift = max(*self.inside.nodes, *mapping.values()) + 1 + shift = max(*self.graph.nodes, *mapping.values()) + 1 mapping_sequential = { - node: i for i, node in enumerate(sorted(other.inside.nodes - mapping.keys()), start=shift) + node: i for i, node in enumerate(sorted(other.graph.nodes - mapping.keys()), start=shift) } # assigns new labels to nodes in other not specified in mapping mapping_complete = {**mapping, **mapping_sequential} - g2_shifted = nx.relabel_nodes(other.inside, mapping_complete) - g = nx.compose(self.inside, g2_shifted) + g2_shifted = nx.relabel_nodes(other.graph, mapping_complete) + g = nx.compose(self.graph, g2_shifted) - merged = set(mapping_complete.values()) & self.inside.nodes + merged = set(mapping_complete.values()) & self.graph.nodes def merge_ports(p1: Iterable[int], p2: Iterable[int]) -> list[int]: p2_mapped = [mapping_complete[node] for node in p2] @@ -180,10 +396,14 @@ def merge_ports(p1: Iterable[int], p2: Iterable[int]) -> list[int]: part2 = [node for node in p2_mapped if node not in merged] return part1 + part2 - inputs = merge_ports(self.inputs, other.inputs) - outputs = merge_ports(self.outputs, other.outputs) + inputs = merge_ports(self.input_nodes, other.input_nodes) + outputs = merge_ports(self.output_nodes, other.output_nodes) measurements_shifted = {mapping_complete[i]: meas for i, meas in other.measurements.items()} measurements = {**self.measurements, **measurements_shifted} - return OpenGraph(g, measurements, inputs, outputs), mapping_complete + return OpenGraph(g, inputs, outputs, measurements), mapping_complete + + +class OpenGraphError(Exception): + """Exception subclass to handle open graphs errors.""" diff --git a/graphix/pattern.py b/graphix/pattern.py index d136989f..3adf02b6 100644 --- a/graphix/pattern.py +++ b/graphix/pattern.py @@ -23,7 +23,8 @@ from graphix.fundamentals import Axis, Plane, Sign from graphix.gflow import find_flow, find_gflow, get_layers from graphix.graphsim import GraphState -from graphix.measurements import Outcome, PauliMeasurement, toggle_outcome +from graphix.measurements import Measurement, Outcome, PauliMeasurement, toggle_outcome +from graphix.opengraph import OpenGraph from graphix.pretty_print import OutputFormat, pattern_to_str from graphix.simulator import PatternSimulator from graphix.states import BasicStates @@ -1147,6 +1148,37 @@ def extract_isolated_nodes(self) -> set[int]: graph = self.extract_graph() return {node for node, d in graph.degree if d == 0} + def extract_opengraph(self) -> OpenGraph[Measurement]: + """Extract the underlying resource-state open graph from the pattern. + + Returns + ------- + OpenGraph[Measurement] + + Notes + ----- + This operation loses all the information on the Clifford commands. + """ + nodes = set(self.input_nodes) + edges: set[tuple[int, int]] = set() + measurements: dict[int, Measurement] = {} + + for cmd in self.__seq: + if cmd.kind == CommandKind.N: + nodes.add(cmd.node) + elif cmd.kind == CommandKind.E: + u, v = cmd.nodes + if u > v: + u, v = v, u + edges.symmetric_difference_update({(u, v)}) + elif cmd.kind == CommandKind.M: + measurements[cmd.node] = Measurement(cmd.angle, cmd.plane) + + graph = nx.Graph(edges) + graph.add_nodes_from(nodes) + + return OpenGraph(graph, self.input_nodes, self.output_nodes, measurements) + def get_vops(self, conj: bool = False, include_identity: bool = False) -> dict[int, Clifford]: """Get local-Clifford decorations from measurement or Clifford commands. diff --git a/graphix/pyzx.py b/graphix/pyzx.py index e3d35548..33cf68da 100644 --- a/graphix/pyzx.py +++ b/graphix/pyzx.py @@ -31,7 +31,8 @@ def _fraction_of_angle(angle: ExpressionOrFloat) -> Fraction: return Fraction(angle) -def to_pyzx_graph(og: OpenGraph) -> BaseGraph[int, tuple[int, int]]: +# TODO: Adapt to new OpenGraph API +def to_pyzx_graph(og: OpenGraph[Measurement]) -> BaseGraph[int, tuple[int, int]]: """Return a :mod:`pyzx` graph corresponding to the open graph. Example @@ -42,7 +43,7 @@ def to_pyzx_graph(og: OpenGraph) -> BaseGraph[int, tuple[int, int]]: >>> inputs = [0] >>> outputs = [2] >>> measurements = {0: Measurement(0, Plane.XY), 1: Measurement(1, Plane.YZ)} - >>> og = OpenGraph(g, measurements, inputs, outputs) + >>> og = OpenGraph(g, inputs, outputs, measurements) >>> reconstructed_pyzx_graph = to_pyzx_graph(og) """ if zx.__version__ != "0.9.0": @@ -61,11 +62,11 @@ def add_vertices(n: int, ty: VertexType) -> list[VertexType]: return verts # Add input boundary nodes - in_verts = add_vertices(len(og.inputs), VertexType.BOUNDARY) + in_verts = add_vertices(len(og.input_nodes), VertexType.BOUNDARY) g.set_inputs(tuple(in_verts)) # Add nodes for internal Z spiders - not including the phase gadgets - body_verts = add_vertices(len(og.inside), VertexType.Z) + body_verts = add_vertices(len(og.graph), VertexType.Z) # Add nodes for the phase gadgets. In OpenGraph we don't store the # effect as a separate node, it is instead just stored in the @@ -73,23 +74,23 @@ def add_vertices(n: int, ty: VertexType) -> list[VertexType]: x_meas = [i for i, m in og.measurements.items() if m.plane == Plane.YZ] x_meas_verts = add_vertices(len(x_meas), VertexType.Z) - out_verts = add_vertices(len(og.outputs), VertexType.BOUNDARY) + out_verts = add_vertices(len(og.output_nodes), VertexType.BOUNDARY) g.set_outputs(tuple(out_verts)) # Maps a node's ID in the Open Graph to it's corresponding node ID in # the PyZX graph and vice versa. - map_to_og = dict(zip(body_verts, og.inside.nodes())) + map_to_og = dict(zip(body_verts, og.graph.nodes())) map_to_pyzx = {v: i for i, v in map_to_og.items()} # Open Graph's don't have boundary nodes, so we need to connect the # input and output Z spiders to their corresponding boundary nodes in # pyzx. - for pyzx_index, og_index in zip(in_verts, og.inputs): + for pyzx_index, og_index in zip(in_verts, og.input_nodes): g.add_edge((pyzx_index, map_to_pyzx[og_index])) - for pyzx_index, og_index in zip(out_verts, og.outputs): + for pyzx_index, og_index in zip(out_verts, og.output_nodes): g.add_edge((pyzx_index, map_to_pyzx[og_index])) - og_edges = og.inside.edges() + og_edges = og.graph.edges() pyzx_edges = ((map_to_pyzx[a], map_to_pyzx[b]) for a, b in og_edges) g.add_edges(pyzx_edges, EdgeType.HADAMARD) @@ -114,7 +115,7 @@ def _checked_float(x: FractionLike) -> float: return float(x) -def from_pyzx_graph(g: BaseGraph[int, tuple[int, int]]) -> OpenGraph: +def from_pyzx_graph(g: BaseGraph[int, tuple[int, int]]) -> OpenGraph[Measurement]: """Construct an :class:`OpenGraph` from a :mod:`pyzx` graph. This method may add additional nodes to the graph so that it adheres @@ -193,4 +194,4 @@ def from_pyzx_graph(g: BaseGraph[int, tuple[int, int]]) -> OpenGraph: # expects a float measurements[v] = Measurement(-_checked_float(g.phase(v)), Plane.XY) - return OpenGraph(g_nx, measurements, inputs, outputs) + return OpenGraph(g_nx, inputs, outputs, measurements) diff --git a/tests/test_find_pflow.py b/tests/test_find_pflow.py deleted file mode 100644 index 18541813..00000000 --- a/tests/test_find_pflow.py +++ /dev/null @@ -1,682 +0,0 @@ -from __future__ import annotations - -from typing import TYPE_CHECKING, NamedTuple - -import networkx as nx -import numpy as np -import pytest - -from graphix._linalg import MatGF2 -from graphix.find_pflow import ( - OpenGraphIndex, - _compute_pflow_matrices, - _compute_reduced_adj, - _compute_topological_generations, - _find_pflow_simple, - find_pflow, -) -from graphix.fundamentals import Plane -from graphix.generator import _pflow2pattern -from graphix.measurements import Measurement -from graphix.opengraph import OpenGraph -from graphix.parameter import Placeholder -from graphix.random_objects import rand_circuit -from graphix.states import PlanarState -from tests.conftest import fx_rng - -if TYPE_CHECKING: - from numpy.random import Generator - from pytest_benchmark import BenchmarkFixture - - -class OpenGraphTestCase(NamedTuple): - ogi: OpenGraphIndex - radj: MatGF2 | None - flow_demand_mat: MatGF2 | None - order_demand_mat: MatGF2 | None - has_pflow: bool - - -class DAGTestCase(NamedTuple): - adj_mat: MatGF2 - generations: list[list[int]] | None - - -def get_og_rndcircuit(depth: int, n_qubits: int, n_inputs: int | None = None) -> OpenGraph: - """Return an open graph from a random circuit. - - Parameters - ---------- - depth : int - Circuit depth of the random circuits for generating open graphs. - n_qubits : int - Number of qubits in the random circuits for generating open graphs. It controls the number of outputs. - n_inputs : int | None - Optional (default to `None`). Maximum number of inputs in the returned open graph. The returned open graph is the open graph generated from the random circuit where `n_qubits - n_inputs` nodes have been removed from the input-nodes set. This operation does not change the flow properties of the graph. - - Returns - ------- - OpenGraph - Open graph with causal flow. - """ - circuit = rand_circuit(n_qubits, depth, fx_rng._fixture_function()) - pattern = circuit.transpile().pattern - graph = pattern.extract_graph() - - angles = pattern.get_angles() - planes = pattern.get_meas_plane() - meas = {node: Measurement(angle, planes[node]) for node, angle in angles.items()} - - og = OpenGraph( - inside=graph, - inputs=pattern.input_nodes, - outputs=pattern.output_nodes, - measurements=meas, - ) - - if n_inputs is not None: - ni_remove = max(0, n_qubits - n_inputs) - for i in range(ni_remove): - og.inputs.remove(i) - - return og - - -def get_og_dense(ni: int, no: int, m: int) -> OpenGraph: - """Return a dense open graph with causal, gflow and pflow. - - Parameters - ---------- - ni : int - Number of input nodes (must be equal or smaller than `no` ). - no : int - Number of output nodes (must be larger than 1). - m : int - Number of total nodes (it must satisfy `m - 2*no > 0`). - - Returns - ------- - OpenGraph - Open graph with causal and gflow. - - Notes - ----- - Adapted from Fig. 1 in Houshmand et al., Phys. Rev. A, 98 (2018) (arXiv:1705.01535) - """ - if no <= 1: - raise ValueError("Number of outputs must be larger than 1 (no > 1).") - if m - 2 * no <= 0: - raise ValueError("Total number of nodes must be larger than twice the number of outputs (m - 2no > 0).") - - inputs = list(range(no)) # we remove inputs afterwards - outputs = list(range(no, 2 * no)) - edges = [(i, o) for i, o in zip(inputs[:-2], outputs[:-2])] - edges.extend((node, node + 1) for node in range(2 * no - 1, m - 1)) - edges.append((inputs[-2], m - 1)) - - graph: nx.Graph[int] = nx.Graph() - graph.add_nodes_from(range(m)) - graph.add_edges_from(edges) - graph_c = nx.complement(graph) - - meas = {node: Measurement(Placeholder("Angle"), Plane.XY) for node in range(m) if node not in set(outputs)} - - og = OpenGraph( - inside=graph_c, - inputs=inputs, - outputs=outputs, - measurements=meas, - ) # This open graph corresponds to the example in the reference. Now we remove nodes from the set of inputs, since this operation preserves the flow properties. - - ni_remove = max(0, len(og.inputs) - ni) - for i in og.inputs[ni_remove:]: - og.inputs.remove(i) - return og - - -def prepare_test_og() -> list[OpenGraphTestCase]: - test_cases: list[OpenGraphTestCase] = [] - - # Trivial open graph with pflow and nI = nO - def get_og_0() -> OpenGraph: - """Return an open graph with Pauli flow and equal number of outputs and inputs. - - The returned graph has the following structure: - - [0]-1-(2) - """ - graph: nx.Graph[int] = nx.Graph([(0, 1), (1, 2)]) - inputs = [0] - outputs = [2] - meas = { - 0: Measurement(0.1, Plane.XY), # XY - 1: Measurement(0.5, Plane.YZ), # Y - } - return OpenGraph(inside=graph, inputs=inputs, outputs=outputs, measurements=meas) - - test_cases.append( - OpenGraphTestCase( - ogi=OpenGraphIndex(get_og_0()), - radj=MatGF2([[1, 0], [0, 1]]), - flow_demand_mat=MatGF2([[1, 0], [1, 1]]), - order_demand_mat=MatGF2([[0, 0], [0, 0]]), - has_pflow=True, - ) - ) - - # Non-trivial open graph without pflow and nI = nO - def get_og_1() -> OpenGraph: - """Return an open graph without Pauli flow and equal number of outputs and inputs. - - The returned graph has the following structure: - - [0]-2-4-(6) - | | - [1]-3-5-(7) - """ - graph: nx.Graph[int] = nx.Graph([(0, 2), (1, 3), (2, 3), (2, 4), (3, 5), (4, 5), (4, 6), (5, 7)]) - inputs = [1, 0] - outputs = [6, 7] - meas = { - 0: Measurement(0.1, Plane.XY), # XY - 1: Measurement(0.1, Plane.XZ), # XZ - 2: Measurement(0.5, Plane.XZ), # X - 3: Measurement(0.5, Plane.YZ), # Y - 4: Measurement(0.5, Plane.YZ), # Y - 5: Measurement(0.1, Plane.YZ), # YZ - } - return OpenGraph(inside=graph, inputs=inputs, outputs=outputs, measurements=meas) - - test_cases.append( - OpenGraphTestCase( - ogi=OpenGraphIndex(get_og_1()), - radj=MatGF2( - [ - [1, 0, 0, 0, 0, 0], - [0, 1, 0, 0, 0, 0], - [0, 1, 1, 0, 0, 0], - [1, 0, 0, 1, 0, 0], - [1, 0, 0, 1, 1, 0], - [0, 1, 1, 0, 0, 1], - ] - ), - flow_demand_mat=MatGF2( - [ - [1, 0, 0, 0, 0, 0], - [0, 0, 0, 0, 0, 0], - [0, 1, 1, 0, 0, 0], - [1, 1, 0, 1, 0, 0], - [1, 0, 1, 1, 1, 0], - [0, 0, 0, 1, 0, 0], - ] - ), - order_demand_mat=MatGF2( - [ - [0, 0, 0, 0, 0, 0], - [0, 1, 0, 0, 0, 0], - [0, 0, 0, 0, 0, 0], - [0, 0, 0, 0, 0, 0], - [0, 0, 0, 0, 0, 0], - [0, 1, 1, 0, 0, 1], - ] - ), - has_pflow=False, - ) - ) - - # Non-trivial open graph with pflow and nI = nO - def get_og_2() -> OpenGraph: - """Return an open graph with Pauli flow and equal number of outputs and inputs. - - The returned graph has the following structure: - - [0]-2-4-(6) - | | - [1]-3-5-(7) - """ - graph: nx.Graph[int] = nx.Graph([(0, 2), (1, 3), (2, 3), (2, 4), (3, 5), (4, 5), (4, 6), (5, 7)]) - inputs = [0, 1] - outputs = [6, 7] - meas = { - 0: Measurement(0.1, Plane.XY), # XY - 1: Measurement(0.1, Plane.XY), # XY - 2: Measurement(0.0, Plane.XY), # X - 3: Measurement(0.1, Plane.XY), # XY - 4: Measurement(0.0, Plane.XY), # X - 5: Measurement(0.5, Plane.XY), # Y - } - return OpenGraph(inside=graph, inputs=inputs, outputs=outputs, measurements=meas) - - test_cases.append( - OpenGraphTestCase( - ogi=OpenGraphIndex(get_og_2()), - radj=MatGF2( - [ - [1, 0, 0, 0, 0, 0], - [0, 1, 0, 0, 0, 0], - [0, 1, 1, 0, 0, 0], - [1, 0, 0, 1, 0, 0], - [1, 0, 0, 1, 1, 0], - [0, 1, 1, 0, 0, 1], - ] - ), - flow_demand_mat=MatGF2( - [ - [1, 0, 0, 0, 0, 0], - [0, 1, 0, 0, 0, 0], - [0, 1, 1, 0, 0, 0], - [1, 0, 0, 1, 0, 0], - [1, 0, 0, 1, 1, 0], - [0, 1, 1, 1, 0, 1], - ] - ), - order_demand_mat=MatGF2( - [ - [0, 0, 0, 0, 0, 0], - [0, 0, 0, 0, 0, 0], - [0, 0, 0, 0, 0, 0], - [0, 1, 0, 0, 0, 0], - [0, 0, 0, 0, 0, 0], - [0, 0, 0, 0, 0, 0], - ] - ), - has_pflow=True, - ) - ) - - # Non-trivial open graph with pflow and nI != nO - def get_og_3() -> OpenGraph: - """Return an open graph with Pauli flow and unequal number of outputs and inputs. - - Example from Fig. 1 in Mitosek and Backens, 2024 (arXiv:2410.23439). - """ - graph: nx.Graph[int] = nx.Graph( - [(0, 2), (2, 4), (3, 4), (4, 6), (1, 4), (1, 6), (2, 3), (3, 5), (2, 6), (3, 6)] - ) - inputs = [0] - outputs = [5, 6] - meas = { - 0: Measurement(0.1, Plane.XY), # XY - 1: Measurement(0.1, Plane.XZ), # XZ - 2: Measurement(0.5, Plane.YZ), # Y - 3: Measurement(0.1, Plane.XY), # XY - 4: Measurement(0, Plane.XZ), # Z - } - - return OpenGraph(inside=graph, inputs=inputs, outputs=outputs, measurements=meas) - - test_cases.append( - OpenGraphTestCase( - ogi=OpenGraphIndex(get_og_3()), - radj=MatGF2( - [[0, 1, 0, 0, 0, 0], [0, 0, 0, 1, 0, 1], [0, 0, 1, 1, 0, 1], [0, 1, 0, 1, 1, 1], [1, 1, 1, 0, 0, 1]] - ), - flow_demand_mat=MatGF2( - [[0, 1, 0, 0, 0, 0], [1, 0, 0, 0, 0, 0], [0, 1, 1, 1, 0, 1], [0, 1, 0, 1, 1, 1], [0, 0, 0, 1, 0, 0]] - ), - order_demand_mat=MatGF2( - [[0, 0, 0, 0, 0, 0], [1, 0, 0, 1, 0, 1], [0, 0, 0, 0, 0, 0], [0, 0, 1, 0, 0, 0], [0, 0, 0, 0, 0, 0]] - ), - has_pflow=True, - ) - ) - - # The following tests check the final result only, not the intermediate steps. - - # Non-trivial open graph with pflow and nI != nO - def get_og_4() -> OpenGraph: - """Return an open graph with Pauli flow and unequal number of outputs and inputs.""" - graph: nx.Graph[int] = nx.Graph([(0, 2), (1, 3), (2, 3), (2, 6), (3, 4), (4, 7), (4, 5), (7, 8)]) - inputs = [0, 1] - outputs = [5, 6, 8] - meas = { - 0: Measurement(0.1, Plane.XY), - 1: Measurement(0.1, Plane.XY), - 2: Measurement(0.0, Plane.XY), - 3: Measurement(0, Plane.XY), - 4: Measurement(0.5, Plane.XY), - 7: Measurement(0, Plane.XY), - } - - return OpenGraph(inside=graph, inputs=inputs, outputs=outputs, measurements=meas) - - test_cases.append( - OpenGraphTestCase( - ogi=OpenGraphIndex(get_og_4()), - radj=None, - flow_demand_mat=None, - order_demand_mat=None, - has_pflow=True, - ) - ) - - # Non-trivial open graph with pflow and nI != nO - def get_og_5() -> OpenGraph: - """Return an open graph with Pauli flow and unequal number of outputs and inputs.""" - graph: nx.Graph[int] = nx.Graph([(0, 2), (1, 2), (2, 3), (3, 4)]) - inputs = [0, 1] - outputs = [1, 3, 4] - meas = {0: Measurement(0.1, Plane.XY), 2: Measurement(0.5, Plane.YZ)} - - return OpenGraph(inside=graph, inputs=inputs, outputs=outputs, measurements=meas) - - test_cases.append( - OpenGraphTestCase( - ogi=OpenGraphIndex(get_og_5()), - radj=None, - flow_demand_mat=None, - order_demand_mat=None, - has_pflow=True, - ) - ) - - # Non-trivial open graph with pflow and nI != nO - def get_og_6() -> OpenGraph: - """Return an open graph with Pauli flow and unequal number of outputs and inputs.""" - graph: nx.Graph[int] = nx.Graph([(0, 1), (0, 3), (1, 4), (3, 4), (2, 3), (2, 5), (3, 6), (4, 7)]) - inputs = [1] - outputs = [6, 2, 7] - meas = { - 0: Measurement(0.1, Plane.XZ), - 1: Measurement(0.1, Plane.XY), - 3: Measurement(0, Plane.XY), - 4: Measurement(0.1, Plane.XY), - 5: Measurement(0.1, Plane.YZ), - } - - return OpenGraph(inside=graph, inputs=inputs, outputs=outputs, measurements=meas) - - test_cases.append( - OpenGraphTestCase( - ogi=OpenGraphIndex(get_og_6()), - radj=None, - flow_demand_mat=None, - order_demand_mat=None, - has_pflow=True, - ) - ) - - # Disconnected open graph with pflow and nI != nO - def get_og_7() -> OpenGraph: - """Return an open graph with Pauli flow and unequal number of outputs and inputs.""" - graph: nx.Graph[int] = nx.Graph([(0, 1), (0, 2), (2, 3), (1, 3), (4, 6)]) - inputs: list[int] = [] - outputs = [1, 3, 4] - meas = {0: Measurement(0.5, Plane.XZ), 2: Measurement(0, Plane.YZ), 6: Measurement(0.2, Plane.XY)} - - return OpenGraph(inside=graph, inputs=inputs, outputs=outputs, measurements=meas) - - test_cases.append( - OpenGraphTestCase( - ogi=OpenGraphIndex(get_og_7()), - radj=None, - flow_demand_mat=None, - order_demand_mat=None, - has_pflow=True, - ) - ) - - # Non-trivial open graph without pflow and nI != nO - def get_og_8() -> OpenGraph: - """Return an open graph without Pauli flow and unequal number of outputs and inputs.""" - graph: nx.Graph[int] = nx.Graph( - [(0, 1), (0, 3), (1, 4), (3, 4), (2, 3), (2, 5), (3, 6), (4, 7), (5, 6), (6, 7)] - ) - inputs = [1] - outputs = [6, 2, 7] - meas = { - 0: Measurement(0.1, Plane.XZ), - 1: Measurement(0.1, Plane.XY), - 3: Measurement(0, Plane.XY), - 4: Measurement(0.1, Plane.XY), - 5: Measurement(0.1, Plane.XY), - } - - return OpenGraph(inside=graph, inputs=inputs, outputs=outputs, measurements=meas) - - test_cases.append( - OpenGraphTestCase( - ogi=OpenGraphIndex(get_og_8()), - radj=None, - flow_demand_mat=None, - order_demand_mat=None, - has_pflow=False, - ) - ) - - # Disconnected open graph without pflow and nI != nO - def get_og_9() -> OpenGraph: - """Return an open graph without Pauli flow and unequal number of outputs and inputs.""" - graph: nx.Graph[int] = nx.Graph([(0, 1), (0, 2), (2, 3), (1, 3), (4, 6)]) - inputs = [0] - outputs = [1, 3, 4] - meas = {0: Measurement(0.1, Plane.XZ), 2: Measurement(0, Plane.YZ), 6: Measurement(0.2, Plane.XY)} - - return OpenGraph(inside=graph, inputs=inputs, outputs=outputs, measurements=meas) - - test_cases.append( - OpenGraphTestCase( - ogi=OpenGraphIndex(get_og_9()), - radj=None, - flow_demand_mat=None, - order_demand_mat=None, - has_pflow=False, - ) - ) - - # Non-trivial open graph without pflow and nI != nO - def get_og_10() -> OpenGraph: - """Return a graph constructed by adding a disconnected input to graph_6. The resulting graph does not have pflow.""" - graph: nx.Graph[int] = nx.Graph([(0, 1), (0, 3), (1, 4), (3, 4), (2, 3), (2, 5), (3, 6), (4, 7)]) - graph.add_node(8) - inputs = [1, 8] - outputs = [6, 2, 7] - meas = { - 0: Measurement(0.1, Plane.XZ), - 1: Measurement(0.1, Plane.XY), - 3: Measurement(0, Plane.XY), - 4: Measurement(0.1, Plane.XY), - 5: Measurement(0.1, Plane.YZ), - 8: Measurement(0.1, Plane.XY), - } - - return OpenGraph(inside=graph, inputs=inputs, outputs=outputs, measurements=meas) - - test_cases.append( - OpenGraphTestCase( - ogi=OpenGraphIndex(get_og_10()), - radj=None, - flow_demand_mat=None, - order_demand_mat=None, - has_pflow=False, - ) - ) - - # Open graph with only Pauli measurements, without pflow and nI != nO - def get_og_11() -> OpenGraph: - """Return an open graph without Pauli flow and unequal number of outputs and inputs.""" - graph: nx.Graph[int] = nx.Graph([(0, 2), (1, 3), (2, 3), (2, 6), (3, 4), (4, 7), (4, 5), (7, 8)]) - inputs = [0, 1] - outputs = [5, 6, 8] - meas = { - 0: Measurement(0, Plane.XY), - 1: Measurement(0, Plane.XY), - 2: Measurement(0, Plane.XZ), - 3: Measurement(0, Plane.XY), - 4: Measurement(0.5, Plane.XY), - 7: Measurement(0, Plane.YZ), - } - - return OpenGraph(inside=graph, inputs=inputs, outputs=outputs, measurements=meas) - - test_cases.append( - OpenGraphTestCase( - ogi=OpenGraphIndex(get_og_11()), - radj=None, - flow_demand_mat=None, - order_demand_mat=None, - has_pflow=False, - ) - ) - - # Open graph with only Pauli measurements, with pflow and nI != nO - def get_og_12() -> OpenGraph: - """Return an open graph with Pauli flow and unequal number of outputs and inputs. Even though all nodes are Pauli-measured, open graph has flow because none of them are inputs.""" - graph: nx.Graph[int] = nx.Graph([(0, 2), (1, 3), (2, 3), (2, 6), (3, 4), (4, 7), (4, 5), (7, 8)]) - outputs = [5, 6, 8] - meas = { - 0: Measurement(0, Plane.XZ), - 1: Measurement(0, Plane.XZ), - 2: Measurement(0, Plane.XZ), - 3: Measurement(0, Plane.XZ), - 4: Measurement(0, Plane.XZ), - 7: Measurement(0, Plane.XZ), - } - - return OpenGraph(inside=graph, inputs=[], outputs=outputs, measurements=meas) - - test_cases.append( - OpenGraphTestCase( - ogi=OpenGraphIndex(get_og_12()), - radj=None, - flow_demand_mat=None, - order_demand_mat=None, - has_pflow=True, - ) - ) - - return test_cases - - -def prepare_benchmark_og() -> list[OpenGraphTestCase]: - test_cases: list[OpenGraphTestCase] = [] - - # Open graph from random circuit - test_cases.extend( - ( - OpenGraphTestCase( - ogi=OpenGraphIndex(get_og_rndcircuit(depth=20, n_qubits=7, n_inputs=1)), - radj=None, - flow_demand_mat=None, - order_demand_mat=None, - has_pflow=True, - ), - OpenGraphTestCase( - ogi=OpenGraphIndex(get_og_dense(ni=3, no=6, m=400)), - radj=None, - flow_demand_mat=None, - order_demand_mat=None, - has_pflow=True, - ), - ) - ) - return test_cases - - -def prepare_test_dag() -> list[DAGTestCase]: - test_cases: list[DAGTestCase] = [] - - # Simple DAG - test_cases.extend( - ( # Simple DAG - DAGTestCase( - adj_mat=MatGF2([[0, 0, 0, 0], [1, 0, 0, 0], [1, 0, 0, 0], [0, 1, 1, 0]]), generations=[[0], [1, 2], [3]] - ), - # Graph with loop - DAGTestCase(adj_mat=MatGF2([[0, 0, 0, 0], [1, 0, 0, 1], [1, 0, 0, 0], [0, 1, 1, 0]]), generations=None), - # Disconnected graph - DAGTestCase( - adj_mat=MatGF2([[0, 0, 0, 0, 0], [1, 0, 0, 0, 0], [0, 0, 0, 0, 0], [0, 0, 1, 0, 0], [0, 0, 1, 0, 0]]), - generations=[[0, 2], [1, 3, 4]], - ), - ) - ) - - return test_cases - - -class TestPflow: - @pytest.mark.parametrize("test_case", prepare_test_og()) - def test_compute_reduced_adj(self, test_case: OpenGraphTestCase) -> None: - if test_case.radj is not None: - ogi = test_case.ogi - radj = _compute_reduced_adj(ogi) - assert np.all(radj == test_case.radj) - - @pytest.mark.parametrize("test_case", prepare_test_og()) - def test_compute_pflow_matrices(self, test_case: OpenGraphTestCase) -> None: - if test_case.flow_demand_mat is not None and test_case.order_demand_mat is not None: - ogi = test_case.ogi - flow_demand_matrix, order_demand_matrix = _compute_pflow_matrices(ogi) - - assert np.all(flow_demand_matrix == test_case.flow_demand_mat) - assert np.all(order_demand_matrix == test_case.order_demand_mat) - - @pytest.mark.parametrize("test_case", prepare_test_og()) - def test_find_pflow_simple(self, test_case: OpenGraphTestCase) -> None: - if test_case.flow_demand_mat is not None: - ogi = test_case.ogi - if len(ogi.og.outputs) == len(ogi.og.inputs): - pflow_algebraic = _find_pflow_simple(ogi) - - if not test_case.has_pflow: - assert pflow_algebraic is None - else: - assert pflow_algebraic is not None - correction_matrix, _ = pflow_algebraic - ident = MatGF2(np.eye(len(ogi.non_outputs), dtype=np.uint8)) - assert np.all( - (test_case.flow_demand_mat @ correction_matrix) % 2 == ident - ) # Test with numpy matrix product. - - @pytest.mark.parametrize("test_case", prepare_test_og()) - def test_find_pflow_determinism(self, test_case: OpenGraphTestCase, fx_rng: Generator) -> None: - og = test_case.ogi.og - - pflow = find_pflow(og) - - if not test_case.has_pflow: - assert pflow is None - else: - assert pflow is not None - - pattern = _pflow2pattern( - graph=og.inside, - inputs=set(og.inputs), - meas_planes={i: m.plane for i, m in og.measurements.items()}, - angles={i: m.angle for i, m in og.measurements.items()}, - p=pflow[0], - l_k=pflow[1], - ) - pattern.reorder_output_nodes(og.outputs) - - alpha = 2 * np.pi * fx_rng.random() - state_ref = pattern.simulate_pattern(input_state=PlanarState(Plane.XY, alpha)) - - n_shots = 5 - results = [] - for _ in range(n_shots): - state = pattern.simulate_pattern(input_state=PlanarState(Plane.XY, alpha)) - results.append(np.abs(np.dot(state.flatten().conjugate(), state_ref.flatten()))) - - avg = sum(results) / n_shots - - assert avg == pytest.approx(1) - - @pytest.mark.parametrize("test_case", prepare_benchmark_og()) - def test_benchmark_pflow(self, benchmark: BenchmarkFixture, test_case: OpenGraphTestCase) -> None: - og = test_case.ogi.og - - pflow = benchmark(find_pflow, og) - - if not test_case.has_pflow: - assert pflow is None - else: - assert pflow is not None - - @pytest.mark.parametrize("test_case", prepare_test_dag()) - def test_compute_topological_generations(self, test_case: DAGTestCase) -> None: - adj_mat = test_case.adj_mat - generations_ref = test_case.generations - - assert generations_ref == _compute_topological_generations(adj_mat) diff --git a/tests/test_flow_core.py b/tests/test_flow_core.py new file mode 100644 index 00000000..9b15a9a7 --- /dev/null +++ b/tests/test_flow_core.py @@ -0,0 +1,585 @@ +from __future__ import annotations + +from typing import TYPE_CHECKING, NamedTuple + +import networkx as nx +import numpy as np +import pytest + +from graphix.command import E, M, N, X, Z +from graphix.flow.core import CausalFlow, GFlow, PauliFlow, XZCorrections +from graphix.fundamentals import AbstractMeasurement, AbstractPlanarMeasurement, Axis, Plane +from graphix.measurements import Measurement +from graphix.opengraph import OpenGraph +from graphix.pattern import Pattern +from graphix.states import PlanarState + +if TYPE_CHECKING: + from collections.abc import Mapping + + from numpy.random import Generator + + +def generate_causal_flow_0() -> CausalFlow[Plane]: + """Generate causal flow on linear open graph. + + Open graph structure: + + [0]-1-2-(3) + + Causal flow: + c(0) = 1, c(1) = 2, c(2) = 3 + {3} > {2} > {1} > {0} + """ + og = OpenGraph( + graph=nx.Graph([(0, 1), (1, 2), (2, 3)]), + input_nodes=[0], + output_nodes=[3], + measurements=dict.fromkeys(range(3), Plane.XY), + ) + return CausalFlow( + og=og, + correction_function={0: {1}, 1: {2}, 2: {3}}, + partial_order_layers=[{3}, {2}, {1}, {0}], + ) + + +def generate_causal_flow_1() -> CausalFlow[Measurement]: + """Generate causal flow on H-shaped open graph. + + Open graph structure: + + [0]-2-(4) + | + [1]-3-(5) + + Causal flow: + c(0) = 2, c(1) = 3, c(2) = 4, c(3) = 5 + {4, 5} > {2, 3} > {0, 1} + """ + og = OpenGraph( + graph=nx.Graph([(0, 2), (2, 3), (1, 3), (2, 4), (3, 5)]), + input_nodes=[0, 1], + output_nodes=[4, 5], + measurements=dict.fromkeys(range(4), Measurement(angle=0, plane=Plane.XY)), + ) + return CausalFlow( + og=og, + correction_function={0: {2}, 1: {3}, 2: {4}, 3: {5}}, + partial_order_layers=[{4, 5}, {2, 3}, {0, 1}], + ) + + +def generate_gflow_0() -> GFlow[Measurement]: + """Generate gflow on H-shaped open graph. + + Open graph structure: + + [0]-2-(4) + | + [1]-3-(5) + + GFlow: + g(0) = {2, 5}, g(1) = {3, 4}, g(2) = {4}, g(3) = {5} + {4, 5} > {0, 1, 2, 3} + + Notes + ----- + This is the same open graph as in `:func: generate_causal_flow_1` but now we consider a gflow which has lower depth than the causal flow. + """ + og = OpenGraph( + graph=nx.Graph([(0, 2), (2, 3), (1, 3), (2, 4), (3, 5)]), + input_nodes=[0, 1], + output_nodes=[4, 5], + measurements=dict.fromkeys(range(4), Measurement(angle=0, plane=Plane.XY)), + ) + return GFlow( + og=og, + correction_function={0: {2, 5}, 1: {3, 4}, 2: {4}, 3: {5}}, + partial_order_layers=[{4, 5}, {0, 1, 2, 3}], + ) + + +def generate_gflow_1() -> GFlow[Plane]: + r"""Generate gflow on open graph without causal flow. + + Open graph structure: + + 1 + \ + (4)-[0]-(3) + / + 2 + + GFlow: + g(0) = {3}, g(1) = {1}, g(2) = {2, 3, 4} + {3, 4} > {1} > {0, 2} + """ + og = OpenGraph( + graph=nx.Graph([(0, 3), (0, 4), (1, 4), (2, 4)]), + input_nodes=[0], + output_nodes=[3, 4], + measurements={0: Plane.XY, 1: Plane.YZ, 2: Plane.XZ}, + ) + return GFlow( + og=og, + correction_function={0: {3}, 1: {1}, 2: {2, 3, 4}}, + partial_order_layers=[{3, 4}, {1}, {0, 2}], + ) + + +def generate_gflow_2() -> GFlow[Plane]: + r"""Generate gflow on open graph without causal flow. + + Open graph structure: + + [0]-(3) + X + [1]-(4) + X + [2]-(5) + + GFlow: + g(0) = {4, 5}, g(1) = {3, 4, 5}, g(2) = {3, 4} + {3, 4, 5} > {0, 1, 2} + """ + og = OpenGraph( + graph=nx.Graph([(0, 3), (0, 4), (1, 3), (1, 4), (1, 5), (2, 4), (2, 5)]), + input_nodes=[0, 1, 2], + output_nodes=[3, 4, 5], + measurements=dict.fromkeys(range(3), Plane.XY), + ) + return GFlow( + og=og, + correction_function={0: {4, 5}, 1: {3, 4, 5}, 2: {3, 4}}, + partial_order_layers=[{3, 4}, {1}, {0, 2}], + ) + + +def generate_pauli_flow_0() -> PauliFlow[Axis]: + """Generate Pauli flow on linear open graph. + + Open graph structure: + + [0]-1-2-(3) + + Pauli flow: + p(0) = {1, 3}, p(1) = {2}, p(2) = {3} + {3} > {0, 1, 2} + """ + og = OpenGraph( + graph=nx.Graph([(0, 1), (1, 2), (2, 3)]), + input_nodes=[0], + output_nodes=[3], + measurements=dict.fromkeys(range(3), Axis.X), + ) + return PauliFlow( + og=og, + correction_function={0: {1, 3}, 1: {2}, 2: {3}}, + partial_order_layers=[{3}, {0, 1, 2}], + ) + + +def generate_pauli_flow_1() -> PauliFlow[Measurement]: + """Generate Pauli flow on double-H-shaped open graph. + + Open graph structure: + + [0]-2-4-(6) + | | + [1]-3-5-(7) + + Pauli flow: + p(0) = {2, 5, 7}, p(1) = {3, 4}, p(2) = {4, 7}, p(3) = {5, 6, 7}, + p(4) = {6}, p(5) = 7 + {6, 7} > {3} > {0, 1, 2, 4, 5} + """ + og = OpenGraph( + graph=nx.Graph([(0, 2), (1, 3), (2, 3), (2, 4), (3, 5), (4, 5), (4, 6), (5, 7)]), + input_nodes=[0, 1], + output_nodes=[6, 7], + measurements={ + 0: Measurement(0.1, Plane.XY), # XY + 1: Measurement(0.1, Plane.XY), # XY + 2: Measurement(0.0, Plane.XY), # X + 3: Measurement(0.1, Plane.XY), # XY + 4: Measurement(0.0, Plane.XY), # X + 5: Measurement(0.5, Plane.XY), # Y + }, + ) + return PauliFlow( + og=og, + correction_function={0: {2, 5, 7}, 1: {3, 4}, 2: {4, 7}, 3: {5, 6, 7}, 4: {6}, 5: {7}}, + partial_order_layers=[{6, 7}, {3}, {0, 1, 2, 4, 5}], + ) + + +class XZCorrectionsTestCase(NamedTuple): + flow: CausalFlow[AbstractPlanarMeasurement] | GFlow[AbstractPlanarMeasurement] | PauliFlow[AbstractMeasurement] + x_corr: Mapping[int, set[int]] + z_corr: Mapping[int, set[int]] + pattern: Pattern | None + # Patterns can only be extracted from `Measurement`-type objects. If `flow` is of parametric type, we set `pattern = None`. + + +def prepare_test_xzcorrections() -> list[XZCorrectionsTestCase]: + test_cases: list[XZCorrectionsTestCase] = [] + + test_cases.extend( + ( + XZCorrectionsTestCase( + flow=generate_causal_flow_0(), x_corr={0: {1}, 1: {2}, 2: {3}}, z_corr={0: {2}, 1: {3}}, pattern=None + ), + XZCorrectionsTestCase( + flow=generate_causal_flow_1(), + x_corr={0: {2}, 1: {3}, 2: {4}, 3: {5}}, + z_corr={0: {3, 4}, 1: {2, 5}}, + pattern=Pattern( + input_nodes=[0, 1], + cmds=[ + N(2), + N(3), + N(4), + N(5), + E((0, 2)), + E((2, 3)), + E((2, 4)), + E((3, 1)), + E((3, 5)), + M(0), + Z(3, {0}), + Z(4, {0}), + X(2, {0}), + M(1), + Z(2, {1}), + Z(5, {1}), + X(3, {1}), + M(2), + X(4, {2}), + M(3), + X(5, {3}), + ], + output_nodes=[4, 5], + ), + ), + XZCorrectionsTestCase( + flow=generate_gflow_0(), + x_corr={0: {2, 5}, 1: {3, 4}, 2: {4}, 3: {5}}, + z_corr={0: {4}, 1: {5}}, + pattern=Pattern( + input_nodes=[0, 1], + cmds=[ + N(2), + N(3), + N(4), + N(5), + E((0, 2)), + E((2, 3)), + E((2, 4)), + E((3, 1)), + E((3, 5)), + M(0), + Z(3, {0}), + Z(4, {0}), + X(2, {0}), + M(1), + Z(2, {1}), + Z(5, {1}), + X(3, {1}), + M(2), + X(4, {2}), + M(3), + X(5, {3}), + ], + output_nodes=[4, 5], + ), + ), + XZCorrectionsTestCase( + flow=generate_gflow_1(), + x_corr={0: {3}, 2: {3, 4}}, + z_corr={1: {4}, 2: {1, 4}}, + pattern=None, + ), + XZCorrectionsTestCase( + flow=generate_gflow_2(), x_corr={0: {4, 5}, 1: {3, 4, 5}, 2: {3, 4}}, z_corr={}, pattern=None + ), + XZCorrectionsTestCase( + flow=generate_pauli_flow_0(), + x_corr={0: {3}, 2: {3}}, + z_corr={1: {3}}, + pattern=None, + ), + XZCorrectionsTestCase( + flow=generate_pauli_flow_1(), + x_corr={0: {7}, 1: {3}, 2: {7}, 3: {6, 7}, 4: {6}, 5: {7}}, + z_corr={0: {7}, 1: {6}, 2: {6}, 3: {7}}, + pattern=Pattern( + input_nodes=[0, 1], + cmds=[ + N(2), + N(3), + N(4), + N(5), + N(6), + N(7), + E((0, 2)), + E((2, 3)), + E((2, 4)), + E((1, 3)), + E((3, 5)), + E((4, 5)), + E((4, 6)), + E((5, 7)), + M(0, angle=0.1), + Z(3, {0}), + Z(4, {0}), + X(2, {0}), + M(1, angle=0.1), + Z(2, {1}), + Z(5, {1}), + X(3, {1}), + M(2), + Z(5, {2}), + Z(6, {2}), + X(4, {2}), + M(3, angle=0.1), + Z(4, {3}), + Z(7, {3}), + X(5, {3}), + M(4), + X(6, {4}), + M(5, angle=0.5), + X(7, {5}), + ], + output_nodes=[6, 7], + ), + ), + ) + ) + + return test_cases + + +class TestFlowPatternConversion: + """Bundle for unit tests of the flow to XZ-corrections to pattern methods. + + The module `tests.test_opengraph.py` provides an additional (more comprehensive) suite of unit tests on this transformation. + """ + + @pytest.mark.parametrize("test_case", prepare_test_xzcorrections()) + def test_flow_to_corrections(self, test_case: XZCorrectionsTestCase) -> None: + flow = test_case.flow + corrections = flow.to_corrections() + assert corrections.z_corrections == test_case.z_corr + assert corrections.x_corrections == test_case.x_corr + + @pytest.mark.parametrize("test_case", prepare_test_xzcorrections()) + def test_corrections_to_pattern(self, test_case: XZCorrectionsTestCase, fx_rng: Generator) -> None: + if test_case.pattern is not None: + pattern = test_case.flow.to_corrections().to_pattern() # type: ignore[misc] + n_shots = 2 + results = [] + + for plane in {Plane.XY, Plane.XZ, Plane.YZ}: + alpha = 2 * np.pi * fx_rng.random() + state_ref = test_case.pattern.simulate_pattern(input_state=PlanarState(plane, alpha)) + + for _ in range(n_shots): + state = pattern.simulate_pattern(input_state=PlanarState(plane, alpha)) + results.append(np.abs(np.dot(state.flatten().conjugate(), state_ref.flatten()))) + + avg = sum(results) / (n_shots * 3) + + assert avg == pytest.approx(1) + + +class TestXZCorrections: + """Bundle for unit tests of :class:`XZCorrections`.""" + + # See `:func: generate_causal_flow_0` + def test_order_0(self) -> None: + corrections = generate_causal_flow_0().to_corrections() + + assert corrections.generate_total_measurement_order() == [0, 1, 2] + assert corrections.is_compatible([0, 1, 2]) # Correct order + assert not corrections.is_compatible([1, 0, 2]) # Wrong order + assert not corrections.is_compatible([1, 2]) # Incomplete order + assert not corrections.is_compatible([0, 1, 2, 3]) # Contains outputs + + assert nx.utils.graphs_equal(corrections.extract_dag(), nx.DiGraph([(0, 1), (0, 2), (1, 2), (2, 3), (1, 3)])) + + # See `:func: generate_causal_flow_1` + def test_order_1(self) -> None: + og = OpenGraph( + graph=nx.Graph([(0, 2), (2, 3), (1, 3), (2, 4), (3, 5)]), + input_nodes=[0, 1], + output_nodes=[4, 5], + measurements=dict.fromkeys(range(4), Measurement(angle=0, plane=Plane.XY)), + ) + + corrections = XZCorrections.from_measured_nodes_mapping( + og=og, x_corrections={0: {2}, 1: {3}, 2: {4}, 3: {5}}, z_corrections={0: {3, 4}, 1: {2, 5}} + ) + + assert corrections.is_compatible([0, 1, 2, 3]) + assert corrections.is_compatible([1, 0, 2, 3]) + assert corrections.is_compatible([1, 0, 3, 2]) + assert not corrections.is_compatible([0, 2, 1, 3]) # Wrong order + assert not corrections.is_compatible([1, 0, 3]) # Incomplete order + assert not corrections.is_compatible([0, 1, 1, 2, 3]) # Duplicates + assert not corrections.is_compatible([0, 1, 2, 3, 4, 5]) # Contains outputs + + assert nx.utils.graphs_equal( + corrections.extract_dag(), nx.DiGraph([(0, 2), (0, 3), (0, 4), (1, 2), (1, 3), (1, 5), (2, 4), (3, 5)]) + ) + + # Incomplete corrections + def test_order_2(self) -> None: + og = OpenGraph( + graph=nx.Graph([(0, 1), (1, 2), (1, 3)]), + input_nodes=[0], + output_nodes=[2, 3], + measurements=dict.fromkeys(range(2), Measurement(angle=0, plane=Plane.XY)), + ) + + corrections = XZCorrections.from_measured_nodes_mapping(og=og, x_corrections={1: {0}}) + + assert corrections.partial_order_layers == (frozenset({2, 3}), frozenset({0}), frozenset({1})) + assert corrections.is_compatible([1, 0]) + assert not corrections.is_compatible([0, 1]) # Wrong order + assert not corrections.is_compatible([0]) # Incomplete order + assert not corrections.is_compatible([0, 0, 1]) # Duplicates + assert not corrections.is_compatible([1, 0, 2, 3]) # Contains outputs + + assert nx.utils.graphs_equal(corrections.extract_dag(), nx.DiGraph([(1, 0)])) + + # OG without outputs + def test_order_3(self) -> None: + og = OpenGraph( + graph=nx.Graph([(0, 1), (1, 2)]), + input_nodes=[0], + output_nodes=[], + measurements=dict.fromkeys(range(3), Measurement(angle=0, plane=Plane.XY)), + ) + + corrections = XZCorrections.from_measured_nodes_mapping( + og=og, x_corrections={0: {1, 2}}, z_corrections={0: {1}} + ) + + assert corrections.partial_order_layers == (frozenset({1, 2}), frozenset({0})) + assert corrections.is_compatible([0, 1, 2]) + assert not corrections.is_compatible([2, 0, 1]) # Wrong order + assert not corrections.is_compatible([0, 1]) # Incomplete order + assert corrections.generate_total_measurement_order() in ([0, 1, 2], [0, 2, 1]) + assert nx.utils.graphs_equal(corrections.extract_dag(), nx.DiGraph([(0, 1), (0, 2)])) + + # Only output nodes + def test_from_measured_nodes_mapping_0(self) -> None: + og: OpenGraph[Plane] = OpenGraph( + graph=nx.Graph([(0, 1)]), + input_nodes=[], + output_nodes=[0, 1], + measurements={}, + ) + + corrections = XZCorrections.from_measured_nodes_mapping(og=og) + assert corrections.x_corrections == {} + assert corrections.z_corrections == {} + assert corrections.partial_order_layers == (frozenset({0, 1}),) + + # Empty corrections + def test_from_measured_nodes_mapping_1(self) -> None: + og: OpenGraph[Plane] = OpenGraph( + graph=nx.Graph([(0, 1)]), + input_nodes=[], + output_nodes=[1], + measurements={0: Plane.XY}, + ) + + corrections = XZCorrections.from_measured_nodes_mapping(og=og) + assert corrections.x_corrections == {} + assert corrections.z_corrections == {} + assert corrections.partial_order_layers == (frozenset({1}), frozenset({0})) + + def test_from_measured_nodes_mapping_2(self) -> None: + og = OpenGraph( + graph=nx.Graph([(0, 1), (2, 3)]), + input_nodes=[], + output_nodes=[3], + measurements=dict.fromkeys([0, 1, 2], Measurement(angle=0, plane=Plane.XY)), + ) + x_corrections = {0: {1}, 2: {3}} + + corrections = XZCorrections.from_measured_nodes_mapping(og=og, x_corrections=x_corrections) + + assert all(corrections.partial_order_layers) # No empty sets + assert all( + sum(1 for layer in corrections.partial_order_layers if node in layer) == 1 for node in og.graph.nodes + ) + + # All output nodes in corrections + def test_from_measured_nodes_mapping_3(self) -> None: + og = OpenGraph( + graph=nx.Graph([(0, 1), (2, 3)]), + input_nodes=[], + output_nodes=[1, 3], + measurements=dict.fromkeys([0, 2], Measurement(angle=0, plane=Plane.XY)), + ) + x_corrections = {0: {1}, 2: {3}} + + corrections = XZCorrections.from_measured_nodes_mapping(og=og, x_corrections=x_corrections) + + assert corrections.partial_order_layers == (frozenset({1, 3}), frozenset({0, 2})) + + # Some output nodes in corrections + def test_from_measured_nodes_mapping_4(self) -> None: + og = OpenGraph( + graph=nx.Graph([(0, 1), (2, 3)]), + input_nodes=[], + output_nodes=[1, 3], + measurements=dict.fromkeys([0, 2], Measurement(angle=0, plane=Plane.XY)), + ) + x_corrections = {2: {3}} + + corrections = XZCorrections.from_measured_nodes_mapping(og=og, x_corrections=x_corrections) + + assert corrections.partial_order_layers == (frozenset({1, 3}), frozenset({0, 2})) + + # No output nodes in corrections + def test_from_measured_nodes_mapping_5(self) -> None: + og = OpenGraph( + graph=nx.Graph([(0, 1), (2, 3)]), + input_nodes=[0], + output_nodes=[1, 3], + measurements=dict.fromkeys([0, 2], Measurement(angle=0, plane=Plane.XY)), + ) + x_corrections = {0: {2}} + + corrections = XZCorrections.from_measured_nodes_mapping(og=og, x_corrections=x_corrections) + + assert corrections.partial_order_layers == (frozenset({1, 3}), frozenset({2}), frozenset({0})) + + # Test exceptions + def test_from_measured_nodes_mapping_exceptions(self) -> None: + og = OpenGraph( + graph=nx.Graph([(0, 1), (1, 2), (2, 3)]), + input_nodes=[0], + output_nodes=[3], + measurements=dict.fromkeys(range(3), Measurement(angle=0, plane=Plane.XY)), + ) + with pytest.raises(ValueError, match=r"Keys of input X-corrections contain non-measured nodes."): + XZCorrections.from_measured_nodes_mapping(og=og, x_corrections={3: {1, 2}}) + + with pytest.raises(ValueError, match=r"Keys of input Z-corrections contain non-measured nodes."): + XZCorrections.from_measured_nodes_mapping(og=og, z_corrections={3: {1, 2}}) + + with pytest.raises( + ValueError, + match=r"Input XZ-corrections are not runnable since the induced directed graph contains closed loops.", + ): + XZCorrections.from_measured_nodes_mapping(og=og, x_corrections={0: {1}, 1: {2}}, z_corrections={2: {0}}) + + with pytest.raises( + ValueError, match=r"Values of input mapping contain labels which are not nodes of the input open graph." + ): + XZCorrections.from_measured_nodes_mapping(og=og, x_corrections={0: {4}}) diff --git a/tests/test_flow_find_gpflow.py b/tests/test_flow_find_gpflow.py new file mode 100644 index 00000000..704008db --- /dev/null +++ b/tests/test_flow_find_gpflow.py @@ -0,0 +1,280 @@ +"""Unit tests for the algebraic flow finding algorithm (for generalised or Pauli flow). + +This module tests the following: + - Computation of the reduced adjacency matrix. + - Computation of the flow-demand and order-demand matrices. We check this routine for open graphs with Pauli-angle measurements, intepreted as planes or as axes. + - Computation of the correction matrix. + - Computation of topological generations on small DAGs. + +The second part of the flow-finding algorithm (namely, verifying if the correction matrix is compatible with a DAG) is not done in this test module. For a complete test on the flow-finding algorithms see `tests.test_opengraph.py`. +""" + +from __future__ import annotations + +from typing import TYPE_CHECKING, NamedTuple + +import networkx as nx +import numpy as np +import pytest + +from graphix._linalg import MatGF2 +from graphix.flow._find_gpflow import ( + AlgebraicOpenGraph, + PlanarAlgebraicOpenGraph, + _compute_topological_generations, + compute_correction_matrix, +) +from graphix.fundamentals import Axis, Plane +from graphix.measurements import Measurement +from graphix.opengraph import OpenGraph + +if TYPE_CHECKING: + from graphix.fundamentals import AbstractMeasurement, AbstractPlanarMeasurement + + +class AlgebraicOpenGraphTestCase(NamedTuple): + aog: AlgebraicOpenGraph[AbstractMeasurement] | PlanarAlgebraicOpenGraph[AbstractPlanarMeasurement] + radj: MatGF2 + flow_demand_mat: MatGF2 + order_demand_mat: MatGF2 + has_corr_mat: bool + + +class DAGTestCase(NamedTuple): + adj_mat: MatGF2 + generations: list[list[int]] | None + + +def prepare_test_og() -> list[AlgebraicOpenGraphTestCase]: + test_cases: list[AlgebraicOpenGraphTestCase] = [] + + # Trivial open graph with pflow and nI = nO + def get_og_0() -> OpenGraph[Plane | Axis]: + """Return an open graph with Pauli flow and equal number of outputs and inputs. + + The returned graph has the following structure: + + [0]-1-(2) + """ + return OpenGraph( + graph=nx.Graph([(0, 1), (1, 2)]), input_nodes=[0], output_nodes=[2], measurements={0: Plane.XY, 1: Axis.Y} + ) + + test_cases.append( + AlgebraicOpenGraphTestCase( + aog=AlgebraicOpenGraph(get_og_0()), + radj=MatGF2([[1, 0], [0, 1]]), + flow_demand_mat=MatGF2([[1, 0], [1, 1]]), + order_demand_mat=MatGF2([[0, 0], [0, 0]]), + has_corr_mat=True, + ) + ) + + # Non-trivial open graph with pflow and nI = nO + def get_og_1() -> OpenGraph[Measurement]: + """Return an open graph with Pauli flow and equal number of outputs and inputs. + + The returned graph has the following structure: + + [0]-2-4-(6) + | | + [1]-3-5-(7) + """ + graph: nx.Graph[int] = nx.Graph([(0, 2), (1, 3), (2, 3), (2, 4), (3, 5), (4, 5), (4, 6), (5, 7)]) + inputs = [0, 1] + outputs = [6, 7] + meas = { + 0: Measurement(0.1, Plane.XY), # XY + 1: Measurement(0.1, Plane.XY), # XY + 2: Measurement(0.0, Plane.XY), # X + 3: Measurement(0.1, Plane.XY), # XY + 4: Measurement(0.0, Plane.XY), # X + 5: Measurement(0.5, Plane.XY), # Y + } + return OpenGraph(graph=graph, input_nodes=inputs, output_nodes=outputs, measurements=meas) + + test_cases.extend( + ( + AlgebraicOpenGraphTestCase( + aog=AlgebraicOpenGraph(get_og_1()), + radj=MatGF2( + [ + [1, 0, 0, 0, 0, 0], + [0, 1, 0, 0, 0, 0], + [0, 1, 1, 0, 0, 0], + [1, 0, 0, 1, 0, 0], + [1, 0, 0, 1, 1, 0], + [0, 1, 1, 0, 0, 1], + ] + ), + flow_demand_mat=MatGF2( + [ + [1, 0, 0, 0, 0, 0], + [0, 1, 0, 0, 0, 0], + [0, 1, 1, 0, 0, 0], + [1, 0, 0, 1, 0, 0], + [1, 0, 0, 1, 1, 0], + [0, 1, 1, 1, 0, 1], + ] + ), + order_demand_mat=MatGF2( + [ + [0, 0, 0, 0, 0, 0], + [0, 0, 0, 0, 0, 0], + [0, 0, 0, 0, 0, 0], + [0, 1, 0, 0, 0, 0], + [0, 0, 0, 0, 0, 0], + [0, 0, 0, 0, 0, 0], + ] + ), + has_corr_mat=True, + ), + # Same open graph but we interpret the measurements on Pauli axes as planar measurements, therefore, there flow-demand and order demand matrices are different. + AlgebraicOpenGraphTestCase( + aog=PlanarAlgebraicOpenGraph(get_og_1()), + radj=MatGF2( + [ + [1, 0, 0, 0, 0, 0], + [0, 1, 0, 0, 0, 0], + [0, 1, 1, 0, 0, 0], + [1, 0, 0, 1, 0, 0], + [1, 0, 0, 1, 1, 0], + [0, 1, 1, 0, 0, 1], + ] + ), + flow_demand_mat=MatGF2( + [ + [1, 0, 0, 0, 0, 0], + [0, 1, 0, 0, 0, 0], + [0, 1, 1, 0, 0, 0], + [1, 0, 0, 1, 0, 0], + [1, 0, 0, 1, 1, 0], + [0, 1, 1, 0, 0, 1], + ] + ), + order_demand_mat=MatGF2( + [ + [0, 0, 0, 0, 0, 0], + [0, 0, 0, 0, 0, 0], + [1, 0, 0, 0, 0, 0], + [0, 1, 0, 0, 0, 0], + [0, 0, 1, 0, 0, 0], + [0, 0, 0, 1, 0, 0], + ] + ), + has_corr_mat=True, + ), + ) + ) + + # Non-trivial open graph with pflow and nI != nO + def get_og_2() -> OpenGraph[Measurement]: + """Return an open graph with Pauli flow and unequal number of outputs and inputs. + + Example from Fig. 1 in Mitosek and Backens, 2024 (arXiv:2410.23439). + """ + graph: nx.Graph[int] = nx.Graph( + [(0, 2), (2, 4), (3, 4), (4, 6), (1, 4), (1, 6), (2, 3), (3, 5), (2, 6), (3, 6)] + ) + inputs = [0] + outputs = [5, 6] + meas = { + 0: Measurement(0.1, Plane.XY), # XY + 1: Measurement(0.1, Plane.XZ), # XZ + 2: Measurement(0.5, Plane.YZ), # Y + 3: Measurement(0.1, Plane.XY), # XY + 4: Measurement(0, Plane.XZ), # Z + } + + return OpenGraph(graph=graph, input_nodes=inputs, output_nodes=outputs, measurements=meas) + + test_cases.extend( + ( + AlgebraicOpenGraphTestCase( + aog=AlgebraicOpenGraph(get_og_2()), + radj=MatGF2( + [[0, 1, 0, 0, 0, 0], [0, 0, 0, 1, 0, 1], [0, 0, 1, 1, 0, 1], [0, 1, 0, 1, 1, 1], [1, 1, 1, 0, 0, 1]] + ), + flow_demand_mat=MatGF2( + [[0, 1, 0, 0, 0, 0], [1, 0, 0, 0, 0, 0], [0, 1, 1, 1, 0, 1], [0, 1, 0, 1, 1, 1], [0, 0, 0, 1, 0, 0]] + ), + order_demand_mat=MatGF2( + [[0, 0, 0, 0, 0, 0], [1, 0, 0, 1, 0, 1], [0, 0, 0, 0, 0, 0], [0, 0, 1, 0, 0, 0], [0, 0, 0, 0, 0, 0]] + ), + has_corr_mat=True, + ), + # Same open graph but we interpret the measurements on Pauli axes as planar measurements, therefore, there flow-demand and order demand matrices are different. + # The new flow-demand matrix is not invertible, therefore the open graph does not have gflow. + AlgebraicOpenGraphTestCase( + aog=PlanarAlgebraicOpenGraph(get_og_2()), + radj=MatGF2( + [[0, 1, 0, 0, 0, 0], [0, 0, 0, 1, 0, 1], [0, 0, 1, 1, 0, 1], [0, 1, 0, 1, 1, 1], [1, 1, 1, 0, 0, 1]] + ), + flow_demand_mat=MatGF2( + [[0, 1, 0, 0, 0, 0], [1, 0, 0, 0, 0, 0], [0, 1, 0, 0, 0, 0], [0, 1, 0, 1, 1, 1], [0, 0, 0, 1, 0, 0]] + ), + order_demand_mat=MatGF2( + [[0, 0, 0, 0, 0, 0], [1, 0, 0, 1, 0, 1], [0, 0, 1, 1, 0, 1], [0, 0, 1, 0, 0, 0], [1, 1, 1, 1, 0, 1]] + ), + has_corr_mat=False, + ), + ) + ) + return test_cases + + +def prepare_test_dag() -> list[DAGTestCase]: + test_cases: list[DAGTestCase] = [] + + # Simple DAG + test_cases.extend( + ( # Simple DAG + DAGTestCase( + adj_mat=MatGF2([[0, 0, 0, 0], [1, 0, 0, 0], [1, 0, 0, 0], [0, 1, 1, 0]]), generations=[[0], [1, 2], [3]] + ), + # Graph with loop + DAGTestCase(adj_mat=MatGF2([[0, 0, 0, 0], [1, 0, 0, 1], [1, 0, 0, 0], [0, 1, 1, 0]]), generations=None), + # Disconnected graph + DAGTestCase( + adj_mat=MatGF2([[0, 0, 0, 0, 0], [1, 0, 0, 0, 0], [0, 0, 0, 0, 0], [0, 0, 1, 0, 0], [0, 0, 1, 0, 0]]), + generations=[[0, 2], [1, 3, 4]], + ), + ) + ) + + return test_cases + + +class TestAlgebraicFlow: + @pytest.mark.parametrize("test_case", prepare_test_og()) + def test_compute_reduced_adj(self, test_case: AlgebraicOpenGraphTestCase) -> None: + aog = test_case.aog + radj = aog._compute_reduced_adj() + assert np.all(radj == test_case.radj) + + @pytest.mark.parametrize("test_case", prepare_test_og()) + def test_og_matrices(self, test_case: AlgebraicOpenGraphTestCase) -> None: + aog = test_case.aog + assert np.all(aog.flow_demand_matrix == test_case.flow_demand_mat) + assert np.all(aog.order_demand_matrix == test_case.order_demand_mat) + + @pytest.mark.parametrize("test_case", prepare_test_og()) + def test_correction_matrix(self, test_case: AlgebraicOpenGraphTestCase) -> None: + aog = test_case.aog + corr_matrix = compute_correction_matrix(aog) + + ident = MatGF2(np.eye(len(aog.non_outputs), dtype=np.uint8)) + if test_case.has_corr_mat: + assert corr_matrix is not None + assert np.all( + (test_case.flow_demand_mat @ corr_matrix.c_matrix) % 2 == ident + ) # Test with numpy matrix product. + else: + assert corr_matrix is None + + @pytest.mark.parametrize("test_case", prepare_test_dag()) + def test_compute_topological_generations(self, test_case: DAGTestCase) -> None: + adj_mat = test_case.adj_mat + generations_ref = test_case.generations + + assert generations_ref == _compute_topological_generations(adj_mat) diff --git a/tests/test_generator.py b/tests/test_generator.py deleted file mode 100644 index c9f8989d..00000000 --- a/tests/test_generator.py +++ /dev/null @@ -1,156 +0,0 @@ -from __future__ import annotations - -from typing import TYPE_CHECKING - -import networkx as nx -import numpy as np -import pytest - -from graphix.fundamentals import Plane -from graphix.generator import generate_from_graph -from graphix.gflow import find_gflow, find_pauliflow, pauliflow_from_pattern -from graphix.measurements import Measurement -from graphix.opengraph import OpenGraph -from graphix.random_objects import rand_gate - -if TYPE_CHECKING: - from collections.abc import Callable - - from numpy.random import Generator - - from graphix import Pattern - - -def example_flow(rng: Generator) -> Pattern: - graph: nx.Graph[int] = nx.Graph([(0, 3), (1, 4), (2, 5), (1, 3), (2, 4), (3, 6), (4, 7), (5, 8)]) - inputs = [1, 0, 2] # non-trivial order to check order is conserved. - outputs = [7, 6, 8] - angles = dict(zip(range(6), (2 * rng.random(6)).tolist())) - meas_planes = dict.fromkeys(range(6), Plane.XY) - - pattern = generate_from_graph(graph, angles, inputs, outputs, meas_planes=meas_planes) - pattern.standardize() - - assert pattern.input_nodes == inputs - assert pattern.output_nodes == outputs - return pattern - - -def example_gflow(rng: Generator) -> Pattern: - graph: nx.Graph[int] = nx.Graph([(1, 2), (2, 3), (3, 4), (4, 5), (5, 6), (3, 6), (1, 6)]) - inputs = [3, 1, 5] - outputs = [4, 2, 6] - angles = dict(zip([1, 3, 5], (2 * rng.random(3)).tolist())) - meas_planes = dict.fromkeys([1, 3, 5], Plane.XY) - - pattern = generate_from_graph(graph, angles, inputs, outputs, meas_planes=meas_planes) - pattern.standardize() - - assert pattern.input_nodes == inputs - assert pattern.output_nodes == outputs - return pattern - - -def example_graph_pflow(rng: Generator) -> OpenGraph: - """Create a graph which has pflow but no gflow. - - Parameters - ---------- - rng : :class:`numpy.random.Generator` - See graphix.tests.conftest.py - - Returns - ------- - OpenGraph: :class:`graphix.opengraph.OpenGraph` - """ - graph: nx.Graph[int] = nx.Graph( - [(0, 2), (1, 4), (2, 3), (3, 4), (2, 5), (3, 6), (4, 7), (5, 6), (6, 7), (5, 8), (7, 9)] - ) - inputs = [1, 0] - outputs = [9, 8] - - # Heuristic mixture of Pauli and non-Pauli angles ensuring there's no gflow but there's pflow. - meas_angles: dict[int, float] = { - **dict.fromkeys(range(4), 0), - **dict(zip(range(4, 8), (2 * rng.random(4)).tolist())), - } - meas_planes = dict.fromkeys(range(8), Plane.XY) - meas = {i: Measurement(angle, plane) for (i, angle), plane in zip(meas_angles.items(), meas_planes.values())} - - gf, _ = find_gflow(graph=graph, iset=set(inputs), oset=set(outputs), meas_planes=meas_planes) - pf, _ = find_pauliflow( - graph=graph, iset=set(inputs), oset=set(outputs), meas_planes=meas_planes, meas_angles=meas_angles - ) - - assert gf is None # example graph doesn't have gflow - assert pf is not None # example graph has Pauli flow - - return OpenGraph(inside=graph, inputs=inputs, outputs=outputs, measurements=meas) - - -def example_pflow(rng: Generator) -> Pattern: - og = example_graph_pflow(rng) - pattern = og.to_pattern() - pattern.standardize() - assert og.inputs == pattern.input_nodes - assert og.outputs == pattern.output_nodes - return pattern - - -class TestGenerator: - @pytest.mark.parametrize("example", [example_flow, example_gflow, example_pflow]) - def test_pattern_generation_determinism(self, example: Callable[[Generator], Pattern], fx_rng: Generator) -> None: - pattern = example(fx_rng) - pattern.minimize_space() - - repeats = 3 # for testing the determinism of a pattern - results = [pattern.simulate_pattern(rng=fx_rng) for _ in range(repeats)] - - for i in range(1, 3): - inner_product = np.dot(results[0].flatten(), results[i].flatten().conjugate()) - assert abs(inner_product) == pytest.approx(1) - - def test_pattern_generation_flow(self, fx_rng: Generator) -> None: - nqubits = 3 - depth = 2 - pairs = [(0, 1), (1, 2)] - circuit = rand_gate(nqubits, depth, pairs, fx_rng) - # transpile into graph - pattern = circuit.transpile().pattern - pattern.standardize() - pattern.shift_signals() - # get the graph and generate pattern again with flow algorithm - graph = pattern.extract_graph() - input_list = [0, 1, 2] - angles: dict[int, float] = {} - for cmd in pattern.extract_measurement_commands(): - assert isinstance(cmd.angle, float) - angles[cmd.node] = float(cmd.angle) - meas_planes = pattern.get_meas_plane() - pattern2 = generate_from_graph(graph, angles, input_list, pattern.output_nodes, meas_planes) - # check that the new one runs and returns correct result - pattern2.standardize() - pattern2.shift_signals() - pattern2.minimize_space() - state = circuit.simulate_statevector().statevec - state_mbqc = pattern2.simulate_pattern(rng=fx_rng) - assert np.abs(np.dot(state_mbqc.flatten().conjugate(), state.flatten())) == pytest.approx(1) - - def test_pattern_generation_no_internal_nodes(self) -> None: - g: nx.Graph[int] = nx.Graph() - g.add_edges_from([(0, 1), (1, 2)]) - pattern = generate_from_graph(g, {}, {0, 1, 2}, {0, 1, 2}, {}) - graph = pattern.extract_graph() - graph_ref = nx.Graph(((0, 1), (1, 2))) - assert nx.utils.graphs_equal(graph, graph_ref) - - def test_pattern_generation_pflow(self, fx_rng: Generator) -> None: - og = example_graph_pflow(fx_rng) - pattern = og.to_pattern() - - graph_generated_pattern = pattern.extract_graph() - assert nx.is_isomorphic(og.inside, graph_generated_pattern) - - pattern.standardize() - pf_generated_pattern, _ = pauliflow_from_pattern(pattern) - assert pf_generated_pattern is not None diff --git a/tests/test_opengraph.py b/tests/test_opengraph.py index 58844fd5..114f1a3a 100644 --- a/tests/test_opengraph.py +++ b/tests/test_opengraph.py @@ -1,50 +1,641 @@ +"""Unit tests for the class `:class: graphix.opengraph.OpenGraph`. + +This module tests the full conversion Open Graph -> Flow -> XZ-corrections -> Pattern for all three classes of flow. +Output correctness is verified by checking if the resulting pattern is deterministic (when the flow exists). +""" + from __future__ import annotations +from typing import TYPE_CHECKING, NamedTuple + import networkx as nx import numpy as np import pytest -from graphix import Pattern, command +from graphix.command import E from graphix.fundamentals import Plane from graphix.measurements import Measurement -from graphix.opengraph import OpenGraph +from graphix.opengraph import OpenGraph, OpenGraphError +from graphix.pattern import Pattern +from graphix.random_objects import rand_circuit +from graphix.states import PlanarState +if TYPE_CHECKING: + from typing import Callable -# Tests whether an open graph can be converted to and from a pattern and be -# successfully reconstructed. -def test_open_graph_to_pattern() -> None: - g: nx.Graph[int] - g = nx.Graph([(0, 1), (1, 2)]) - inputs = [0] - outputs = [2] - meas = {0: Measurement(0, Plane.XY), 1: Measurement(0, Plane.XY)} - og = OpenGraph(g, meas, inputs, outputs) + from numpy.random import Generator - pattern = og.to_pattern() - og_reconstructed = OpenGraph.from_pattern(pattern) - assert og.isclose(og_reconstructed) +class OpenGraphFlowTestCase(NamedTuple): + og: OpenGraph[Measurement] + has_cflow: bool + has_gflow: bool + has_pflow: bool + + +OPEN_GRAPH_FLOW_TEST_CASES: list[OpenGraphFlowTestCase] = [] + + +def register_open_graph_flow_test_case( + test_case: Callable[[], OpenGraphFlowTestCase], +) -> Callable[[], OpenGraphFlowTestCase]: + OPEN_GRAPH_FLOW_TEST_CASES.append(test_case()) + return test_case + + +@register_open_graph_flow_test_case +def _og_0() -> OpenGraphFlowTestCase: + """Generate open graph. + + Structure: + + [(0)]-[(1)] + """ + meas: dict[int, Measurement] = {} + og = OpenGraph( + graph=nx.Graph([(0, 1)]), + input_nodes=[0, 1], + output_nodes=[0, 1], + measurements=meas, + ) + return OpenGraphFlowTestCase(og, has_cflow=True, has_gflow=True, has_pflow=True) + + +@register_open_graph_flow_test_case +def _og_1() -> OpenGraphFlowTestCase: + """Generate open graph. + + Structure: + + [0]-1-20-30-4-(5) + """ + og = OpenGraph( + graph=nx.Graph([(0, 1), (1, 20), (20, 30), (30, 4), (4, 5)]), + input_nodes=[0], + output_nodes=[5], + measurements={ + 0: Measurement(0.1, Plane.XY), + 1: Measurement(0.2, Plane.XY), + 20: Measurement(0.3, Plane.XY), + 30: Measurement(0.4, Plane.XY), + 4: Measurement(0.5, Plane.XY), + }, + ) + return OpenGraphFlowTestCase(og, has_cflow=True, has_gflow=True, has_pflow=True) + + +@register_open_graph_flow_test_case +def _og_2() -> OpenGraphFlowTestCase: + """Generate open graph. + + Structure: + + [1]-3-(5) + | + [2]-4-(6) + """ + og = OpenGraph( + graph=nx.Graph([(1, 3), (2, 4), (3, 4), (3, 5), (4, 6)]), + input_nodes=[1, 2], + output_nodes=[5, 6], + measurements={ + 1: Measurement(0.1, Plane.XY), + 2: Measurement(0.2, Plane.XY), + 3: Measurement(0.3, Plane.XY), + 4: Measurement(0.4, Plane.XY), + }, + ) + return OpenGraphFlowTestCase(og, has_cflow=True, has_gflow=True, has_pflow=True) + + +@register_open_graph_flow_test_case +def _og_3() -> OpenGraphFlowTestCase: + r"""Generate open graph. + + Structure: + + [1]-(4) + \ / + X____ + / | + [2]-(5) | + \ / | + X | + / \ / + [3]-(6) + """ + og = OpenGraph( + graph=nx.Graph([(1, 4), (1, 6), (2, 4), (2, 5), (2, 6), (3, 5), (3, 6)]), + input_nodes=[1, 2, 3], + output_nodes=[4, 5, 6], + measurements={ + 1: Measurement(0.1, Plane.XY), + 2: Measurement(0.2, Plane.XY), + 3: Measurement(0.3, Plane.XY), + }, + ) + return OpenGraphFlowTestCase(og, has_cflow=False, has_gflow=True, has_pflow=True) + + +@register_open_graph_flow_test_case +def _og_4() -> OpenGraphFlowTestCase: + r"""Generate open graph. + + Structure: + + [0]-[1] + /| | + (4)| | + \| | + 2--(5)-3 + """ + og = OpenGraph( + graph=nx.Graph([(0, 1), (0, 2), (0, 4), (1, 5), (2, 4), (2, 5), (3, 5)]), + input_nodes=[0, 1], + output_nodes=[4, 5], + measurements={ + 0: Measurement(0.1, Plane.XY), + 1: Measurement(0.1, Plane.XY), + 2: Measurement(0.2, Plane.XZ), + 3: Measurement(0.3, Plane.YZ), + }, + ) + return OpenGraphFlowTestCase(og, has_cflow=False, has_gflow=True, has_pflow=True) + + +@register_open_graph_flow_test_case +def _og_5() -> OpenGraphFlowTestCase: + r"""Generate open graph. + + Structure: + + [1]-(3) + \ / + X + / \ + [2]-(4) + """ + og = OpenGraph( + graph=nx.Graph([(1, 3), (1, 4), (2, 3), (2, 4)]), + input_nodes=[1, 2], + output_nodes=[3, 4], + measurements={ + 1: Measurement(0.1, Plane.XY), + 2: Measurement(0.1, Plane.XY), + }, + ) + return OpenGraphFlowTestCase(og, has_cflow=False, has_gflow=False, has_pflow=False) + + +@register_open_graph_flow_test_case +def _og_6() -> OpenGraphFlowTestCase: + r"""Generate open graph. + + Structure: + + [0] + | + 1-2-3 + | + (4) + """ + og = OpenGraph( + graph=nx.Graph([(0, 1), (1, 2), (1, 4), (2, 3)]), + input_nodes=[0], + output_nodes=[4], + measurements={ + 0: Measurement(0.1, Plane.XY), # XY + 1: Measurement(0, Plane.XY), # X + 2: Measurement(0.1, Plane.XY), # XY + 3: Measurement(0, Plane.XY), # X + }, + ) + return OpenGraphFlowTestCase(og, has_cflow=False, has_gflow=False, has_pflow=True) + + +@register_open_graph_flow_test_case +def _og_7() -> OpenGraphFlowTestCase: + r"""Generate open graph. + + Structure: + + [0]-1-(2) + """ + og = OpenGraph( + graph=nx.Graph([(0, 1), (1, 2)]), + input_nodes=[0], + output_nodes=[2], + measurements={ + 0: Measurement(0.1, Plane.XY), # XY + 1: Measurement(0.5, Plane.YZ), # Y + }, + ) + return OpenGraphFlowTestCase(og, has_cflow=False, has_gflow=False, has_pflow=True) + + +@register_open_graph_flow_test_case +def _og_8() -> OpenGraphFlowTestCase: + r"""Generate open graph. + + Structure: + + [0]-2-4-(6) + | | + [1]-3-5-(7) + """ + og = OpenGraph( + graph=nx.Graph([(0, 2), (1, 3), (2, 3), (2, 4), (3, 5), (4, 5), (4, 6), (5, 7)]), + input_nodes=[1, 0], + output_nodes=[6, 7], + measurements={ + 0: Measurement(0.1, Plane.XY), # XY + 1: Measurement(0.1, Plane.XZ), # XZ + 2: Measurement(0.5, Plane.XZ), # X + 3: Measurement(0.5, Plane.YZ), # Y + 4: Measurement(0.5, Plane.YZ), # Y + 5: Measurement(0.1, Plane.YZ), # YZ + }, + ) + return OpenGraphFlowTestCase(og, has_cflow=False, has_gflow=False, has_pflow=False) + + +@register_open_graph_flow_test_case +def _og_9() -> OpenGraphFlowTestCase: + r"""Generate open graph. + + Structure: + + [0]-2-4-(6) + | | + [1]-3-5-(7) + """ + og = OpenGraph( + graph=nx.Graph([(0, 2), (1, 3), (2, 3), (2, 4), (3, 5), (4, 5), (4, 6), (5, 7)]), + input_nodes=[0, 1], + output_nodes=[6, 7], + measurements={ + 0: Measurement(0.1, Plane.XY), # XY + 1: Measurement(0.1, Plane.XY), # XY + 2: Measurement(0.0, Plane.XY), # X + 3: Measurement(0.1, Plane.XY), # XY + 4: Measurement(0.0, Plane.XY), # X + 5: Measurement(0.5, Plane.XY), # Y + }, + ) + return OpenGraphFlowTestCase(og, has_cflow=True, has_gflow=True, has_pflow=True) + + +@register_open_graph_flow_test_case +def _og_10() -> OpenGraphFlowTestCase: + r"""Generate open graph. + + Structure: + + 3-(5) + /| \ + [0]-2-4-(6) + | | /| + | 1 | + |____| + + Notes + ----- + Example from Fig. 1 in Mitosek and Backens, 2024 (arXiv:2410.23439). + """ + og = OpenGraph( + graph=nx.Graph([(0, 2), (2, 4), (3, 4), (4, 6), (1, 4), (1, 6), (2, 3), (3, 5), (2, 6), (3, 6)]), + input_nodes=[0], + output_nodes=[5, 6], + measurements={ + 0: Measurement(0.1, Plane.XY), # XY + 1: Measurement(0.1, Plane.XZ), # XZ + 2: Measurement(0.5, Plane.YZ), # Y + 3: Measurement(0.1, Plane.XY), # XY + 4: Measurement(0, Plane.XZ), # Z + }, + ) + return OpenGraphFlowTestCase(og, has_cflow=False, has_gflow=False, has_pflow=True) + + +@register_open_graph_flow_test_case +def _og_11() -> OpenGraphFlowTestCase: + r"""Generate open graph. + + Structure: + + [0]-2--3--4-7-(8) + | | | + (6)[1](5) + """ + og = OpenGraph( + graph=nx.Graph([(0, 2), (1, 3), (2, 3), (2, 6), (3, 4), (4, 7), (4, 5), (7, 8)]), + input_nodes=[0, 1], + output_nodes=[5, 6, 8], + measurements={ + 0: Measurement(0.1, Plane.XY), + 1: Measurement(0.1, Plane.XY), + 2: Measurement(0.0, Plane.XY), + 3: Measurement(0, Plane.XY), + 4: Measurement(0.5, Plane.XY), + 7: Measurement(0, Plane.XY), + }, + ) + return OpenGraphFlowTestCase(og, has_cflow=True, has_gflow=True, has_pflow=True) + + +@register_open_graph_flow_test_case +def _og_12() -> OpenGraphFlowTestCase: + r"""Generate open graph. + + Structure: + + [0]-2-(3)-(4) + | + [(1)] + """ + og = OpenGraph( + graph=nx.Graph([(0, 2), (1, 2), (2, 3), (3, 4)]), + input_nodes=[0, 1], + output_nodes=[1, 3, 4], + measurements={0: Measurement(0.1, Plane.XY), 2: Measurement(0.5, Plane.YZ)}, + ) + return OpenGraphFlowTestCase(og, has_cflow=False, has_gflow=False, has_pflow=True) + + +@register_open_graph_flow_test_case +def _og_13() -> OpenGraphFlowTestCase: + r"""Generate open graph. + + Structure: + + 0-[1] + | | + 5-(2)-3--4-(7) + | + (6) + """ + og = OpenGraph( + graph=nx.Graph([(0, 1), (0, 3), (1, 4), (3, 4), (2, 3), (2, 5), (3, 6), (4, 7)]), + input_nodes=[1], + output_nodes=[6, 2, 7], + measurements={ + 0: Measurement(0.1, Plane.XZ), + 1: Measurement(0.1, Plane.XY), + 3: Measurement(0, Plane.XY), + 4: Measurement(0.1, Plane.XY), + 5: Measurement(0.1, Plane.YZ), + }, + ) + return OpenGraphFlowTestCase(og, has_cflow=False, has_gflow=True, has_pflow=True) - # 0 -- 1 -- 2 - # | - # 3 -- 4 -- 5 - g = nx.Graph([(0, 1), (1, 2), (1, 4), (3, 4), (4, 5)]) - inputs = [0, 3] - outputs = [2, 5] - meas = { - 0: Measurement(0, Plane.XY), - 1: Measurement(1.0, Plane.XY), - 3: Measurement(0.5, Plane.YZ), - 4: Measurement(1.0, Plane.XY), - } - og = OpenGraph(g, meas, inputs, outputs) +@register_open_graph_flow_test_case +def _og_14() -> OpenGraphFlowTestCase: + r"""Generate open graph. - pattern = og.to_pattern() - og_reconstructed = OpenGraph.from_pattern(pattern) + Structure: + + 0-(1) (4)-6 + | | + 2-(3) + """ + og = OpenGraph( + graph=nx.Graph([(0, 1), (0, 2), (2, 3), (1, 3), (4, 6)]), + input_nodes=[], + output_nodes=[1, 3, 4], + measurements={0: Measurement(0.5, Plane.XZ), 2: Measurement(0, Plane.YZ), 6: Measurement(0.2, Plane.XY)}, + ) + return OpenGraphFlowTestCase(og, has_cflow=False, has_gflow=True, has_pflow=True) + + +@register_open_graph_flow_test_case +def _og_15() -> OpenGraphFlowTestCase: + r"""Generate open graph. + + Structure: + + [1]--0 + | | + 4---3-(2) + | | | + (7)-(6)-5 + """ + og = OpenGraph( + graph=nx.Graph([(0, 1), (0, 3), (1, 4), (3, 4), (2, 3), (2, 5), (3, 6), (4, 7), (5, 6), (6, 7)]), + input_nodes=[1], + output_nodes=[6, 2, 7], + measurements={ + 0: Measurement(0.1, Plane.XZ), + 1: Measurement(0.1, Plane.XY), + 3: Measurement(0, Plane.XY), + 4: Measurement(0.1, Plane.XY), + 5: Measurement(0.1, Plane.XY), + }, + ) + return OpenGraphFlowTestCase(og, has_cflow=False, has_gflow=False, has_pflow=False) - assert og.isclose(og_reconstructed) +@register_open_graph_flow_test_case +def _og_16() -> OpenGraphFlowTestCase: + r"""Generate open graph. + + Structure: + + [0]-(1) (4)-6 + | | + 2--(3) + """ + og = OpenGraph( + graph=nx.Graph([(0, 1), (0, 2), (2, 3), (1, 3), (4, 6)]), + input_nodes=[0], + output_nodes=[1, 3, 4], + measurements={0: Measurement(0.1, Plane.XZ), 2: Measurement(0, Plane.YZ), 6: Measurement(0.2, Plane.XY)}, + ) + return OpenGraphFlowTestCase(og, has_cflow=False, has_gflow=False, has_pflow=False) + + +@register_open_graph_flow_test_case +def _og_17() -> OpenGraphFlowTestCase: + r"""Generate open graph. + + Structure: + + 0-[1] [8] + | | + 5-(2)-3--4-(7) + | + (6) + + Notes + ----- + Graph is constructed by adding a disconnected input to OG 13. + """ + graph: nx.Graph[int] = nx.Graph([(0, 1), (0, 3), (1, 4), (3, 4), (2, 3), (2, 5), (3, 6), (4, 7)]) + graph.add_node(8) + og = OpenGraph( + graph=graph, + input_nodes=[1, 8], + output_nodes=[6, 2, 7], + measurements={ + 0: Measurement(0.1, Plane.XZ), + 1: Measurement(0.1, Plane.XY), + 3: Measurement(0, Plane.XY), + 4: Measurement(0.1, Plane.XY), + 5: Measurement(0.1, Plane.YZ), + 8: Measurement(0.1, Plane.XY), + }, + ) + return OpenGraphFlowTestCase(og, has_cflow=False, has_gflow=False, has_pflow=False) + + +@register_open_graph_flow_test_case +def _og_18() -> OpenGraphFlowTestCase: + r"""Generate open graph. + + Structure: + + [0]-2--3--4-7-(8) + | | | + (6)[1](5) + """ + og = OpenGraph( + graph=nx.Graph([(0, 2), (1, 3), (2, 3), (2, 6), (3, 4), (4, 7), (4, 5), (7, 8)]), + input_nodes=[0, 1], + output_nodes=[5, 6, 8], + measurements={ + 0: Measurement(0, Plane.XY), + 1: Measurement(0, Plane.XY), + 2: Measurement(0, Plane.XZ), + 3: Measurement(0, Plane.XY), + 4: Measurement(0.5, Plane.XY), + 7: Measurement(0, Plane.YZ), + }, + ) + return OpenGraphFlowTestCase(og, has_cflow=False, has_gflow=False, has_pflow=False) + + +@register_open_graph_flow_test_case +def _og_19() -> OpenGraphFlowTestCase: + r"""Generate open graph. + + Structure: + + 0-2--3--4-7-(8) + | | | + (6) 1 (5) + + Notes + ----- + Even though any node is measured in the XY plane, OG has Pauli flow because none of them are inputs. + """ + og = OpenGraph( + graph=nx.Graph([(0, 2), (1, 3), (2, 3), (2, 6), (3, 4), (4, 7), (4, 5), (7, 8)]), + input_nodes=[], + output_nodes=[5, 6, 8], + measurements={ + 0: Measurement(0, Plane.XZ), + 1: Measurement(0, Plane.XZ), + 2: Measurement(0, Plane.XZ), + 3: Measurement(0, Plane.XZ), + 4: Measurement(0, Plane.XZ), + 7: Measurement(0, Plane.XZ), + }, + ) + return OpenGraphFlowTestCase(og, has_cflow=False, has_gflow=False, has_pflow=True) + + +def check_determinism(pattern: Pattern, fx_rng: Generator, n_shots: int = 3) -> bool: + """Verify if the input pattern is deterministic.""" + for plane in {Plane.XY, Plane.XZ, Plane.YZ}: + alpha = 2 * np.pi * fx_rng.random() + state_ref = pattern.simulate_pattern(input_state=PlanarState(plane, alpha)) + + for _ in range(n_shots): + state = pattern.simulate_pattern(input_state=PlanarState(plane, alpha)) + result = np.abs(np.dot(state.flatten().conjugate(), state_ref.flatten())) + + if result: + continue + return False + + return True + + +class TestOpenGraph: + def test_odd_neighbors(self) -> None: + graph: nx.Graph[int] = nx.Graph([(0, 1), (0, 2), (1, 3), (1, 2), (2, 3), (1, 4)]) + og = OpenGraph(graph=graph, input_nodes=[0], output_nodes=[1, 2, 3, 4], measurements={0: Plane.XY}) + + assert og.odd_neighbors([0]) == {1, 2} + assert og.odd_neighbors([0, 1]) == {0, 1, 3, 4} + assert og.odd_neighbors([1, 2, 3]) == {4} + assert og.odd_neighbors([]) == set() + + def test_neighbors(self) -> None: + graph: nx.Graph[int] = nx.Graph([(0, 1), (0, 2), (1, 3), (1, 2), (2, 3), (1, 4)]) + og = OpenGraph(graph=graph, input_nodes=[0], output_nodes=[1, 2, 3, 4], measurements={0: Plane.XY}) + + assert og.neighbors([4]) == {1} + assert og.neighbors([0, 1]) == {0, 1, 2, 3, 4} + assert og.neighbors([1, 2, 3]) == {0, 1, 2, 3, 4} + assert og.neighbors([]) == set() + + @pytest.mark.parametrize("test_case", OPEN_GRAPH_FLOW_TEST_CASES) + def test_cflow(self, test_case: OpenGraphFlowTestCase, fx_rng: Generator) -> None: + og = test_case.og + + if test_case.has_cflow: + pattern = og.extract_causal_flow().to_corrections().to_pattern() + assert check_determinism(pattern, fx_rng) + else: + with pytest.raises(OpenGraphError, match=r"The open graph does not have a causal flow."): + og.extract_causal_flow() + + @pytest.mark.parametrize("test_case", OPEN_GRAPH_FLOW_TEST_CASES) + def test_gflow(self, test_case: OpenGraphFlowTestCase, fx_rng: Generator) -> None: + og = test_case.og + + if test_case.has_gflow: + pattern = og.extract_gflow().to_corrections().to_pattern() + assert check_determinism(pattern, fx_rng) + else: + with pytest.raises(OpenGraphError, match=r"The open graph does not have a gflow."): + og.extract_gflow() + + @pytest.mark.parametrize("test_case", OPEN_GRAPH_FLOW_TEST_CASES) + def test_pflow(self, test_case: OpenGraphFlowTestCase, fx_rng: Generator) -> None: + og = test_case.og + + if test_case.has_pflow: + pattern = og.extract_pauli_flow().to_corrections().to_pattern() + assert check_determinism(pattern, fx_rng) + else: + with pytest.raises(OpenGraphError, match=r"The open graph does not have a Pauli flow."): + og.extract_pauli_flow() + + def test_double_entanglement(self) -> None: + pattern = Pattern(input_nodes=[0, 1], cmds=[E((0, 1)), E((0, 1))]) + pattern2 = pattern.extract_opengraph().to_pattern() + state = pattern.simulate_pattern() + state2 = pattern2.simulate_pattern() + assert np.abs(np.dot(state.flatten().conjugate(), state2.flatten())) == pytest.approx(1) + + def test_from_to_pattern(self, fx_rng: Generator) -> None: + n_qubits = 2 + depth = 2 + circuit = rand_circuit(n_qubits, depth, fx_rng) + pattern_ref = circuit.transpile().pattern + pattern = pattern_ref.extract_opengraph().to_pattern() + + for plane in {Plane.XY, Plane.XZ, Plane.YZ}: + alpha = 2 * np.pi * fx_rng.random() + state_ref = pattern_ref.simulate_pattern(input_state=PlanarState(plane, alpha)) + state = pattern.simulate_pattern(input_state=PlanarState(plane, alpha)) + assert np.abs(np.dot(state.flatten().conjugate(), state_ref.flatten())) == pytest.approx(1) + + +# TODO: Add test `OpenGraph.is_close` +# TODO: rewrite as parametric tests # Tests composition of two graphs @@ -68,7 +659,7 @@ def test_compose_1() -> None: inputs = [1] outputs = [2] meas = {i: Measurement(0, Plane.XY) for i in g.nodes - set(outputs)} - og_1 = OpenGraph(g, meas, inputs, outputs) + og_1 = OpenGraph(g, inputs, outputs, meas) mapping = {1: 100, 2: 200} @@ -76,11 +667,11 @@ def test_compose_1() -> None: expected_graph: nx.Graph[int] expected_graph = nx.Graph([(1, 2), (100, 200)]) - assert nx.is_isomorphic(og.inside, expected_graph) - assert og.inputs == [1, 100] - assert og.outputs == [2, 200] + assert nx.is_isomorphic(og.graph, expected_graph) + assert og.input_nodes == [1, 100] + assert og.output_nodes == [2, 200] - outputs_c = {i for i in og.inside.nodes if i not in og.outputs} + outputs_c = {i for i in og.graph.nodes if i not in og.output_nodes} assert og.measurements.keys() == outputs_c assert mapping.keys() <= mapping_complete.keys() assert set(mapping.values()) <= set(mapping_complete.values()) @@ -110,13 +701,13 @@ def test_compose_2() -> None: inputs = [0, 3] outputs = [13, 23] meas = {i: Measurement(0, Plane.XY) for i in g.nodes - set(outputs)} - og_1 = OpenGraph(g, meas, inputs, outputs) + og_1 = OpenGraph(g, inputs, outputs, meas) g = nx.Graph([(6, 7), (6, 17), (17, 1), (7, 4), (17, 4), (4, 2)]) inputs = [6, 7] outputs = [1, 2] meas = {i: Measurement(0, Plane.XY) for i in g.nodes - set(outputs)} - og_2 = OpenGraph(g, meas, inputs, outputs) + og_2 = OpenGraph(g, inputs, outputs, meas) mapping = {6: 23, 7: 13, 1: 100, 2: 200} @@ -126,11 +717,11 @@ def test_compose_2() -> None: expected_graph = nx.Graph( [(0, 17), (17, 23), (17, 4), (3, 4), (4, 13), (23, 13), (23, 1), (13, 2), (1, 2), (1, 100), (2, 200)] ) - assert nx.is_isomorphic(og.inside, expected_graph) - assert og.inputs == [0, 3] - assert og.outputs == [100, 200] + assert nx.is_isomorphic(og.graph, expected_graph) + assert og.input_nodes == [0, 3] + assert og.output_nodes == [100, 200] - outputs_c = {i for i in og.inside.nodes if i not in og.outputs} + outputs_c = {i for i in og.graph.nodes if i not in og.output_nodes} assert og.measurements.keys() == outputs_c assert mapping.keys() <= mapping_complete.keys() assert set(mapping.values()) <= set(mapping_complete.values()) @@ -154,7 +745,7 @@ def test_compose_3() -> None: inputs = [0, 3] outputs = [13, 23] meas = {i: Measurement(0, Plane.XY) for i in g.nodes - set(outputs)} - og_1 = OpenGraph(g, meas, inputs, outputs) + og_1 = OpenGraph(g, inputs, outputs, meas) mapping = {i: i for i in g.nodes} @@ -187,13 +778,13 @@ def test_compose_4() -> None: inputs = [17, 18] outputs = [3, 17] meas = {i: Measurement(0, Plane.XY) for i in g.nodes - set(outputs)} - og_1 = OpenGraph(g, meas, inputs, outputs) + og_1 = OpenGraph(g, inputs, outputs, meas) g = nx.Graph([(1, 2), (2, 3)]) inputs = [1] outputs = [3] meas = {i: Measurement(0, Plane.XY) for i in g.nodes - set(outputs)} - og_2 = OpenGraph(g, meas, inputs, outputs) + og_2 = OpenGraph(g, inputs, outputs, meas) mapping = {1: 17, 3: 300} @@ -201,11 +792,14 @@ def test_compose_4() -> None: expected_graph: nx.Graph[int] expected_graph = nx.Graph([(18, 17), (17, 3), (17, 2), (2, 300)]) - assert nx.is_isomorphic(og.inside, expected_graph) - assert og.inputs == [17, 18] # the input character of node 17 is kept because node 1 (in G2) is an input - assert og.outputs == [3, 300] # the output character of node 17 is lost because node 1 (in G2) is not an output - - outputs_c = {i for i in og.inside.nodes if i not in og.outputs} + assert nx.is_isomorphic(og.graph, expected_graph) + assert og.input_nodes == [17, 18] # the input character of node 17 is kept because node 1 (in G2) is an input + assert og.output_nodes == [ + 3, + 300, + ] # the output character of node 17 is lost because node 1 (in G2) is not an output + + outputs_c = {i for i in og.graph.nodes if i not in og.output_nodes} assert og.measurements.keys() == outputs_c assert mapping.keys() <= mapping_complete.keys() assert set(mapping.values()) <= set(mapping_complete.values()) @@ -233,13 +827,13 @@ def test_compose_5() -> None: inputs = [1, 3] outputs = [2] meas = {i: Measurement(0, Plane.XY) for i in g.nodes - set(outputs)} - og_1 = OpenGraph(g, meas, inputs, outputs) + og_1 = OpenGraph(g, inputs, outputs, meas) g = nx.Graph([(3, 4)]) inputs = [3] outputs = [4] meas = {i: Measurement(0, Plane.XY) for i in g.nodes - set(outputs)} - og_2 = OpenGraph(g, meas, inputs, outputs) + og_2 = OpenGraph(g, inputs, outputs, meas) mapping = {4: 1, 3: 300} @@ -247,19 +841,11 @@ def test_compose_5() -> None: expected_graph: nx.Graph[int] expected_graph = nx.Graph([(1, 2), (1, 3), (1, 300)]) - assert nx.is_isomorphic(og.inside, expected_graph) - assert og.inputs == [3, 300] - assert og.outputs == [2] + assert nx.is_isomorphic(og.graph, expected_graph) + assert og.input_nodes == [3, 300] + assert og.output_nodes == [2] - outputs_c = {i for i in og.inside.nodes if i not in og.outputs} + outputs_c = {i for i in og.graph.nodes if i not in og.output_nodes} assert og.measurements.keys() == outputs_c assert mapping.keys() <= mapping_complete.keys() assert set(mapping.values()) <= set(mapping_complete.values()) - - -def test_double_entanglement() -> None: - pattern = Pattern(input_nodes=[0, 1], cmds=[command.E((0, 1)), command.E((0, 1))]) - pattern2 = OpenGraph.from_pattern(pattern).to_pattern() - state = pattern.simulate_pattern() - state2 = pattern2.simulate_pattern() - assert np.abs(np.dot(state.flatten().conjugate(), state2.flatten())) == pytest.approx(1) diff --git a/tests/test_pyzx.py b/tests/test_pyzx.py index 53e8644c..a468cca0 100644 --- a/tests/test_pyzx.py +++ b/tests/test_pyzx.py @@ -8,7 +8,6 @@ import pytest from numpy.random import PCG64, Generator -from graphix.opengraph import OpenGraph from graphix.random_objects import rand_circuit from graphix.transpiler import Circuit @@ -81,7 +80,7 @@ def test_random_circuit(fx_bg: PCG64, jumps: int) -> None: depth = 5 circuit = rand_circuit(nqubits, depth, rng) pattern = circuit.transpile().pattern - opengraph = OpenGraph.from_pattern(pattern) + opengraph = pattern.extract_opengraph() zx_graph = to_pyzx_graph(opengraph) opengraph2 = from_pyzx_graph(zx_graph) pattern2 = opengraph2.to_pattern() @@ -113,7 +112,7 @@ def test_full_reduce_toffoli() -> None: c = Circuit(3) c.ccx(0, 1, 2) p = c.transpile().pattern - og = OpenGraph.from_pattern(p) + og = p.extract_opengraph() pyg = to_pyzx_graph(og) pyg.normalize() pyg_copy = deepcopy(pyg) diff --git a/tests/test_visualization.py b/tests/test_visualization.py index a5604565..d74ed359 100644 --- a/tests/test_visualization.py +++ b/tests/test_visualization.py @@ -6,11 +6,14 @@ from typing import TYPE_CHECKING import matplotlib.pyplot as plt +import networkx as nx import pytest from graphix import Circuit, Pattern, command, gflow, visualization +from graphix.fundamentals import Plane +from graphix.measurements import Measurement +from graphix.opengraph import OpenGraph, OpenGraphError from graphix.visualization import GraphVisualizer -from tests.test_generator import example_flow, example_gflow, example_pflow if TYPE_CHECKING: from collections.abc import Callable @@ -19,6 +22,75 @@ from numpy.random import Generator +def example_flow(rng: Generator) -> Pattern: + graph: nx.Graph[int] = nx.Graph([(0, 3), (1, 4), (2, 5), (1, 3), (2, 4), (3, 6), (4, 7), (5, 8)]) + inputs = [1, 0, 2] # non-trivial order to check order is conserved. + outputs = [7, 6, 8] + angles = (2 * rng.random(6)).tolist() + measurements = {node: Measurement(angle, Plane.XY) for node, angle in enumerate(angles)} + + pattern = OpenGraph(graph=graph, input_nodes=inputs, output_nodes=outputs, measurements=measurements).to_pattern() + pattern.standardize() + + assert pattern.input_nodes == inputs + assert pattern.output_nodes == outputs + return pattern + + +def example_gflow(rng: Generator) -> Pattern: + graph: nx.Graph[int] = nx.Graph([(1, 2), (2, 3), (3, 4), (4, 5), (5, 6), (3, 6), (1, 6)]) + inputs = [3, 1, 5] + outputs = [4, 2, 6] + angles = dict(zip([1, 3, 5], (2 * rng.random(3)).tolist())) + measurements = {node: Measurement(angle, Plane.XY) for node, angle in angles.items()} + + pattern = OpenGraph(graph=graph, input_nodes=inputs, output_nodes=outputs, measurements=measurements).to_pattern() + pattern.standardize() + + assert pattern.input_nodes == inputs + assert pattern.output_nodes == outputs + return pattern + + +def example_pflow(rng: Generator) -> Pattern: + """Create a graph which has pflow but no gflow. + + Parameters + ---------- + rng : :class:`numpy.random.Generator` + See graphix.tests.conftest.py + + Returns + ------- + Pattern: :class:`graphix.pattern.Pattern` + """ + graph: nx.Graph[int] = nx.Graph( + [(0, 2), (1, 4), (2, 3), (3, 4), (2, 5), (3, 6), (4, 7), (5, 6), (6, 7), (5, 8), (7, 9)] + ) + inputs = [1, 0] + outputs = [9, 8] + + # Heuristic mixture of Pauli and non-Pauli angles ensuring there's no gflow but there's pflow. + meas_angles: dict[int, float] = { + **dict.fromkeys(range(4), 0), + **dict(zip(range(4, 8), (2 * rng.random(4)).tolist())), + } + measurements = {i: Measurement(angle, Plane.XY) for i, angle in meas_angles.items()} + + og = OpenGraph(graph=graph, input_nodes=inputs, output_nodes=outputs, measurements=measurements) + + try: + og.extract_gflow() # example graph doesn't have gflow + except OpenGraphError: + og.extract_pauli_flow() # example graph has Pauli flow + + pattern = og.to_pattern() + pattern.standardize() + assert og.input_nodes == pattern.input_nodes + assert og.output_nodes == pattern.output_nodes + return pattern + + def test_get_pos_from_flow() -> None: circuit = Circuit(1) circuit.h(0)