Variational Quantum Eigensolver (VQE) with Measurement-Based Quantum Computing (MBQC)

In this example, we solve a simple VQE problem using a measurement-based quantum computing (MBQC) approach. The Hamiltonian for the system is given by:

\[H = Z_0 Z_1 + X_0 + X_1\]

where \(Z\) and \(X\) are the Pauli-Z and Pauli-X matrices, respectively.

This Hamiltonian corresponds to a simple model system often used in quantum computing to demonstrate algorithms like VQE. The goal is to find the ground state energy of this Hamiltonian.

We will build a parameterized quantum circuit and optimize its parameters to minimize the expectation value of the Hamiltonian, effectively finding the ground state energy.

from __future__ import annotations

import itertools
import sys
from timeit import timeit
from typing import TYPE_CHECKING

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.simulator import PatternSimulator

if TYPE_CHECKING:
    from collections.abc import Iterable

    from scipy.optimize import OptimizeResult

    from graphix.fundamentals import ParameterizedAngle
    from graphix.pattern import Pattern
    from graphix.sim.tensornet import MBQCTensorNet

Z = np.array([[1, 0], [0, -1]])
X = np.array([[0, 1], [1, 0]])

Define the Hamiltonian for the VQE problem (Example: H = Z0Z1 + X0 + X1)

def create_hamiltonian() -> npt.NDArray[np.float64]:
    return np.kron(Z, Z) + np.kron(X, np.eye(2)) + np.kron(np.eye(2), X)


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: int, params: Iterable[ParameterizedAngle]) -> Circuit:
    circuit = Circuit(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: int, hamiltonian: npt.NDArray[np.float64]):
        self.n_qubits = n_qubits
        self.hamiltonian = hamiltonian

    # %%
    # Function to build the MBQC pattern
    def build_mbqc_pattern(self, params: Iterable[ParameterizedAngle]) -> Pattern:
        circuit = build_vqe_circuit(self.n_qubits, params)
        pattern = circuit.transpile().pattern
        pattern.standardize()
        pattern.shift_signals()
        pattern.remove_input_nodes()
        pattern.perform_pauli_measurements()  # Perform Pauli measurements
        return pattern

    # %%
    # Function to simulate the MBQC circuit
    def simulate_mbqc(
        self,
        params: Iterable[float],
    ) -> MBQCTensorNet:
        pattern = self.build_mbqc_pattern(params)
        simulator = PatternSimulator(pattern, backend="tensornetwork")
        state = simulator.backend.state
        simulator.run()  # Simulate the MBQC circuit using tensor network
        tn = simulator.backend.state
        tn.default_output_nodes = pattern.output_nodes  # Set the default_output_nodes attribute
        if tn.default_output_nodes is None:
            raise ValueError("Output nodes are not set for tensor network simulation.")
        return state

    # %%
    # Function to compute the energy
    def compute_energy(self, params: Iterable[float]) -> float:
        # Simulate the MBQC circuit using tensor network backend
        tn = self.simulate_mbqc(params)
        # Compute the expectation value using MBQCTensornet.expectation_value
        return tn.expectation_value(self.hamiltonian.astype(np.complex128), qubit_indices=range(self.n_qubits))


class MBQCVQEWithPlaceholders(MBQCVQE):
    def __init__(self, n_qubits: int, hamiltonian: npt.NDArray[np.float64]) -> 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: Iterable[ParameterizedAngle]) -> Pattern:
        return self.pattern.xreplace(dict(zip(self.placeholders, params, strict=True)))

Set parameters for VQE

n_qubits = 2
hamiltonian = create_hamiltonian()

Instantiate the MBQCVQE class

mbqc_vqe: MBQCVQE = MBQCVQEWithPlaceholders(n_qubits, hamiltonian)

Define the cost function

def cost_function(params: Iterable[float]) -> float:
    return mbqc_vqe.compute_energy(params)

Random initial parameters

rng = np.random.default_rng()
initial_params = rng.random(n_qubits * 3)

Perform the optimization using COBYLA

def compute() -> OptimizeResult:
    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}")

/home/docs/checkouts/readthedocs.org/user_builds/graphix/checkouts/stable/examples/mbqc_vqe.py:108: UserWarning: Default random-number generator is used. Results may not be reproducible.
  simulator.run()  # Simulate the MBQC circuit using tensor network
/home/docs/checkouts/readthedocs.org/user_builds/graphix/checkouts/stable/examples/mbqc_vqe.py:108: UserWarning: Default random-number generator is used. Results may not be reproducible.
  simulator.run()  # Simulate the MBQC circuit using tensor network
/home/docs/checkouts/readthedocs.org/user_builds/graphix/checkouts/stable/examples/mbqc_vqe.py:108: UserWarning: Default random-number generator is used. Results may not be reproducible.
  simulator.run()  # Simulate the MBQC circuit using tensor network
Optimized parameters: [0.28788715 1.00002127 1.00001902 1.74236528 0.14754458 1.00002686]
Optimized energy: -2.2360679507249825

Compare with the analytical solution

analytical_solution = -np.sqrt(2) - 1
print(f"Analytical solution: {analytical_solution}")

Analytical solution: -2.414213562373095

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}")

Time with placeholders: 1.282748126000115
Time without placeholders: 1.6936494389997279