From e756691c6b1e9d749b94c293ec69ec4d13d22e07 Mon Sep 17 00:00:00 2001 From: BrunoLiegiBastonLiegi Date: Mon, 14 Oct 2024 17:48:50 +0200 Subject: [PATCH 01/20] feat: implemented expectation for non diagonal observbales --- src/qibo/hamiltonians/expectation.py | 26 ++++++++++++++++++++++++++ 1 file changed, 26 insertions(+) create mode 100644 src/qibo/hamiltonians/expectation.py diff --git a/src/qibo/hamiltonians/expectation.py b/src/qibo/hamiltonians/expectation.py new file mode 100644 index 0000000000..333fdce40a --- /dev/null +++ b/src/qibo/hamiltonians/expectation.py @@ -0,0 +1,26 @@ +from math import prod + +from qibo import Circuit, gates, symbols +from qibo.config import raise_error +from qibo.hamiltonians import SymbolicHamiltonian + + +def Expectation(circuit: Circuit, observable: SymbolicHamiltonian) -> float: + if len(circuit.measurements) > 0: + raise_error(RuntimeError) + exp_val = 0.0 + for term in observable.terms: + Z_observable = SymbolicHamiltonian( + prod([symbols.Z(q) for q in term.target_qubits]), + nqubits=circuit.nqubits, + backend=observable.backend, + ) + measurements = [ + gates.M(factor.target_qubit, basis=factor.gate.__class__) + for factor in term.factors + ] + circ_copy = circuit.copy(True) + [circ_copy.add(m) for m in measurements] + freq = observable.backend.execute_circuit(circ_copy).frequencies() + exp_val += Z_observable.expectation_from_samples(freq) + return exp_val From 2ced67252723ce68ddc80bcd482fe7298af86899 Mon Sep 17 00:00:00 2001 From: BrunoLiegiBastonLiegi Date: Tue, 15 Oct 2024 08:58:09 +0200 Subject: [PATCH 02/20] feat: making expectation_from_samples take input circuits --- src/qibo/hamiltonians/hamiltonians.py | 69 +++++++++++++++++++++++++-- 1 file changed, 64 insertions(+), 5 deletions(-) diff --git a/src/qibo/hamiltonians/hamiltonians.py b/src/qibo/hamiltonians/hamiltonians.py index 75fd10c7ba..21627eb0ee 100644 --- a/src/qibo/hamiltonians/hamiltonians.py +++ b/src/qibo/hamiltonians/hamiltonians.py @@ -1,7 +1,8 @@ """Module defining Hamiltonian classes.""" from itertools import chain -from typing import Optional +from math import prod +from typing import Optional, Union import numpy as np import sympy @@ -552,17 +553,75 @@ def calculate_dense(self): def expectation(self, state, normalize=False): return Hamiltonian.expectation(self, state, normalize) - def expectation_from_samples(self, freq, qubit_map=None): + def expectation_from_samples( + self, freq: Union[dict, "Circuit"], qubit_map: dict = None, nshots: int = None + ): + """ + Calculate the expectation value starting from the samples. This even works + for observables not completely diagonal in the computational basis, but only + diagonal at a term level in a defined basis. Namely, for an observable of the + form :math:``H = \\sum_i H_i``, where each :math:``H_i`` consists in a `n`-qubits + pauli operator :math:`P_0 \\otimes P_1 \\otimes \\cdots \\otimes P_n`, the + expectation value is computed by rotating the input circuit in the suitable + basis for each term :math:``H_i`` thus extracting the `term-wise` expectations + that are then summed to build the global expectation value. + + Args: + freq (dict | :class:`qibo.models.Circuit`): either the ``dict`` of the + frequencies of the samples or a :class:`qibo.models.Circuit` where to + extract the samples from. A :class:`qibo.models.Circuit` is needed in + case the observable is not diagonal in the computational basis. + qubit_map (dict): Qubit map. + nshots (int, optional): Number of shots for the expecation value + calculation. This is used only when a :class:`qibo.models.Circuit` + is passed as input. + + Returns: + (float) the computed expectation value. + """ + from qibo import Circuit, gates + terms = self.terms for term in terms: + if len(term.factors) != len(set(term.factors)): + raise_error(NotImplementedError, "Z^k is not implemented since Z^2=I.") # pylint: disable=E1101 for factor in term.factors: - if not isinstance(factor, Z): + if not isinstance(factor, Z) and not isinstance(freq, Circuit): raise_error( NotImplementedError, "Observable is not a Z Pauli string." ) - if len(term.factors) != len(set(term.factors)): - raise_error(NotImplementedError, "Z^k is not implemented since Z^2=I.") + + if isinstance(freq, Circuit): + if nshots is None: + nshots = 1000 + rotated_circuits = [] + for term in terms: + Z_observable = SymbolicHamiltonian( + prod([Z(q) for q in term.target_qubits]), + nqubits=freq.nqubits, + backend=self.backend, + ) + measurements = [ + gates.M(factor.target_qubit, basis=factor.gate.__class__) + for factor in term.factors + ] + circ_copy = freq.copy(True) + [circ_copy.add(m) for m in measurements] + rotated_circuits.append(circ_copy) + frequencies = [ + result.frequencies() + for result in self.backend.execute_circuits( + rotated_circuits, nshots=nshots + ) + ] + return sum( + [ + Z_observable.expectation_from_samples(f, qubit_map) + for f in frequencies + ] + ) + keys = list(freq.keys()) counts = self.backend.cast(list(freq.values()), self.backend.precision) / sum( freq.values() From ed8b641aa2859534def1f16c3264a3b9b02cb4e7 Mon Sep 17 00:00:00 2001 From: BrunoLiegiBastonLiegi Date: Tue, 15 Oct 2024 09:21:13 +0200 Subject: [PATCH 03/20] fix: added coefficients --- src/qibo/hamiltonians/hamiltonians.py | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/src/qibo/hamiltonians/hamiltonians.py b/src/qibo/hamiltonians/hamiltonians.py index 21627eb0ee..d2c72a4212 100644 --- a/src/qibo/hamiltonians/hamiltonians.py +++ b/src/qibo/hamiltonians/hamiltonians.py @@ -596,7 +596,9 @@ def expectation_from_samples( if nshots is None: nshots = 1000 rotated_circuits = [] + coefficients = [] for term in terms: + coefficients.append(term.coefficient) Z_observable = SymbolicHamiltonian( prod([Z(q) for q in term.target_qubits]), nqubits=freq.nqubits, @@ -617,8 +619,8 @@ def expectation_from_samples( ] return sum( [ - Z_observable.expectation_from_samples(f, qubit_map) - for f in frequencies + c * Z_observable.expectation_from_samples(f, qubit_map) + for c, f in zip(coefficients, frequencies) ] ) From db41268913df509c22614a00765c4aa4c38aa5ad Mon Sep 17 00:00:00 2001 From: BrunoLiegiBastonLiegi Date: Tue, 15 Oct 2024 09:25:43 +0200 Subject: [PATCH 04/20] feat: removed expecation.py --- src/qibo/hamiltonians/expectation.py | 26 -------------------------- 1 file changed, 26 deletions(-) delete mode 100644 src/qibo/hamiltonians/expectation.py diff --git a/src/qibo/hamiltonians/expectation.py b/src/qibo/hamiltonians/expectation.py deleted file mode 100644 index 333fdce40a..0000000000 --- a/src/qibo/hamiltonians/expectation.py +++ /dev/null @@ -1,26 +0,0 @@ -from math import prod - -from qibo import Circuit, gates, symbols -from qibo.config import raise_error -from qibo.hamiltonians import SymbolicHamiltonian - - -def Expectation(circuit: Circuit, observable: SymbolicHamiltonian) -> float: - if len(circuit.measurements) > 0: - raise_error(RuntimeError) - exp_val = 0.0 - for term in observable.terms: - Z_observable = SymbolicHamiltonian( - prod([symbols.Z(q) for q in term.target_qubits]), - nqubits=circuit.nqubits, - backend=observable.backend, - ) - measurements = [ - gates.M(factor.target_qubit, basis=factor.gate.__class__) - for factor in term.factors - ] - circ_copy = circuit.copy(True) - [circ_copy.add(m) for m in measurements] - freq = observable.backend.execute_circuit(circ_copy).frequencies() - exp_val += Z_observable.expectation_from_samples(freq) - return exp_val From 6fd2bf143a879c1f0e15bb94e3e8f581e1787540 Mon Sep 17 00:00:00 2001 From: BrunoLiegiBastonLiegi Date: Wed, 16 Oct 2024 12:11:16 +0200 Subject: [PATCH 05/20] test: working on a test --- tests/test_hamiltonians.py | 59 ++++++++++++++++++++++++++++---------- 1 file changed, 44 insertions(+), 15 deletions(-) diff --git a/tests/test_hamiltonians.py b/tests/test_hamiltonians.py index 5dce386185..cabaa3de79 100644 --- a/tests/test_hamiltonians.py +++ b/tests/test_hamiltonians.py @@ -1,11 +1,14 @@ """Test methods in `qibo/core/hamiltonians.py`.""" +from math import prod + import numpy as np import pytest from qibo import Circuit, gates, hamiltonians +from qibo.hamiltonians.hamiltonians import SymbolicHamiltonian from qibo.quantum_info.random_ensembles import random_density_matrix, random_statevector -from qibo.symbols import I, Z +from qibo.symbols import I, X, Y, Z from .utils import random_sparse_matrix @@ -261,23 +264,49 @@ def test_hamiltonian_expectation_errors(backend): h.expectation("test") -def test_hamiltonian_expectation_from_samples(backend): +@pytest.mark.parametrize( + "observable", + [2 * Z(0) * Z(1) + Z(0) * Z(2), X(0) * Z(1) + Y(0) * X(2) / 2 - Z(0) * Y(1)], +) +def test_hamiltonian_expectation_from_samples(backend, observable): """Test Hamiltonian expectation value calculation.""" backend.set_seed(12) - obs0 = 2 * Z(0) * Z(1) + Z(0) * Z(2) - obs1 = 2 * Z(0) * Z(1) + Z(0) * Z(2) * I(3) - h0 = hamiltonians.SymbolicHamiltonian(obs0, backend=backend) - h1 = hamiltonians.SymbolicHamiltonian(obs1, backend=backend) - matrix = backend.to_numpy(h0.matrix) - c = Circuit(4) - c.add(gates.RX(0, np.random.rand())) - c.add(gates.RX(1, np.random.rand())) - c.add(gates.RX(2, np.random.rand())) - c.add(gates.RX(3, np.random.rand())) - c.add(gates.M(0, 1, 2, 3)) + + nqubits = 3 nshots = 10**5 - # result = c(nshots=nshots) - result = backend.execute_circuit(c, nshots=nshots) + c = Circuit(nqubits) + for q in range(nqubits): + c.add(gates.RX(0, np.random.rand())) + + H = hamiltonians.SymbolicHamiltonian(observable, nqubits=nqubits, backend=backend) + matrix = backend.to_numpy(H.matrix) + diagonal = ( + backend.np.count_nonzero(matrix - backend.np.diag(backend.np.diagonal(matrix))) + == 0 + ) + if not diagonal: + matrices = [] + for term in H.terms: + non_Z = [ + factor + for factor in term.factors + if not isinstance(factor.gate.__class__, gates.Z) + ] + diagonalizator = SymbolicHamiltonian( + prod(non_Z), nqubits=nqubits, backend=backend + ).matrix + matrices.append( + backend.np.transpose(backend.np.conj(diagonalizator)) + @ matrix + @ diagonalizator + ) + exp_from_samples = H.expectation_from_samples(c) + state = backend.execute_circuit(c).state() + exp = H.expectation(state) + else: + matrices = [matrix] + c.add(gates.M(*range(nqubits))) + freq = result.frequencies(binary=True) Obs0 = hamiltonians.Hamiltonian( From e4b35918861a247c8909c6556bf4ac2f46e7cfe5 Mon Sep 17 00:00:00 2001 From: BrunoLiegiBastonLiegi Date: Wed, 16 Oct 2024 14:16:41 +0200 Subject: [PATCH 06/20] test: finished to modify the test --- tests/test_hamiltonians.py | 28 +++++++++++----------------- 1 file changed, 11 insertions(+), 17 deletions(-) diff --git a/tests/test_hamiltonians.py b/tests/test_hamiltonians.py index cabaa3de79..f57e7e8aea 100644 --- a/tests/test_hamiltonians.py +++ b/tests/test_hamiltonians.py @@ -6,7 +6,7 @@ import pytest from qibo import Circuit, gates, hamiltonians -from qibo.hamiltonians.hamiltonians import SymbolicHamiltonian +from qibo.hamiltonians.hamiltonians import Hamiltonian, SymbolicHamiltonian from qibo.quantum_info.random_ensembles import random_density_matrix, random_statevector from qibo.symbols import I, X, Y, Z @@ -284,8 +284,9 @@ def test_hamiltonian_expectation_from_samples(backend, observable): backend.np.count_nonzero(matrix - backend.np.diag(backend.np.diagonal(matrix))) == 0 ) + final_state = backend.execute_circuit(c.copy(True)).state() if not diagonal: - matrices = [] + exp = 0.0 for term in H.terms: non_Z = [ factor @@ -295,29 +296,22 @@ def test_hamiltonian_expectation_from_samples(backend, observable): diagonalizator = SymbolicHamiltonian( prod(non_Z), nqubits=nqubits, backend=backend ).matrix - matrices.append( + diag_matrix = ( backend.np.transpose(backend.np.conj(diagonalizator)) @ matrix @ diagonalizator ) + exp += Hamiltonian(nqubits, diag_matrix, backend=backend).expectation( + final_state + ) exp_from_samples = H.expectation_from_samples(c) - state = backend.execute_circuit(c).state() - exp = H.expectation(state) else: - matrices = [matrix] c.add(gates.M(*range(nqubits))) + freq = backend.execute_circuit(c).frequencies() + exp_from_samples = H.expectation_from_samples(freq) + exp = H.expectation(final_state) - freq = result.frequencies(binary=True) - - Obs0 = hamiltonians.Hamiltonian( - 3, matrix, backend=backend - ).expectation_from_samples(freq, qubit_map=None) - Obs0 = backend.cast(Obs0, dtype=Obs0.dtype) - - Obs1 = h1.expectation(result.state()) - Obs1 = backend.cast(Obs1, dtype=Obs1.dtype) - - backend.assert_allclose(Obs0, Obs1, atol=10 / np.sqrt(nshots)) + backend.assert_allclose(exp, exp_from_samples, atol=10 / np.sqrt(nshots)) def test_hamiltonian_expectation_from_samples_errors(backend): From 761eaa8341525b111caada599983db590eb3266b Mon Sep 17 00:00:00 2001 From: BrunoLiegiBastonLiegi Date: Wed, 16 Oct 2024 15:41:57 +0200 Subject: [PATCH 07/20] fix: fixing the test --- tests/test_hamiltonians.py | 32 ++++++-------------------------- 1 file changed, 6 insertions(+), 26 deletions(-) diff --git a/tests/test_hamiltonians.py b/tests/test_hamiltonians.py index f57e7e8aea..1f02f642e7 100644 --- a/tests/test_hamiltonians.py +++ b/tests/test_hamiltonians.py @@ -1,7 +1,5 @@ """Test methods in `qibo/core/hamiltonians.py`.""" -from math import prod - import numpy as np import pytest @@ -273,45 +271,27 @@ def test_hamiltonian_expectation_from_samples(backend, observable): backend.set_seed(12) nqubits = 3 - nshots = 10**5 + nshots = 4 * 10**6 c = Circuit(nqubits) for q in range(nqubits): c.add(gates.RX(0, np.random.rand())) H = hamiltonians.SymbolicHamiltonian(observable, nqubits=nqubits, backend=backend) - matrix = backend.to_numpy(H.matrix) + matrix = H.matrix diagonal = ( backend.np.count_nonzero(matrix - backend.np.diag(backend.np.diagonal(matrix))) == 0 ) final_state = backend.execute_circuit(c.copy(True)).state() + exp = H.expectation(final_state) if not diagonal: - exp = 0.0 - for term in H.terms: - non_Z = [ - factor - for factor in term.factors - if not isinstance(factor.gate.__class__, gates.Z) - ] - diagonalizator = SymbolicHamiltonian( - prod(non_Z), nqubits=nqubits, backend=backend - ).matrix - diag_matrix = ( - backend.np.transpose(backend.np.conj(diagonalizator)) - @ matrix - @ diagonalizator - ) - exp += Hamiltonian(nqubits, diag_matrix, backend=backend).expectation( - final_state - ) - exp_from_samples = H.expectation_from_samples(c) + exp_from_samples = H.expectation_from_samples(c, nshots=nshots) else: c.add(gates.M(*range(nqubits))) - freq = backend.execute_circuit(c).frequencies() + freq = backend.execute_circuit(c, nshots=nshots).frequencies() exp_from_samples = H.expectation_from_samples(freq) - exp = H.expectation(final_state) - backend.assert_allclose(exp, exp_from_samples, atol=10 / np.sqrt(nshots)) + backend.assert_allclose(exp, exp_from_samples, atol=1e-2) def test_hamiltonian_expectation_from_samples_errors(backend): From 4cdfe6ba503892b28a0050d2753b97c97e4b9b7b Mon Sep 17 00:00:00 2001 From: BrunoLiegiBastonLiegi Date: Wed, 16 Oct 2024 16:49:14 +0200 Subject: [PATCH 08/20] feat: splitted the circuit and frequency cases --- src/qibo/hamiltonians/hamiltonians.py | 144 +++++++++++++++----------- 1 file changed, 86 insertions(+), 58 deletions(-) diff --git a/src/qibo/hamiltonians/hamiltonians.py b/src/qibo/hamiltonians/hamiltonians.py index d2c72a4212..4a4a15783e 100644 --- a/src/qibo/hamiltonians/hamiltonians.py +++ b/src/qibo/hamiltonians/hamiltonians.py @@ -553,83 +553,77 @@ def calculate_dense(self): def expectation(self, state, normalize=False): return Hamiltonian.expectation(self, state, normalize) - def expectation_from_samples( - self, freq: Union[dict, "Circuit"], qubit_map: dict = None, nshots: int = None - ): + def _exp_from_circuit(self, circuit: dict, qubit_map: dict, nshots: int = None): """ - Calculate the expectation value starting from the samples. This even works - for observables not completely diagonal in the computational basis, but only - diagonal at a term level in a defined basis. Namely, for an observable of the - form :math:``H = \\sum_i H_i``, where each :math:``H_i`` consists in a `n`-qubits - pauli operator :math:`P_0 \\otimes P_1 \\otimes \\cdots \\otimes P_n`, the - expectation value is computed by rotating the input circuit in the suitable - basis for each term :math:``H_i`` thus extracting the `term-wise` expectations - that are then summed to build the global expectation value. + Calculate the expecation value from a circuit. + This allows for non diagonal observables. Each term of the observable is + treated separately, by measuring in the correct basis and re-executing the + circuit. Args: - freq (dict | :class:`qibo.models.Circuit`): either the ``dict`` of the - frequencies of the samples or a :class:`qibo.models.Circuit` where to - extract the samples from. A :class:`qibo.models.Circuit` is needed in - case the observable is not diagonal in the computational basis. - qubit_map (dict): Qubit map. - nshots (int, optional): Number of shots for the expecation value - calculation. This is used only when a :class:`qibo.models.Circuit` - is passed as input. + freq (dict): input frequencies. + qubit_map (dict): qubit map. Returns: - (float) the computed expectation value. + (float): the calculated expectation value """ from qibo import Circuit, gates - terms = self.terms - for term in terms: - if len(term.factors) != len(set(term.factors)): - raise_error(NotImplementedError, "Z^k is not implemented since Z^2=I.") + if nshots is None: + nshots = 1000 + rotated_circuits = [] + coefficients = [] + for term in self.terms: + coefficients.append(term.coefficient) + Z_observable = SymbolicHamiltonian( + prod([Z(q) for q in term.target_qubits]), + nqubits=circuit.nqubits, + backend=self.backend, + ) + measurements = [ + gates.M(factor.target_qubit, basis=factor.gate.__class__) + for factor in term.factors + ] + circ_copy = circuit.copy(True) + [circ_copy.add(m) for m in measurements] + rotated_circuits.append(circ_copy) + frequencies = [ + result.frequencies() + for result in self.backend.execute_circuits(rotated_circuits, nshots=nshots) + ] + return sum( + [ + c * Z_observable._exp_from_freq(f, qubit_map) + for c, f in zip(coefficients, frequencies) + ] + ) + + def _exp_from_freq(self, freq: dict, qubit_map: dict): + """ + Calculate the expecation value from some frequencies. + The observable has to be diagonal in the computational basis. + + Args: + freq (dict): input frequencies. + qubit_map (dict): qubit map. + + Returns: + (float): the calculated expectation value + """ + for term in self.terms: # pylint: disable=E1101 for factor in term.factors: - if not isinstance(factor, Z) and not isinstance(freq, Circuit): + if not isinstance(factor, Z): raise_error( NotImplementedError, "Observable is not a Z Pauli string." ) - if isinstance(freq, Circuit): - if nshots is None: - nshots = 1000 - rotated_circuits = [] - coefficients = [] - for term in terms: - coefficients.append(term.coefficient) - Z_observable = SymbolicHamiltonian( - prod([Z(q) for q in term.target_qubits]), - nqubits=freq.nqubits, - backend=self.backend, - ) - measurements = [ - gates.M(factor.target_qubit, basis=factor.gate.__class__) - for factor in term.factors - ] - circ_copy = freq.copy(True) - [circ_copy.add(m) for m in measurements] - rotated_circuits.append(circ_copy) - frequencies = [ - result.frequencies() - for result in self.backend.execute_circuits( - rotated_circuits, nshots=nshots - ) - ] - return sum( - [ - c * Z_observable.expectation_from_samples(f, qubit_map) - for c, f in zip(coefficients, frequencies) - ] - ) - keys = list(freq.keys()) counts = self.backend.cast(list(freq.values()), self.backend.precision) / sum( freq.values() ) qubits = [] - for term in terms: + for term in self.terms: qubits_term = [] for k in term.target_qubits: qubits_term.append(k) @@ -649,6 +643,40 @@ def expectation_from_samples( expval += expval_q * self.terms[j].coefficient.real return expval + self.constant.real + def expectation_from_samples( + self, data: Union[dict, "Circuit"], qubit_map: dict = None, nshots: int = None + ): + """ + Calculate the expectation value starting from the samples. This even works + for observables not completely diagonal in the computational basis, but only + diagonal at a term level in a defined basis. Namely, for an observable of the + form :math:``H = \\sum_i H_i``, where each :math:``H_i`` consists in a `n`-qubits + pauli operator :math:`P_0 \\otimes P_1 \\otimes \\cdots \\otimes P_n`, the + expectation value is computed by rotating the input circuit in the suitable + basis for each term :math:``H_i`` thus extracting the `term-wise` expectations + that are then summed to build the global expectation value. + + Args: + freq (dict | :class:`qibo.models.Circuit`): either the ``dict`` of the + frequencies of the samples or a :class:`qibo.models.Circuit` where to + extract the samples from. A :class:`qibo.models.Circuit` is needed in + case the observable is not diagonal in the computational basis. + qubit_map (dict): qubit map. + nshots (int, optional): number of shots for the expecation value + calculation. This is used only when a :class:`qibo.models.Circuit` + is passed as input. + + Returns: + (float) the computed expectation value. + """ + for term in self.terms: + if len(term.factors) != len(set(term.factors)): + raise_error(NotImplementedError, "Z^k is not implemented since Z^2=I.") + + if isinstance(data, dict): + return self._exp_from_freq(data, qubit_map) + return self._exp_from_circuit(data, qubit_map, nshots) + def __add__(self, o): if isinstance(o, self.__class__): if self.nqubits != o.nqubits: From 63b1ce702e5249ea5476b2a256158d4e08a3a1b8 Mon Sep 17 00:00:00 2001 From: BrunoLiegiBastonLiegi Date: Wed, 16 Oct 2024 16:55:08 +0200 Subject: [PATCH 09/20] doc: improved Hamiltonian.exp_from_samples error message --- src/qibo/hamiltonians/hamiltonians.py | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/src/qibo/hamiltonians/hamiltonians.py b/src/qibo/hamiltonians/hamiltonians.py index 4a4a15783e..a7328041b2 100644 --- a/src/qibo/hamiltonians/hamiltonians.py +++ b/src/qibo/hamiltonians/hamiltonians.py @@ -153,7 +153,10 @@ def expectation_from_samples(self, freq, qubit_map=None): ) != 0 ): - raise_error(NotImplementedError, "Observable is not diagonal.") + raise_error( + NotImplementedError, + "Observable is not diagonal. Expectation of non diagonal observables starting from samples is currently supported for `qibo.hamiltonians.hamiltonians.SymbolicHamiltonian` only.", + ) keys = list(freq.keys()) if qubit_map is None: qubit_map = list(range(int(np.log2(len(obs))))) @@ -565,7 +568,7 @@ def _exp_from_circuit(self, circuit: dict, qubit_map: dict, nshots: int = None): qubit_map (dict): qubit map. Returns: - (float): the calculated expectation value + (float): the calculated expectation value. """ from qibo import Circuit, gates @@ -608,7 +611,7 @@ def _exp_from_freq(self, freq: dict, qubit_map: dict): qubit_map (dict): qubit map. Returns: - (float): the calculated expectation value + (float): the calculated expectation value. """ for term in self.terms: # pylint: disable=E1101 From 2b34b5cd0ec2bbdda3e3840a79def3a3eddec822 Mon Sep 17 00:00:00 2001 From: BrunoLiegiBastonLiegi Date: Wed, 16 Oct 2024 17:32:30 +0200 Subject: [PATCH 10/20] feat: added support for power elevation + removed outdated test --- src/qibo/hamiltonians/hamiltonians.py | 5 +---- tests/test_hamiltonians.py | 5 ++++- tests/test_hamiltonians_symbolic.py | 10 ---------- 3 files changed, 5 insertions(+), 15 deletions(-) diff --git a/src/qibo/hamiltonians/hamiltonians.py b/src/qibo/hamiltonians/hamiltonians.py index a7328041b2..df1a85b4d6 100644 --- a/src/qibo/hamiltonians/hamiltonians.py +++ b/src/qibo/hamiltonians/hamiltonians.py @@ -585,7 +585,7 @@ def _exp_from_circuit(self, circuit: dict, qubit_map: dict, nshots: int = None): ) measurements = [ gates.M(factor.target_qubit, basis=factor.gate.__class__) - for factor in term.factors + for factor in set(term.factors) ] circ_copy = circuit.copy(True) [circ_copy.add(m) for m in measurements] @@ -672,9 +672,6 @@ def expectation_from_samples( Returns: (float) the computed expectation value. """ - for term in self.terms: - if len(term.factors) != len(set(term.factors)): - raise_error(NotImplementedError, "Z^k is not implemented since Z^2=I.") if isinstance(data, dict): return self._exp_from_freq(data, qubit_map) diff --git a/tests/test_hamiltonians.py b/tests/test_hamiltonians.py index 1f02f642e7..4bdf6cbd03 100644 --- a/tests/test_hamiltonians.py +++ b/tests/test_hamiltonians.py @@ -264,7 +264,10 @@ def test_hamiltonian_expectation_errors(backend): @pytest.mark.parametrize( "observable", - [2 * Z(0) * Z(1) + Z(0) * Z(2), X(0) * Z(1) + Y(0) * X(2) / 2 - Z(0) * Y(1)], + [ + 2 * Z(0) * Z(1) ** 2 + Z(0) * Z(2), + X(0) * Z(1) + Y(0) * X(2) / 2 - Z(0) * Y(1) ** 3, + ], ) def test_hamiltonian_expectation_from_samples(backend, observable): """Test Hamiltonian expectation value calculation.""" diff --git a/tests/test_hamiltonians_symbolic.py b/tests/test_hamiltonians_symbolic.py index 54d8dac2e3..d1d8485703 100644 --- a/tests/test_hamiltonians_symbolic.py +++ b/tests/test_hamiltonians_symbolic.py @@ -330,16 +330,6 @@ def test_hamiltonian_expectation_from_samples(backend): backend.assert_allclose(ev0, ev1, atol=20 / np.sqrt(nshots)) -def test_hamiltonian_expectation_from_samples_errors(backend): - obs = [Z(0) * Y(1), Z(0) * Z(1) ** 3] - h1 = hamiltonians.SymbolicHamiltonian(obs[0], backend=backend) - h2 = hamiltonians.SymbolicHamiltonian(obs[1], backend=backend) - with pytest.raises(NotImplementedError): - h1.expectation_from_samples(None, qubit_map=None) - with pytest.raises(NotImplementedError): - h2.expectation_from_samples(None, qubit_map=None) - - @pytest.mark.parametrize("density_matrix", [False, True]) @pytest.mark.parametrize("calcterms", [False, True]) def test_symbolic_hamiltonian_abstract_symbol_ev(backend, density_matrix, calcterms): From 68d4e55d904adcba65eca13a8677b545b4a0e6e6 Mon Sep 17 00:00:00 2001 From: BrunoLiegiBastonLiegi Date: Wed, 16 Oct 2024 18:07:43 +0200 Subject: [PATCH 11/20] fix: typos in the docstrings --- src/qibo/hamiltonians/hamiltonians.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/qibo/hamiltonians/hamiltonians.py b/src/qibo/hamiltonians/hamiltonians.py index df1a85b4d6..ee43a884d8 100644 --- a/src/qibo/hamiltonians/hamiltonians.py +++ b/src/qibo/hamiltonians/hamiltonians.py @@ -558,13 +558,13 @@ def expectation(self, state, normalize=False): def _exp_from_circuit(self, circuit: dict, qubit_map: dict, nshots: int = None): """ - Calculate the expecation value from a circuit. + Calculate the expectation value from a circuit. This allows for non diagonal observables. Each term of the observable is treated separately, by measuring in the correct basis and re-executing the circuit. Args: - freq (dict): input frequencies. + circuit (dict): input frequencies. qubit_map (dict): qubit map. Returns: @@ -603,7 +603,7 @@ def _exp_from_circuit(self, circuit: dict, qubit_map: dict, nshots: int = None): def _exp_from_freq(self, freq: dict, qubit_map: dict): """ - Calculate the expecation value from some frequencies. + Calculate the expectation value from some frequencies. The observable has to be diagonal in the computational basis. Args: @@ -660,7 +660,7 @@ def expectation_from_samples( that are then summed to build the global expectation value. Args: - freq (dict | :class:`qibo.models.Circuit`): either the ``dict`` of the + data (dict | :class:`qibo.models.Circuit`): either the ``dict`` of the frequencies of the samples or a :class:`qibo.models.Circuit` where to extract the samples from. A :class:`qibo.models.Circuit` is needed in case the observable is not diagonal in the computational basis. From 63fba048bcd6a19ba9159300090b6c45cea12bc7 Mon Sep 17 00:00:00 2001 From: BrunoLiegiBastonLiegi Date: Fri, 18 Oct 2024 15:16:14 +0200 Subject: [PATCH 12/20] fix: minor fixes to docstring and coverage --- src/qibo/hamiltonians/hamiltonians.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/qibo/hamiltonians/hamiltonians.py b/src/qibo/hamiltonians/hamiltonians.py index ee43a884d8..85ec563e68 100644 --- a/src/qibo/hamiltonians/hamiltonians.py +++ b/src/qibo/hamiltonians/hamiltonians.py @@ -556,7 +556,7 @@ def calculate_dense(self): def expectation(self, state, normalize=False): return Hamiltonian.expectation(self, state, normalize) - def _exp_from_circuit(self, circuit: dict, qubit_map: dict, nshots: int = None): + def _exp_from_circuit(self, circuit: Circuit, qubit_map: dict, nshots: int = None): """ Calculate the expectation value from a circuit. This allows for non diagonal observables. Each term of the observable is @@ -564,7 +564,7 @@ def _exp_from_circuit(self, circuit: dict, qubit_map: dict, nshots: int = None): circuit. Args: - circuit (dict): input frequencies. + circuit (Circuit): input frequencies. qubit_map (dict): qubit map. Returns: @@ -572,7 +572,7 @@ def _exp_from_circuit(self, circuit: dict, qubit_map: dict, nshots: int = None): """ from qibo import Circuit, gates - if nshots is None: + if nshots is None: # pragma: no cover nshots = 1000 rotated_circuits = [] coefficients = [] From 84191d9652b210be786563ad4ecdc565125ae937 Mon Sep 17 00:00:00 2001 From: BrunoLiegiBastonLiegi Date: Fri, 18 Oct 2024 15:46:43 +0200 Subject: [PATCH 13/20] doc: added example to advanced code examples --- doc/final_result.npy | Bin 0 -> 937 bytes doc/source/code-examples/advancedexamples.rst | 23 +++++++++++++++--- src/qibo/hamiltonians/hamiltonians.py | 4 ++- 3 files changed, 22 insertions(+), 5 deletions(-) create mode 100644 doc/final_result.npy diff --git a/doc/final_result.npy b/doc/final_result.npy new file mode 100644 index 0000000000000000000000000000000000000000..ace5306f742633936a56a1bd0fa34f4330f32812 GIT binary patch literal 937 zcmbtS%Wl&^6m^=k1=I3=Kh5xP0>MhUKvjiQq)04?nj%X`rYN$;&ZM^R>&}cpQ5uO& zRkG$MS+ZpXcWg&$RdvN(j6L_9dpzfie-1wGeE3w;zG_o*#DveyEyFz6HQR*Zm}qvx()WoLO-m$&hkRC%Vg-pC&@`+?P6N zB}-36B?`5^b_h$ong{hGXxOma?AKZqX(T5Z(_zI~(ll*%6ISOxj^LaP=bQbyQ^Vc? z=!3!F{a>tN>2FK8;LV4yHiPx1vxbzIQIT^NvqTEGID|{iAM8#^LSsg@3{v5wZICbt zrAzs~!2ah>+vpA#%`e}o_LOk89}3C1TL~n26d}z|qljigNsh~SiAP0PNRd6IzDJZ? z)QOI@LK!jvR{||4F5+r{w{Xpd>v#+O#v$AY%J9dVa8reU%ZA&F z@ZZ3lli_O;?t1f{f$~_tx69}074Av6?>Vc8RL6q{17`)di3n1>riVMf>drFmiyr_U z6({b6+|NU47v8{d8i+&arSP~I?N}RDcLPWj kRJK|+yg*yQX+tm3&|9Qpoo%SZ&t9P|;gy$D36^L515r`}CIA2c literal 0 HcmV?d00001 diff --git a/doc/source/code-examples/advancedexamples.rst b/doc/source/code-examples/advancedexamples.rst index 572342f79f..3980586a81 100644 --- a/doc/source/code-examples/advancedexamples.rst +++ b/doc/source/code-examples/advancedexamples.rst @@ -2015,7 +2015,7 @@ to allow calculation of expectation values directly from such samples: .. testcode:: from qibo import Circuit, gates - from qibo.hamiltonians import Z + from qibo.hamiltonians import Z as HZ circuit = Circuit(4) circuit.add(gates.H(i) for i in range(4)) @@ -2024,7 +2024,7 @@ to allow calculation of expectation values directly from such samples: circuit.add(gates.CNOT(2, 3)) circuit.add(gates.M(*range(4))) - hamiltonian = Z(4) + hamiltonian = HZ(4) result = circuit(nshots=1024) expectation_value = hamiltonian.expectation_from_samples(result.frequencies()) @@ -2039,8 +2039,23 @@ This can also be invoked directly from the ``result`` object: expectation_value = result.expectation_from_samples(hamiltonian) -The expectation from samples currently works only for Hamiltonians that are diagonal in -the computational basis. +For Hamiltonians that are not diagonal in the computational basis, or that are sum of terms that cannot be +diagonalised simultaneously, one has to calculate the expectation value starting from the circuit: + +.. testcode:: + + from qibo.symbols import X, Y, Z + from qibo.hamiltonians import SymbolicHamiltonian + + hamiltonian = SymbolicHamiltonian(3 * Z(4) * X(1) ** 2 - (Y(0) * X(3)) / 2, nqubits=4) + expectation_value = hamiltonian.expectation_from_samples(circuit) + +What is happening under the hood in this case, is that the expectation value is calculated for each term +individually by measuring the circuit in the correct (rotated) basis. All the contributions are then +summed to recover the global expectation value. This means, in particular, that several copies of the +circuit are parallely executed, one for each term of the Hamiltonian. Note that, at the moment, no +check is performed to verify whether a subset of the terms could be diagonalised simultaneously, but +rather each term is treated separately every time. .. _tutorials_transpiler: diff --git a/src/qibo/hamiltonians/hamiltonians.py b/src/qibo/hamiltonians/hamiltonians.py index 85ec563e68..afd9db0ecd 100644 --- a/src/qibo/hamiltonians/hamiltonians.py +++ b/src/qibo/hamiltonians/hamiltonians.py @@ -556,7 +556,9 @@ def calculate_dense(self): def expectation(self, state, normalize=False): return Hamiltonian.expectation(self, state, normalize) - def _exp_from_circuit(self, circuit: Circuit, qubit_map: dict, nshots: int = None): + def _exp_from_circuit( + self, circuit: "Circuit", qubit_map: dict, nshots: int = None + ): """ Calculate the expectation value from a circuit. This allows for non diagonal observables. Each term of the observable is From 1ddb2c25480d63692b36fbd58cd8b81989d48d89 Mon Sep 17 00:00:00 2001 From: BrunoLiegiBastonLiegi Date: Tue, 22 Oct 2024 20:29:44 +0200 Subject: [PATCH 14/20] doc: small fix --- doc/final_result.npy | Bin 937 -> 937 bytes doc/source/code-examples/advancedexamples.rst | 12 +++++++++++- 2 files changed, 11 insertions(+), 1 deletion(-) diff --git a/doc/final_result.npy b/doc/final_result.npy index ace5306f742633936a56a1bd0fa34f4330f32812..fcfb33b7e4d2ecc664e35ac4cf5bebc772c5a2f8 100644 GIT binary patch delta 13 UcmZ3I|I|fWH#o-02=cHivR!s diff --git a/doc/source/code-examples/advancedexamples.rst b/doc/source/code-examples/advancedexamples.rst index 3980586a81..1c3b65c8d3 100644 --- a/doc/source/code-examples/advancedexamples.rst +++ b/doc/source/code-examples/advancedexamples.rst @@ -2047,7 +2047,17 @@ diagonalised simultaneously, one has to calculate the expectation value starting from qibo.symbols import X, Y, Z from qibo.hamiltonians import SymbolicHamiltonian - hamiltonian = SymbolicHamiltonian(3 * Z(4) * X(1) ** 2 - (Y(0) * X(3)) / 2, nqubits=4) + # build the circuit as before + circuit = Circuit(4) + circuit.add(gates.H(i) for i in range(4)) + circuit.add(gates.CNOT(0, 1)) + circuit.add(gates.CNOT(1, 2)) + circuit.add(gates.CNOT(2, 3)) + # but don't add any measurement at the end! + # they will be automatically added with the proper basis + # while calculating the expectation value + + hamiltonian = SymbolicHamiltonian(3 * Z(2) * X(1) ** 2 - (Y(0) * X(3)) / 2, nqubits=4) expectation_value = hamiltonian.expectation_from_samples(circuit) What is happening under the hood in this case, is that the expectation value is calculated for each term From e8480b50344e95847e3eb63e24bfdafe246a54bb Mon Sep 17 00:00:00 2001 From: BrunoLiegiBastonLiegi Date: Sat, 26 Oct 2024 10:11:49 +0200 Subject: [PATCH 15/20] fix: various fixes to make powers work properly --- src/qibo/hamiltonians/hamiltonians.py | 67 ++++++++++++++++----------- src/qibo/hamiltonians/terms.py | 14 +++++- tests/test_hamiltonians.py | 6 +-- 3 files changed, 55 insertions(+), 32 deletions(-) diff --git a/src/qibo/hamiltonians/hamiltonians.py b/src/qibo/hamiltonians/hamiltonians.py index afd9db0ecd..4141048a1e 100644 --- a/src/qibo/hamiltonians/hamiltonians.py +++ b/src/qibo/hamiltonians/hamiltonians.py @@ -578,13 +578,20 @@ def _exp_from_circuit( nshots = 1000 rotated_circuits = [] coefficients = [] + Z_observables = [] + q_maps = [] for term in self.terms: + # store coefficient coefficients.append(term.coefficient) - Z_observable = SymbolicHamiltonian( - prod([Z(q) for q in term.target_qubits]), - nqubits=circuit.nqubits, - backend=self.backend, + # build diagonal observable + Z_observables.append( + SymbolicHamiltonian( + prod([Z(q) for q in term.target_qubits]), + nqubits=circuit.nqubits, + backend=self.backend, + ) ) + # prepare the measurement basis and append it to the circuit measurements = [ gates.M(factor.target_qubit, basis=factor.gate.__class__) for factor in set(term.factors) @@ -592,18 +599,29 @@ def _exp_from_circuit( circ_copy = circuit.copy(True) [circ_copy.add(m) for m in measurements] rotated_circuits.append(circ_copy) + # map the original qubits into 0, 1, 2, ... + q_maps.append( + dict( + zip( + [qubit_map[q] for q in term.target_qubits], + range(len(term.target_qubits)), + ) + ) + ) frequencies = [ result.frequencies() for result in self.backend.execute_circuits(rotated_circuits, nshots=nshots) ] return sum( [ - c * Z_observable._exp_from_freq(f, qubit_map) - for c, f in zip(coefficients, frequencies) + coeff * obs._exp_from_freq(freq, qmap) + for coeff, freq, obs, qmap in zip( + coefficients, frequencies, Z_observables, q_maps + ) ] ) - def _exp_from_freq(self, freq: dict, qubit_map: dict): + def _exp_from_freq(self, freq: dict, qubit_map: dict) -> float: """ Calculate the expectation value from some frequencies. The observable has to be diagonal in the computational basis. @@ -627,26 +645,20 @@ def _exp_from_freq(self, freq: dict, qubit_map: dict): counts = self.backend.cast(list(freq.values()), self.backend.precision) / sum( freq.values() ) - qubits = [] + expvals = [] for term in self.terms: - qubits_term = [] - for k in term.target_qubits: - qubits_term.append(k) - qubits.append(qubits_term) - if qubit_map is None: - qubit_map = list(range(len(keys[0]))) - expval = 0 - for j, q in enumerate(qubits): - subk = [] - expval_q = 0 - for i, k in enumerate(keys): - subk = [int(k[qubit_map.index(s)]) for s in q] - expval_k = 1 - if subk.count(1) % 2 == 1: - expval_k = -1 - expval_q += expval_k * counts[i] - expval += expval_q * self.terms[j].coefficient.real - return expval + self.constant.real + qubits = term.target_qubits + expvals.extend( + [ + term.coefficient.real + * (-1) ** [state[qubit_map[q]] for q in qubits].count("1") + for state in keys + ] + ) + expvals = self.backend.cast(expvals, dtype=counts.dtype).reshape( + len(self.terms), len(freq) + ) + return self.backend.np.sum(expvals @ counts.T) + self.constant.real def expectation_from_samples( self, data: Union[dict, "Circuit"], qubit_map: dict = None, nshots: int = None @@ -674,7 +686,8 @@ def expectation_from_samples( Returns: (float) the computed expectation value. """ - + if qubit_map is None: + qubit_map = list(range(self.nqubits)) if isinstance(data, dict): return self._exp_from_freq(data, qubit_map) return self._exp_from_circuit(data, qubit_map, nshots) diff --git a/src/qibo/hamiltonians/terms.py b/src/qibo/hamiltonians/terms.py index f1d2f7cda0..54483d2687 100644 --- a/src/qibo/hamiltonians/terms.py +++ b/src/qibo/hamiltonians/terms.py @@ -1,8 +1,9 @@ import numpy as np import sympy -from qibo import gates +from qibo import gates, symbols from qibo.config import raise_error +from qibo.symbols import I class HamiltonianTerm: @@ -152,7 +153,16 @@ def __init__(self, coefficient, factors=1, symbol_map={}): factor, pow = factor.args assert isinstance(pow, sympy.Integer) assert isinstance(factor, sympy.Symbol) - pow = int(pow) + # if the symbol is a Pauli (i.e. a qibo symbol) and `pow` is even + # the power is the identity, thus the factor vanishes. Otherwise, + # for an odd exponent, it remains unchanged (i.e. `pow`=1) + if not factor in symbol_map: + if not int(pow) % 2: + factor = sympy.N(1) + else: + pow = 1 + else: + pow = int(pow) else: pow = 1 diff --git a/tests/test_hamiltonians.py b/tests/test_hamiltonians.py index 4bdf6cbd03..95129b0963 100644 --- a/tests/test_hamiltonians.py +++ b/tests/test_hamiltonians.py @@ -265,8 +265,8 @@ def test_hamiltonian_expectation_errors(backend): @pytest.mark.parametrize( "observable", [ - 2 * Z(0) * Z(1) ** 2 + Z(0) * Z(2), - X(0) * Z(1) + Y(0) * X(2) / 2 - Z(0) * Y(1) ** 3, + 2 * Z(0) * (1 - Z(1)) ** 2 + Z(0) * Z(2), + X(0) * Z(1) + Y(0) * X(2) / 2 - Z(0) * (1 - Y(1)) ** 3, ], ) def test_hamiltonian_expectation_from_samples(backend, observable): @@ -277,7 +277,7 @@ def test_hamiltonian_expectation_from_samples(backend, observable): nshots = 4 * 10**6 c = Circuit(nqubits) for q in range(nqubits): - c.add(gates.RX(0, np.random.rand())) + c.add(gates.RX(q, np.random.rand())) H = hamiltonians.SymbolicHamiltonian(observable, nqubits=nqubits, backend=backend) matrix = H.matrix From 0dfb4a87be62e2b650c61fbb7bb17507daccabc4 Mon Sep 17 00:00:00 2001 From: BrunoLiegiBastonLiegi Date: Sat, 26 Oct 2024 11:22:30 +0200 Subject: [PATCH 16/20] fix: small fixes for tests --- src/qibo/hamiltonians/terms.py | 4 ++-- tests/test_hamiltonians.py | 4 ++-- tests/test_hamiltonians_terms.py | 6 ++---- 3 files changed, 6 insertions(+), 8 deletions(-) diff --git a/src/qibo/hamiltonians/terms.py b/src/qibo/hamiltonians/terms.py index 54483d2687..c46557abeb 100644 --- a/src/qibo/hamiltonians/terms.py +++ b/src/qibo/hamiltonians/terms.py @@ -3,7 +3,7 @@ from qibo import gates, symbols from qibo.config import raise_error -from qibo.symbols import I +from qibo.symbols import I, X, Y, Z class HamiltonianTerm: @@ -156,7 +156,7 @@ def __init__(self, coefficient, factors=1, symbol_map={}): # if the symbol is a Pauli (i.e. a qibo symbol) and `pow` is even # the power is the identity, thus the factor vanishes. Otherwise, # for an odd exponent, it remains unchanged (i.e. `pow`=1) - if not factor in symbol_map: + if factor.__class__ in (I, X, Y, Z): if not int(pow) % 2: factor = sympy.N(1) else: diff --git a/tests/test_hamiltonians.py b/tests/test_hamiltonians.py index 95129b0963..bcb5e60f30 100644 --- a/tests/test_hamiltonians.py +++ b/tests/test_hamiltonians.py @@ -447,8 +447,8 @@ def construct_hamiltonian(): H1 = construct_hamiltonian() _ = H1.eigenvectors() - backend.assert_allclose(H.exp(0.5), target_matrix) - backend.assert_allclose(H1.exp(0.5), target_matrix) + backend.assert_allclose(H.exp(0.5), target_matrix, atol=1e-6) + backend.assert_allclose(H1.exp(0.5), target_matrix, atol=1e-6) def test_hamiltonian_energy_fluctuation(backend): diff --git a/tests/test_hamiltonians_terms.py b/tests/test_hamiltonians_terms.py index 8cf6f10771..6d3dafdbad 100644 --- a/tests/test_hamiltonians_terms.py +++ b/tests/test_hamiltonians_terms.py @@ -125,11 +125,9 @@ def test_symbolic_term_with_power_creation(backend): expression = X(0) ** 4 * Z(1) ** 2 * X(2) term = terms.SymbolicTerm(2, expression) - assert term.target_qubits == (0, 1, 2) - assert len(term.matrix_map) == 3 + assert term.target_qubits == (2,) + assert len(term.matrix_map) == 1 assert term.coefficient == 2 - backend.assert_allclose(term.matrix_map.get(0), 4 * [matrices.X]) - backend.assert_allclose(term.matrix_map.get(1), 2 * [matrices.Z]) backend.assert_allclose(term.matrix_map.get(2), [matrices.X]) From 9a10b11f85069468aac5e4d6623764b9456db67f Mon Sep 17 00:00:00 2001 From: BrunoLiegiBastonLiegi Date: Sun, 27 Oct 2024 09:59:04 +0100 Subject: [PATCH 17/20] doc: changed example observable --- doc/source/code-examples/advancedexamples.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/doc/source/code-examples/advancedexamples.rst b/doc/source/code-examples/advancedexamples.rst index 1c3b65c8d3..ffc9050060 100644 --- a/doc/source/code-examples/advancedexamples.rst +++ b/doc/source/code-examples/advancedexamples.rst @@ -2057,7 +2057,7 @@ diagonalised simultaneously, one has to calculate the expectation value starting # they will be automatically added with the proper basis # while calculating the expectation value - hamiltonian = SymbolicHamiltonian(3 * Z(2) * X(1) ** 2 - (Y(0) * X(3)) / 2, nqubits=4) + hamiltonian = SymbolicHamiltonian(3 * Z(2) * (1 - X(1)) ** 2 - (Y(0) * X(3)) / 2, nqubits=4) expectation_value = hamiltonian.expectation_from_samples(circuit) What is happening under the hood in this case, is that the expectation value is calculated for each term From e2d1c0511a5fd16e4243933c682b3ff5be90dd74 Mon Sep 17 00:00:00 2001 From: BrunoLiegiBastonLiegi <45011234+BrunoLiegiBastonLiegi@users.noreply.github.com> Date: Tue, 29 Oct 2024 08:04:52 +0100 Subject: [PATCH 18/20] Update src/qibo/hamiltonians/hamiltonians.py Co-authored-by: Stavros Efthymiou <35475381+stavros11@users.noreply.github.com> --- src/qibo/hamiltonians/hamiltonians.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/qibo/hamiltonians/hamiltonians.py b/src/qibo/hamiltonians/hamiltonians.py index 4141048a1e..a1abcbb72e 100644 --- a/src/qibo/hamiltonians/hamiltonians.py +++ b/src/qibo/hamiltonians/hamiltonians.py @@ -597,7 +597,7 @@ def _exp_from_circuit( for factor in set(term.factors) ] circ_copy = circuit.copy(True) - [circ_copy.add(m) for m in measurements] + circ_copy.add(measurements) rotated_circuits.append(circ_copy) # map the original qubits into 0, 1, 2, ... q_maps.append( From ea570bf8a7b52025e49a726a6df5beea62b9273b Mon Sep 17 00:00:00 2001 From: BrunoLiegiBastonLiegi Date: Tue, 29 Oct 2024 09:07:40 +0100 Subject: [PATCH 19/20] feat: changed Z import in docs --- doc/final_result.npy | Bin 937 -> 0 bytes doc/source/code-examples/advancedexamples.rst | 4 ++-- 2 files changed, 2 insertions(+), 2 deletions(-) delete mode 100644 doc/final_result.npy diff --git a/doc/final_result.npy b/doc/final_result.npy deleted file mode 100644 index fcfb33b7e4d2ecc664e35ac4cf5bebc772c5a2f8..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 937 zcmbtS%Wl&^6m^=k1=I3=KaF`ffnX&dp{hbEQY02cO_3!eQxsWaXHt*w>&}cpQ5uO& zRkG$MS+ZpXcWg&$RdvN(j6L_9dpzfie-1uwfB00>zG_ozL`A@@O~cyXvD$_;PPybH z@lrmbs{M{cf?{8UB%`=*w%W$b_`lM>*MlUFvxyy~oZ4|7Nk%v)6W!x9NE0D>9!MRu zlBFl35`|h{JAkEL&4>CSG#prN_G_(*G?J5y>agN2X_~gP0ju*Lhj7M$v(0|ptzmB; z^ub{8{x4Rs^tU9O^XEfYo56b1T|-JtiO4yPX(9!jAHoIq4|b< zlOExF0{fpkZKFF}G{1bS+EbI$JtidOUL}y^QG~Q0jUti>B{?qRHGLvjA!Y6=^<7iR zWgTQxGg-HKpCRB$`~;OAJ?LE Date: Tue, 29 Oct 2024 09:47:44 +0100 Subject: [PATCH 20/20] feat: separated the two expectation functions --- doc/source/code-examples/advancedexamples.rst | 2 +- src/qibo/hamiltonians/hamiltonians.py | 73 +++++++------------ tests/test_hamiltonians.py | 46 ++++++------ 3 files changed, 51 insertions(+), 70 deletions(-) diff --git a/doc/source/code-examples/advancedexamples.rst b/doc/source/code-examples/advancedexamples.rst index 623c6c1ce0..bbd592f820 100644 --- a/doc/source/code-examples/advancedexamples.rst +++ b/doc/source/code-examples/advancedexamples.rst @@ -2058,7 +2058,7 @@ diagonalised simultaneously, one has to calculate the expectation value starting # while calculating the expectation value hamiltonian = SymbolicHamiltonian(3 * Z(2) * (1 - X(1)) ** 2 - (Y(0) * X(3)) / 2, nqubits=4) - expectation_value = hamiltonian.expectation_from_samples(circuit) + expectation_value = hamiltonian.expectation_from_circuit(circuit) What is happening under the hood in this case, is that the expectation value is calculated for each term individually by measuring the circuit in the correct (rotated) basis. All the contributions are then diff --git a/src/qibo/hamiltonians/hamiltonians.py b/src/qibo/hamiltonians/hamiltonians.py index a1abcbb72e..710a7b528c 100644 --- a/src/qibo/hamiltonians/hamiltonians.py +++ b/src/qibo/hamiltonians/hamiltonians.py @@ -2,7 +2,7 @@ from itertools import chain from math import prod -from typing import Optional, Union +from typing import Optional import numpy as np import sympy @@ -556,26 +556,34 @@ def calculate_dense(self): def expectation(self, state, normalize=False): return Hamiltonian.expectation(self, state, normalize) - def _exp_from_circuit( - self, circuit: "Circuit", qubit_map: dict, nshots: int = None - ): + def expectation_from_circuit( + self, circuit: "Circuit", qubit_map: dict = None, nshots: int = 1000 + ) -> float: """ Calculate the expectation value from a circuit. - This allows for non diagonal observables. Each term of the observable is - treated separately, by measuring in the correct basis and re-executing the - circuit. + This even works for observables not completely diagonal in the computational + basis, but only diagonal at a term level in a defined basis. Namely, for + an observable of the form :math:``H = \\sum_i H_i``, where each :math:``H_i`` + consists in a `n`-qubits pauli operator :math:`P_0 \\otimes P_1 \\otimes \\cdots \\otimes P_n`, + the expectation value is computed by rotating the input circuit in the suitable + basis for each term :math:``H_i`` thus extracting the `term-wise` expectations + that are then summed to build the global expectation value. + Each term of the observable is treated separately, by measuring in the correct + basis and re-executing the circuit. Args: - circuit (Circuit): input frequencies. - qubit_map (dict): qubit map. + circuit (Circuit): input circuit. + qubit_map (dict): qubit map, defaults to ``{0: 0, 1: 1, ..., n: n}`` + nshots (int): number of shots, defaults to 1000. Returns: (float): the calculated expectation value. """ from qibo import Circuit, gates - if nshots is None: # pragma: no cover - nshots = 1000 + if qubit_map is None: + qubit_map = list(range(self.nqubits)) + rotated_circuits = [] coefficients = [] Z_observables = [] @@ -614,20 +622,20 @@ def _exp_from_circuit( ] return sum( [ - coeff * obs._exp_from_freq(freq, qmap) + coeff * obs.expectation_from_samples(freq, qmap) for coeff, freq, obs, qmap in zip( coefficients, frequencies, Z_observables, q_maps ) ] ) - def _exp_from_freq(self, freq: dict, qubit_map: dict) -> float: + def expectation_from_samples(self, freq: dict, qubit_map: dict = None) -> float: """ - Calculate the expectation value from some frequencies. + Calculate the expectation value from the samples. The observable has to be diagonal in the computational basis. Args: - freq (dict): input frequencies. + freq (dict): input frequencies of the samples. qubit_map (dict): qubit map. Returns: @@ -641,6 +649,9 @@ def _exp_from_freq(self, freq: dict, qubit_map: dict) -> float: NotImplementedError, "Observable is not a Z Pauli string." ) + if qubit_map is None: + qubit_map = list(range(self.nqubits)) + keys = list(freq.keys()) counts = self.backend.cast(list(freq.values()), self.backend.precision) / sum( freq.values() @@ -660,38 +671,6 @@ def _exp_from_freq(self, freq: dict, qubit_map: dict) -> float: ) return self.backend.np.sum(expvals @ counts.T) + self.constant.real - def expectation_from_samples( - self, data: Union[dict, "Circuit"], qubit_map: dict = None, nshots: int = None - ): - """ - Calculate the expectation value starting from the samples. This even works - for observables not completely diagonal in the computational basis, but only - diagonal at a term level in a defined basis. Namely, for an observable of the - form :math:``H = \\sum_i H_i``, where each :math:``H_i`` consists in a `n`-qubits - pauli operator :math:`P_0 \\otimes P_1 \\otimes \\cdots \\otimes P_n`, the - expectation value is computed by rotating the input circuit in the suitable - basis for each term :math:``H_i`` thus extracting the `term-wise` expectations - that are then summed to build the global expectation value. - - Args: - data (dict | :class:`qibo.models.Circuit`): either the ``dict`` of the - frequencies of the samples or a :class:`qibo.models.Circuit` where to - extract the samples from. A :class:`qibo.models.Circuit` is needed in - case the observable is not diagonal in the computational basis. - qubit_map (dict): qubit map. - nshots (int, optional): number of shots for the expecation value - calculation. This is used only when a :class:`qibo.models.Circuit` - is passed as input. - - Returns: - (float) the computed expectation value. - """ - if qubit_map is None: - qubit_map = list(range(self.nqubits)) - if isinstance(data, dict): - return self._exp_from_freq(data, qubit_map) - return self._exp_from_circuit(data, qubit_map, nshots) - def __add__(self, o): if isinstance(o, self.__class__): if self.nqubits != o.nqubits: diff --git a/tests/test_hamiltonians.py b/tests/test_hamiltonians.py index bcb5e60f30..5d2e5d57eb 100644 --- a/tests/test_hamiltonians.py +++ b/tests/test_hamiltonians.py @@ -262,38 +262,40 @@ def test_hamiltonian_expectation_errors(backend): h.expectation("test") -@pytest.mark.parametrize( - "observable", - [ - 2 * Z(0) * (1 - Z(1)) ** 2 + Z(0) * Z(2), - X(0) * Z(1) + Y(0) * X(2) / 2 - Z(0) * (1 - Y(1)) ** 3, - ], -) -def test_hamiltonian_expectation_from_samples(backend, observable): - """Test Hamiltonian expectation value calculation.""" - backend.set_seed(12) +def non_exact_expectation_test_setup(backend, observable): nqubits = 3 - nshots = 4 * 10**6 c = Circuit(nqubits) for q in range(nqubits): c.add(gates.RX(q, np.random.rand())) H = hamiltonians.SymbolicHamiltonian(observable, nqubits=nqubits, backend=backend) - matrix = H.matrix - diagonal = ( - backend.np.count_nonzero(matrix - backend.np.diag(backend.np.diagonal(matrix))) - == 0 - ) final_state = backend.execute_circuit(c.copy(True)).state() exp = H.expectation(final_state) - if not diagonal: - exp_from_samples = H.expectation_from_samples(c, nshots=nshots) - else: - c.add(gates.M(*range(nqubits))) - freq = backend.execute_circuit(c, nshots=nshots).frequencies() - exp_from_samples = H.expectation_from_samples(freq) + return exp, H, c + + +def test_hamiltonian_expectation_from_samples(backend): + """Test Hamiltonian expectation value calculation.""" + backend.set_seed(12) + + nshots = 4 * 10**6 + observable = 2 * Z(0) * (1 - Z(1)) ** 2 + Z(0) * Z(2) + exp, H, c = non_exact_expectation_test_setup(backend, observable) + c.add(gates.M(*range(c.nqubits))) + freq = backend.execute_circuit(c, nshots=nshots).frequencies() + exp_from_samples = H.expectation_from_samples(freq) + backend.assert_allclose(exp, exp_from_samples, atol=1e-2) + +def test_hamiltonian_expectation_from_circuit(backend): + """Test Hamiltonian expectation value calculation.""" + backend.set_seed(12) + + nshots = 4 * 10**6 + observable = X(0) * Z(1) + Y(0) * X(2) / 2 - Z(0) * (1 - Y(1)) ** 3 + exp, H, c = non_exact_expectation_test_setup(backend, observable) + exp_from_samples = H.expectation_from_circuit(c, nshots=nshots) backend.assert_allclose(exp, exp_from_samples, atol=1e-2)