Fix #45: Parameterized measurement patterns

CHANGELOG.md

@@ -13,6 +13,22 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 - The following changes were made (#155):
   - Added `class Command` and all its child classes that represent all the pattern commands.
   - Added `class Instruction` for the gate network expression in quantum circuit model. Every instruction can be instanciated using this class by passing its name as defined in the Enum `InstructionName`.
+- Parameterized circuits and patterns: angles in instructions and
+  measures can be expressions with parameters created with
+  `parameter.Placeholder` class.  Parameterized circuits can be
+  transpiled and parameterized patterns can be optimized
+  (standardization, minimization, signal shifting and Pauli
+  preprocessing) before being instantiated with the method `subs`.  An
+  additional package,
+  [graphix-symbolic](https://github.com/TeamGraphix/graphix-symbolic),
+  provides parameters that suppor symbolic simulation, and the
+  resulting (symbolic) state vector or density matrix can be
+  instantiated with the method `subs` (probabilities cannot be
+  computed symbolically, so `pr_calc=False` should be passed to
+  simulators for symbolic computation, and an arbitrary path will be
+  computed).
+- Simulator back-ends have an additional optional argument `rng`,
+  to specify the random generator to use during the simulation.
 
 ### Fixed
 

docs/source/data.rst

@@ -70,7 +70,31 @@ This module defines standard data structure for gate seqence (circuit model) use
 
 .. autoclass:: CNOT
 
+:mod:`graphix.parameter` module
++++++++++++++++++++++++++++++++
 
+This module defines parameter objects and parameterized expressions.
+Parameterized expressions can appear in measurement angles in patterns
+and rotation angles in circuits, and they can be substituted with
+actual values.
 
+The module provides generic interfaces for parameters and expressions,
+as well as a simple :class:`Placeholder` class that can be used in
+affine expressions (:class:`AffineExpression`). Affine expressions are
+sufficient for transpiling and pattern optimizations (such as
+standardization, minimization, signal shifting, and Pauli
+preprocessing), but they do not support simulation.
 
+Parameter objects that support symbolic simulation with `sympy` are
+available in a separate package:
+https://github.com/TeamGraphix/graphix-symbolic.
 
+.. currentmodule:: graphix.parameter
+
+.. autoclass:: Expression
+
+.. autoclass:: Parameter
+
+.. autoclass:: AffineExpression
+
+.. autoclass:: Placeholder

examples/MBQCvqe.py

@@ -19,55 +19,77 @@
 the expectation value of the Hamiltonian, effectively finding the ground state energy.
 """
 
+import itertools
+import sys
+from timeit import timeit
+from typing import Iterable
+
 import numpy as np
+import numpy.typing as npt
 from scipy.optimize import minimize
 
 from graphix import Circuit
+from graphix.parameter import Placeholder
+from graphix.pattern import Pattern
 from graphix.simulator import PatternSimulator
+from graphix.transpiler import Angle
 
 
 # %%
 # Define the Hamiltonian for the VQE problem (Example: H = Z0Z1 + X0 + X1)
-def create_hamiltonian():
+def create_hamiltonian() -> npt.NDArray:
     Z = np.array([[1, 0], [0, -1]])
     X = np.array([[0, 1], [1, 0]])
     H = np.kron(Z, Z) + np.kron(X, np.eye(2)) + np.kron(np.eye(2), X)
     return H
 
 
+if sys.version_info >= (3, 12):
+    batched = itertools.batched
+else:
+    # From https://docs.python.org/3/library/itertools.html#itertools.batched
+    def batched(iterable, n):
+        # batched('ABCDEFG', 3) → ABC DEF G
+        if n < 1:
+            raise ValueError("n must be at least one")
+        iterator = iter(iterable)
+        while batch := tuple(itertools.islice(iterator, n)):
+            yield batch
+
+
 # %%
 # Function to build the VQE circuit
-def build_vqe_circuit(n_qubits, params):
+def build_vqe_circuit(n_qubits: int, params: Iterable[Angle]) -> Circuit:
     circuit = Circuit(n_qubits)
-    for i in range(n_qubits):
-        circuit.rx(i, params[i])
-        circuit.ry(i, params[i + n_qubits])
-        circuit.rz(i, params[i + 2 * n_qubits])
+    for i, (x, y, z) in enumerate(batched(params, n=3)):
+        circuit.rx(i, x)
+        circuit.ry(i, y)
+        circuit.rz(i, z)
     for i in range(n_qubits - 1):
         circuit.cnot(i, i + 1)
     return circuit
 
 
 # %%
 class MBQCVQE:
-    def __init__(self, n_qubits, hamiltonian):
+    def __init__(self, n_qubits: int, hamiltonian: npt.NDArray):
         self.n_qubits = n_qubits
         self.hamiltonian = hamiltonian
 
     # %%
     # Function to build the MBQC pattern
-    def build_mbqc_pattern(self, params):
+    def build_mbqc_pattern(self, params: Iterable[Angle]) -> Pattern:
         circuit = build_vqe_circuit(self.n_qubits, params)
         pattern = circuit.transpile().pattern
         pattern.standardize()
         pattern.shift_signals()
+        pattern.perform_pauli_measurements()  # Perform Pauli measurements
         return pattern
 
     # %%
     # Function to simulate the MBQC circuit
-    def simulate_mbqc(self, params, backend="tensornetwork"):
+    def simulate_mbqc(self, params: Iterable[float], backend="tensornetwork"):
         pattern = self.build_mbqc_pattern(params)
-        pattern.perform_pauli_measurements()  # Perform Pauli measurements
         simulator = PatternSimulator(pattern, backend=backend)
         if backend == "tensornetwork":
             tn = simulator.run()  # Simulate the MBQC circuit using tensor network
@@ -81,22 +103,32 @@ def simulate_mbqc(self, params, backend="tensornetwork"):
 
     # %%
     # Function to compute the energy
-    def compute_energy(self, params):
+    def compute_energy(self, params: Iterable[float]):
         # Simulate the MBQC circuit using tensor network backend
         tn = self.simulate_mbqc(params, backend="tensornetwork")
         # Compute the expectation value using MBQCTensornet.expectation_value
         energy = tn.expectation_value(self.hamiltonian, qubit_indices=range(self.n_qubits))
         return energy
 
 
+class MBQCVQEWithPlaceholders(MBQCVQE):
+    def __init__(self, n_qubits: int, hamiltonian) -> None:
+        super().__init__(n_qubits, hamiltonian)
+        self.placeholders = tuple(Placeholder(f"{r}[{q}]") for q in range(n_qubits) for r in ("X", "Y", "Z"))
+        self.pattern = super().build_mbqc_pattern(self.placeholders)
+
+    def build_mbqc_pattern(self, params):
+        return self.pattern.xreplace({placeholder: value for placeholder, value in zip(self.placeholders, params)})
+
+
 # %%
 # Set parameters for VQE
 n_qubits = 2
 hamiltonian = create_hamiltonian()
 
 # %%
 # Instantiate the MBQCVQE class
-mbqc_vqe = MBQCVQE(n_qubits, hamiltonian)
+mbqc_vqe = MBQCVQEWithPlaceholders(n_qubits, hamiltonian)
 
 
 # %%
@@ -109,9 +141,14 @@ def cost_function(params):
 # Random initial parameters
 initial_params = np.random.rand(n_qubits * 3)
 
+
 # %%
 # Perform the optimization using COBYLA
-result = minimize(cost_function, initial_params, method="COBYLA", options={"maxiter": 100})
+def compute():
+    return minimize(cost_function, initial_params, method="COBYLA", options={"maxiter": 100})
+
+
+result = compute()
 
 print(f"Optimized parameters: {result.x}")
 print(f"Optimized energy: {result.fun}")
@@ -120,3 +157,14 @@ def cost_function(params):
 # Compare with the analytical solution
 analytical_solution = -np.sqrt(2) - 1
 print(f"Analytical solution: {analytical_solution}")
+
+# %%
+# Compare performances between using parameterized circuits (with placeholders) or not
+
+mbqc_vqe = MBQCVQEWithPlaceholders(n_qubits, hamiltonian)
+time_with_placeholders = timeit(compute, number=2)
+print(f"Time with placeholders: {time_with_placeholders}")
+
+mbqc_vqe = MBQCVQE(n_qubits, hamiltonian)
+time_without_placeholders = timeit(compute, number=2)
+print(f"Time without placeholders: {time_without_placeholders}")

examples/qnn.py

@@ -22,6 +22,7 @@
 from scipy.optimize import minimize
 from sklearn.datasets import make_circles
 
+from graphix.parameter import Placeholder
 from graphix.transpiler import Circuit
 
 np.random.seed(0)
@@ -164,25 +165,23 @@ def get_expectation_value(self, sv):
         exp_val = np.dot(sv.conj(), exp_val)
         return exp_val.real
 
-    def compute_expectation(self, data_point, params):
+    def compute_expectation(self, pattern, data_point_placeholders, data_point):
         """
         Computes the expectation value of a quantum circuit given a data point and
         parameters.
 
         Args:
+          pattern: Data re-uploading MBQC pattern parameterized by data-point placeholders
+          data_point_placeholders: Data-point placeholders
           data_point: Input to the quantum circuit represented as a 1D numpy array.
-          params: The `params` parameter is a set of parameters that are used to
-        construct a quantum circuit. The specific details of what these parameters
-        represent is described in  `data_reuploading_circuit` method.
 
         Returns:
           the expectation value of a quantum circuit, which is computed using the
         statevector of the output state of the circuit.
         """
-        circuit = self.data_reuploading_circuit(data_point, params)
-        pattern = circuit.transpile().pattern
-        pattern.standardize()
-        pattern.shift_signals()
+        pattern = pattern.xreplace(
+            {placeholder: data for placeholder, data in zip(data_point_placeholders, data_point)}
+        )
         out_state = pattern.simulate_pattern("tensornetwork")
         sv = out_state.to_statevector().flatten()
         return self.get_expectation_value(sv)
@@ -207,7 +206,12 @@ def cost(self, params, x, y):
         Returns:
           the cost value
         """
-        y_pred = [self.compute_expectation(data_point, params) for data_point in x]
+        data_point_placeholders = tuple(Placeholder(f"d[{f}]") for f in range(n_features))
+        circuit = self.data_reuploading_circuit(data_point_placeholders, params)
+        pattern = circuit.transpile().pattern
+        pattern.standardize()
+        pattern.shift_signals()
+        y_pred = [self.compute_expectation(pattern, data_point_placeholders, data_point) for data_point in x]
         cost_val = np.mean(np.abs(y - y_pred))
         self.cost_values.append(cost_val)
         return cost_val

graphix/command.py

@@ -7,9 +7,11 @@
 from typing import TYPE_CHECKING
 
 import numpy as np
-from pydantic import BaseModel
+from pydantic import BaseModel, ConfigDict
 
+# TCH001: ExpressionOrFloat is used in pydantic models
 import graphix.clifford
+from graphix.parameter import ExpressionOrFloat  # noqa: TCH001
 from graphix.pauli import Plane
 
 if TYPE_CHECKING:
@@ -51,10 +53,12 @@ class M(Command):
     Measurement command. By default the plane is set to 'XY', the angle to 0, empty domains and identity vop.
     """
 
+    model_config = ConfigDict(arbitrary_types_allowed=True)
+
     kind: CommandKind = CommandKind.M
     node: Node
     plane: Plane = Plane.XY
-    angle: float = 0.0
+    angle: ExpressionOrFloat = 0.0
     s_domain: set[Node] = set()
     t_domain: set[Node] = set()
 

graphix/instruction.py

@@ -3,9 +3,10 @@
 import abc
 import enum
 
-from pydantic import BaseModel
+from pydantic import BaseModel, ConfigDict
 
-# MEMO: Cannot use TYPE_CHECKING here for pydantic
+# TCH001: ExpressionOrFloat and Plane are used in pydantic models
+from graphix.parameter import ExpressionOrFloat  # noqa: TCH001
 from graphix.pauli import Plane  # noqa: TCH001
 
 
@@ -59,7 +60,9 @@ class RotationInstruction(OneQubitInstruction):
     Rotation instruction base class model.
     """
 
-    angle: float
+    model_config = ConfigDict(arbitrary_types_allowed=True)
+
+    angle: ExpressionOrFloat
 
 
 class OneControlInstruction(OneQubitInstruction):
@@ -180,9 +183,11 @@ class M(OneQubitInstruction):
     M circuit instruction.
     """
 
+    model_config = ConfigDict(arbitrary_types_allowed=True)
+
     kind: InstructionKind = InstructionKind.M
     plane: Plane
-    angle: float
+    angle: ExpressionOrFloat
 
 
 class RX(RotationInstruction):