feat: added expectiation function to compute expectation value of a symbolic hamiltonian using quimb backend. Added function to convers Qibo symbolic hamiltonian in Quimb compatible list of observables. Added a grouping function to combine in a single list local observables acting on the same qubits, to allow Quimb to compute them all in a single contraction (not working for MPS ansatz). Also added argument in setup_backend_specifics to specify an optimizer.

This commit is contained in:
Mattia Robbiano
2025-09-08 18:25:08 +02:00
parent 1fa1730fb3
commit 5e141b200a

View File

@@ -1,12 +1,16 @@
from collections import Counter
import re
from collections import Counter, defaultdict
import numpy as np
import quimb.tensor as qtn
import quimb as qu
from qibo.backends import NumpyBackend
from qibo.config import raise_error
from qibo.result import QuantumState
from qibotn.backends.abstract import QibotnBackend
from qibotn.result import TensorNetworkResult
import warnings
class QuimbBackend(QibotnBackend, NumpyBackend):
@@ -22,7 +26,7 @@ class QuimbBackend(QibotnBackend, NumpyBackend):
def configure_tn_simulation(
self,
ansatz: str = "MPS",
ansatz: str = "any",
max_bond_dimension: int = 10,
n_most_frequent_states: int = 100,
):
@@ -31,7 +35,8 @@ class QuimbBackend(QibotnBackend, NumpyBackend):
Args:
ansatz : str, optional
The tensor network ansatz to use. Currently, only "MPS" is supported. Default is "MPS".
The tensor network ansatz to use. Currently, only "MPS" or "any" is supported. In the second case
the generic Circuit Quimb class is used.
max_bond_dimension : int, optional
The maximum bond dimension for the MPS ansatz. Default is 10.
@@ -43,13 +48,16 @@ class QuimbBackend(QibotnBackend, NumpyBackend):
self.max_bond_dimension = max_bond_dimension
self.n_most_frequent_states = n_most_frequent_states
def setup_backend_specifics(self, qimb_backend="numpy"):
def setup_backend_specifics(self, qimb_backend="numpy", optimizer='auto-hq'):
"""Setup backend specifics.
Args:
qimb_backend: str
The backend to use for the quimb tensor network simulation.
optimizer: str, optional
The optimizer to use for the quimb tensor network simulation.
"""
self.backend = qimb_backend
self.optimizer = optimizer
def execute_circuit(
self,
@@ -57,7 +65,6 @@ class QuimbBackend(QibotnBackend, NumpyBackend):
initial_state=None,
nshots=None,
return_array=False,
**prob_kwargs,
):
"""
Execute a quantum circuit using the specified tensor network ansatz and initial state.
@@ -71,10 +78,6 @@ class QuimbBackend(QibotnBackend, NumpyBackend):
The number of shots for sampling the circuit. If None, no sampling is performed, and the full statevector is used.
return_array : bool, optional
If True, returns the statevector as a dense array. Default is False.
n_most_frequent_states : int, optional
The number of most frequent computational basis states to return. Default is 100.
**prob_kwargs : dict, optional
Additional keyword arguments for probability computation (currently unused).
Returns:
TensorNetworkResult
@@ -131,7 +134,7 @@ class QuimbBackend(QibotnBackend, NumpyBackend):
frequencies = None
measured_probabilities = None
statevector = circ_quimb.to_dense() if return_array else None
statevector = circ_quimb.to_dense(backend=self.backend, optimize=self.optimizer) if return_array else None
return TensorNetworkResult(
nqubits=circuit.nqubits,
backend=self,
@@ -140,3 +143,136 @@ class QuimbBackend(QibotnBackend, NumpyBackend):
prob_type="default",
statevector=statevector,
)
def expectation(self, circuit, observable):
"""Compute the expectation value of a Qibo-friendly ``observable`` on the Tensor Network constructed from a Qibo ``circuit``.
This method takes a Qibo-style symbolic Hamiltonian (e.g., `X(0)*Z(1) + 2.0*Y(2)*Z(0)`)
as the observable, converts it into a Quimb observable and computes its expectation
value using the provided circuit. In case of multiple terms on the same group of qubits, they can be computed in a single contraction.
A grouping procedure is applied to optimize the number of contractions performed.
Args:
circuit: A Qibo quantum circuit object on which the expectation value
is computed.
observable: The observable whose expectation value we want to compute.
This must be provided in the symbolic Hamiltonian form supported by Qibo
(e.g., `X(0)*Y(1)` or `Z(0)*Z(1) + 1.5*Y(2)`).
Returns:
float: The expectation value (real part).
"""
# Map the Qibo observable to Quimb operators and group local operators on the same sites
# for computing them in a single contraction. This does not work with CircuitMPS for some now
# for Quimb 1.11.1
operators_list, sites_list, coeffs_list = self._qiboobs_to_quimbobs(observable)
sites_list_grouped, operators_list_grouped, coeffs_list_grouped = self._group_by_tuples(sites_list, operators_list, coeffs_list)
if self.ansatz == "MPS":
if len(sites_list)-len(sites_list_grouped) > 10:
warnings.warn(
"More than 10 local operators on the same sites are not being grouped as this is not compatible with CircuitMPS. Expected value computation can be more efficient without an MPS ansatz."
)
circ_ansatz = (qtn.circuit.CircuitMPS)
circ = circ_ansatz.from_openqasm2_str(circuit.to_qasm())
expectation_value = 0.0
for ops, sites, coeffs in zip(operators_list, sites_list, coeffs_list):
exp_values = circ.local_expectation(ops, where=sites, backend=self.backend, optimize=self.optimizer)
expectation_value += np.dot(coeffs, exp_values)
return np.real(expectation_value)
else:
circ_ansatz = qtn.circuit.Circuit
circ = circ_ansatz.from_openqasm2_str(circuit.to_qasm())
expectation_value = 0.0
for ops, sites, coeffs in zip(operators_list_grouped, sites_list_grouped, coeffs_list_grouped):
exp_values = circ.local_expectation(ops, where=sites, backend=self.backend, optimize=self.optimizer)
expectation_value += np.dot(coeffs, exp_values)
return np.real(expectation_value)
def _qiboobs_to_quimbobs(
self,
hamiltonian
):
"""
Convert a Qibo SymbolicHamiltonian into a Quimb-compatible decomposition.
Returns three lists:
- operators_list: Quimb operators (tensor products of Pauli matrices).
- sites_list: tuples of qubit indices the operators act on.
- coeffs_list: coefficients for each term.
"""
factor_pattern = re.compile(r"([^\d]+)(\d+)")
operators_list = []
sites_list = []
coeffs_list = []
for term in hamiltonian.terms:
coeff = term.coefficient
term_ops = []
term_sites = []
for factor in term.factors:
match = factor_pattern.match(str(factor))
if not match:
raise ValueError(
f"Factor '{str(factor)}' does not match the expected format."
)
operator_name = match.group(1)
qubit_index = int(match.group(2))
# Build the single-qubit operator
if operator_name not in {"X", "Y", "Z", "I"}:
raise ValueError(f"Unsupported operator {operator_name}")
op = qu.pauli(operator_name)
term_ops.append(op)
term_sites.append(qubit_index)
# Build the tensor product if more than one factor
if term_ops:
full_op = term_ops[0]
for op in term_ops[1:]:
full_op = full_op & op
else:
# Identity term (just coefficient)
full_op = qu.eye(2)
operators_list.append(full_op)
sites_list.append(tuple(term_sites))
coeffs_list.append(coeff)
return operators_list, sites_list, coeffs_list
def _group_by_tuples(self, A, B, C):
"""
Groups the elements of B and C by the unique tuples in A.
Parameters:
A (list of tuples): key tuples (can contain duplicates)
B (list): values aligned with A
C (list): values aligned with A
Returns:
(A_new, B_new, C_new):
A_new: list of unique tuples
B_new: list of lists of grouped values from B
C_new: list of lists of grouped values from C
"""
grouped_B = defaultdict(list)
grouped_C = defaultdict(list)
for a, b, c in zip(A, B, C):
grouped_B[a].append(b)
grouped_C[a].append(c)
A_new = list(grouped_B.keys())
B_new = list(grouped_B.values())
C_new = list(grouped_C.values())
return A_new, B_new, C_new