From 037ae6375e07583b53f396594f1333d29055bfad Mon Sep 17 00:00:00 2001 From: Declan Millar Date: Tue, 27 Sep 2022 14:43:40 +0100 Subject: [PATCH] VQE implementation with estimator primitive (#8702) * Add minimally working VQE with estimator primitive implementation. No gradients, no tests etc. * Revert from dataclass results to original result classes. * Enforce positional and keyword VQE arguments. * Move aux op eval logic to function * Update docstring. * Remove max_evals_grouped. Force to set directly on optimizer. * Remove validate min import. * Make note that eval_observables will be used to eval aux ops. * Add initial vqe tests. * Have VQE inherit from VariationalAlgorithm. * Move energy evaluation to unnested function. * Construct h2_op using SparsePauliOp * Add gradient with primitives support. * Update docstrings * update broadcast handling * update eval_energy output for batching * add incomplete QNSPA test * fix batch evaluation of QNSPSA test * remove vqe callback * move estimator to first arg * remove usused imports * add minimum eigensolvers test init file * add aux ops tests and prepare for new eval_operators * no longer support account for Nones in aux_ops * correct typing for MinimumEigensolver * Compute default initial point using ansatz bounds. * Add NumPyMinimumEigensolverResult * Fix type hints * Fix type hints * Formatting * Do not store NumPyMES result inside the algo. * Provide default values for ansatz and estimator * Formatting * fix old and new batching * Add tests for NumpyMES and import in module. * Use lazy formatting in log messages * Use lazy formatting in log messages * Add back callback to VQE. * minor renaming * raise algorithm error if primitive jobs fail * Add return documentation * Improve var names and docstrings. * Apply suggestions from code review Co-authored-by: ElePT <57907331+ElePT@users.noreply.github.com> Co-authored-by: dlasecki * minor formatting * minor formatting * Ensure evaluate_energy/gradient match Co-authored-by: Max Rossmannek * return bounds logic; fix some docstrings and minor refactor * Force keyword arguments in vqe. * Use estimate_observables function * break up numpy mes tests with subTest * formatting * remove redundant eval_aux_ops * add typehints to docstring attributes * remove usused imports * remove default ansatz * remove usused imports * update typehints * avoid changing the original ansatz * avoid changing the original ansatz * create separate function to build vqe result * Correct aux operator eignvalue type hint Co-authored-by: ElePT <57907331+ElePT@users.noreply.github.com> * add variance to callback * use std_dev in callback rather than variance * formatting * return variance and shots in callback * return full metadata in callback * Move validation functions to algorithms/utils * correct the callback attribute documentation * correct the callback attribute typehint docstring * update VQE docstring * release note and pending-depreciate old algs * update vqe class docstring * Apply suggestions from code review Co-authored-by: ElePT <57907331+ElePT@users.noreply.github.com> * Do not copy ansatz * Note pending depreciation of old algs * fix docstrings and imports * Fix issues with building docs * Include OptimizerResult in VQEResult * Remove trailing whitespace * Fix math notation in docstring * estimate_obervables to return metadata @ElePT +VQE * add example in release note * Update evaluate_observables docstring * Fix observables_evaluator tests. Co-authored-by: ElePT <57907331+ElePT@users.noreply.github.com> * fix trotter_qrte tests and remove depreciation * formatting * remove unused import * Apply suggestions to docstring from code review Co-authored-by: Steve Wood <40241007+woodsp-ibm@users.noreply.github.com> * Update VQE docstring * Remove printing Co-authored-by: ElePT <57907331+ElePT@users.noreply.github.com> * Add parameters to estimate observables * Fix estimate obs. unit test * Update arg description * Update arg description * keep equation part of sentence * dict -> dict[str, Any] * Update qiskit/algorithms/optimizers/spsa.py Apply suggestion Imamichi-san * introduce FilterType and aux_operator_eigenvalues -> aux_operators_evaluated * Correct typehint Co-authored-by: Ikko Hamamura * update vqe docstring and use old typing for type alias Co-authored-by: ElePT <57907331+ElePT@users.noreply.github.com> Co-authored-by: dlasecki Co-authored-by: Max Rossmannek Co-authored-by: Steve Wood <40241007+woodsp-ibm@users.noreply.github.com> Co-authored-by: ElePT Co-authored-by: Ikko Hamamura --- qiskit/algorithms/__init__.py | 16 +- .../minimum_eigensolvers/__init__.py | 48 ++ .../minimum_eigensolver.py | 99 ++++ .../numpy_minimum_eigensolver.py | 106 +++++ qiskit/algorithms/minimum_eigensolvers/vqe.py | 341 ++++++++++++++ qiskit/algorithms/observables_evaluator.py | 81 +--- qiskit/algorithms/optimizers/qnspsa.py | 18 +- qiskit/algorithms/optimizers/spsa.py | 25 +- .../trotterization/trotter_qrte.py | 1 + qiskit/algorithms/utils/__init__.py | 21 + qiskit/algorithms/utils/validate_bounds.py | 44 ++ .../utils/validate_initial_point.py | 69 +++ qiskit/algorithms/variational_algorithm.py | 11 + ...-and-setters-for-vqe-edc753591b368980.yaml | 6 +- ...ct-for-aux-operators-c3c9ad380c208afd.yaml | 6 +- ...r_empty_operator_fix-53ce20e5d2b68fd6.yaml | 16 +- ...-estimator-primitive-7cbcc462ad4dc593.yaml | 34 ++ .../minimum_eigensolvers/__init__.py | 11 + .../test_numpy_minimum_eigensolver.py | 241 ++++++++++ .../minimum_eigensolvers/test_vqe.py | 444 ++++++++++++++++++ .../algorithms/test_observables_evaluator.py | 36 +- .../time_evolvers/test_trotter_qrte.py | 10 +- test/python/algorithms/utils/__init__.py | 11 + .../algorithms/utils/test_validate_bounds.py | 53 +++ .../utils/test_validate_initial_point.py | 50 ++ 25 files changed, 1715 insertions(+), 83 deletions(-) create mode 100644 qiskit/algorithms/minimum_eigensolvers/__init__.py create mode 100644 qiskit/algorithms/minimum_eigensolvers/minimum_eigensolver.py create mode 100644 qiskit/algorithms/minimum_eigensolvers/numpy_minimum_eigensolver.py create mode 100644 qiskit/algorithms/minimum_eigensolvers/vqe.py create mode 100644 qiskit/algorithms/utils/__init__.py create mode 100644 qiskit/algorithms/utils/validate_bounds.py create mode 100644 qiskit/algorithms/utils/validate_initial_point.py create mode 100644 releasenotes/notes/vqe-with-estimator-primitive-7cbcc462ad4dc593.yaml create mode 100644 test/python/algorithms/minimum_eigensolvers/__init__.py create mode 100644 test/python/algorithms/minimum_eigensolvers/test_numpy_minimum_eigensolver.py create mode 100644 test/python/algorithms/minimum_eigensolvers/test_vqe.py create mode 100644 test/python/algorithms/utils/__init__.py create mode 100644 test/python/algorithms/utils/test_validate_bounds.py create mode 100644 test/python/algorithms/utils/test_validate_initial_point.py diff --git a/qiskit/algorithms/__init__.py b/qiskit/algorithms/__init__.py index 322acba21d..55648a040a 100644 --- a/qiskit/algorithms/__init__.py +++ b/qiskit/algorithms/__init__.py @@ -183,10 +183,12 @@ Algorithms to solve linear systems of equations. linear_solvers -Minimum Eigensolvers --------------------- +Minimum Eigen Solvers +--------------------- Algorithms that can find the minimum eigenvalue of an operator. +These algorithms are pending depreciation. One should instead make use of the +Minimum Eigensolver classes in the section below, which leverage Runtime primitives. .. autosummary:: :toctree: ../stubs/ @@ -203,6 +205,16 @@ Algorithms that can find the minimum eigenvalue of an operator. QAOA VQE +Minimum Eigensolvers +-------------------- + +Algorithms that can find the minimum eigenvalue of an operator and leverage primitives. + +.. autosummary:: + :toctree: ../stubs/ + + minimum_eigensolvers + Optimizers ---------- diff --git a/qiskit/algorithms/minimum_eigensolvers/__init__.py b/qiskit/algorithms/minimum_eigensolvers/__init__.py new file mode 100644 index 0000000000..9da052fab6 --- /dev/null +++ b/qiskit/algorithms/minimum_eigensolvers/__init__.py @@ -0,0 +1,48 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. + +""" +============================================================================ +Minimum Eigensolvers Package (:mod:`qiskit.algorithms.minimum_eigensolvers`) +============================================================================ + +.. currentmodule:: qiskit.algorithms.minimum_eigensolvers + +Minimum Eigensolvers +==================== +.. autosummary:: + :toctree: ../stubs/ + + MinimumEigensolver + NumPyMinimumEigensolver + VQE + +.. autosummary:: + :toctree: ../stubs/ + + MinimumEigensolverResult + NumPyMinimumEigensolverResult + VQEResult +""" + +from .minimum_eigensolver import MinimumEigensolver, MinimumEigensolverResult +from .numpy_minimum_eigensolver import NumPyMinimumEigensolver, NumPyMinimumEigensolverResult +from .vqe import VQE, VQEResult + +__all__ = [ + "MinimumEigensolver", + "MinimumEigensolverResult", + "NumPyMinimumEigensolver", + "NumPyMinimumEigensolverResult", + "VQE", + "VQEResult", +] diff --git a/qiskit/algorithms/minimum_eigensolvers/minimum_eigensolver.py b/qiskit/algorithms/minimum_eigensolvers/minimum_eigensolver.py new file mode 100644 index 0000000000..30d3b49ce4 --- /dev/null +++ b/qiskit/algorithms/minimum_eigensolvers/minimum_eigensolver.py @@ -0,0 +1,99 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. + +"""The minimum eigensolver interface and result.""" + +from __future__ import annotations + +from abc import ABC, abstractmethod +from typing import Any + +from qiskit.opflow import PauliSumOp +from qiskit.quantum_info.operators.base_operator import BaseOperator + +from ..algorithm_result import AlgorithmResult +from ..list_or_dict import ListOrDict + + +class MinimumEigensolver(ABC): + """The minimum eigensolver interface. + + Algorithms that can compute a minimum eigenvalue for an operator may implement this interface to + allow different algorithms to be used interchangeably. + """ + + @abstractmethod + def compute_minimum_eigenvalue( + self, + operator: BaseOperator | PauliSumOp, + aux_operators: ListOrDict[BaseOperator | PauliSumOp] | None = None, + ) -> "MinimumEigensolverResult": + """ + Computes the minimum eigenvalue. The ``operator`` and ``aux_operators`` can be supplied here + and if not ``None`` will override any already set into algorithm so it can be reused with + different operators. While an ``operator`` is required by algorithms, ``aux_operators`` are + optional. + + Args: + operator: Qubit operator of the observable. + aux_operators: Optional list of auxiliary operators to be evaluated with the + parameters of the minimum eigenvalue main result and their expectation values + returned. For instance in chemistry these can be dipole operators and total particle + count operators, so we can get values for these at the ground state. + + Returns: + A minimum eigensolver result. + """ + return MinimumEigensolverResult() + + @classmethod + def supports_aux_operators(cls) -> bool: + """Whether computing the expectation value of auxiliary operators is supported. + + If the minimum eigensolver computes an eigenvalue of the main ``operator`` then it can + compute the expectation value of the ``aux_operators`` for that state. Otherwise they will + be ignored. + + Returns: + True if aux_operator expectations can be evaluated, False otherwise + """ + return False + + +class MinimumEigensolverResult(AlgorithmResult): + """Minimum eigensolver result.""" + + def __init__(self) -> None: + super().__init__() + self._eigenvalue = None + self._aux_operators_evaluated = None + + @property + def eigenvalue(self) -> complex | None: + """The computed minimum eigenvalue.""" + return self._eigenvalue + + @eigenvalue.setter + def eigenvalue(self, value: complex) -> None: + self._eigenvalue = value + + @property + def aux_operators_evaluated(self) -> ListOrDict[tuple[complex, dict[str, Any]]] | None: + """The aux operator expectation values. + + These values are in fact tuples formatted as (mean, (variance, shots)). + """ + return self._aux_operators_evaluated + + @aux_operators_evaluated.setter + def aux_operators_evaluated(self, value: ListOrDict[tuple[complex, dict[str, Any]]]) -> None: + self._aux_operators_evaluated = value diff --git a/qiskit/algorithms/minimum_eigensolvers/numpy_minimum_eigensolver.py b/qiskit/algorithms/minimum_eigensolvers/numpy_minimum_eigensolver.py new file mode 100644 index 0000000000..ec11f27eb2 --- /dev/null +++ b/qiskit/algorithms/minimum_eigensolvers/numpy_minimum_eigensolver.py @@ -0,0 +1,106 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. + +"""The NumPy minimum eigensolver algorithm and result.""" + +from __future__ import annotations + +from typing import Callable, List, Union, Optional +import logging +import numpy as np + +from qiskit.opflow import PauliSumOp +from qiskit.quantum_info.operators.base_operator import BaseOperator + +# TODO this path will need updating when VQD is merged +from ..eigen_solvers.numpy_eigen_solver import NumPyEigensolver +from .minimum_eigensolver import MinimumEigensolver, MinimumEigensolverResult +from ..list_or_dict import ListOrDict + +logger = logging.getLogger(__name__) + +# future type annotations not supported in type aliases in 3.8 +FilterType = Callable[[Union[List, np.ndarray], float, Optional[ListOrDict[float]]], bool] + + +class NumPyMinimumEigensolver(MinimumEigensolver): + """ + The NumPy minimum eigensolver algorithm. + """ + + def __init__( + self, + filter_criterion: FilterType | None = None, + ) -> None: + """ + Args: + filter_criterion: callable that allows to filter eigenvalues/eigenstates. The minimum + eigensolver is only searching over feasible states and returns an eigenstate that + has the smallest eigenvalue among feasible states. The callable has the signature + ``filter(eigenstate, eigenvalue, aux_values)`` and must return a boolean to indicate + whether to consider this value or not. If there is no feasible element, the result + can even be empty. + """ + self._eigensolver = NumPyEigensolver(filter_criterion=filter_criterion) + + @property + def filter_criterion( + self, + ) -> FilterType | None: + """Returns the criterion for filtering eigenstates/eigenvalues.""" + return self._eigensolver.filter_criterion + + @filter_criterion.setter + def filter_criterion( + self, + filter_criterion: FilterType, + ) -> None: + self._eigensolver.filter_criterion = filter_criterion + + @classmethod + def supports_aux_operators(cls) -> bool: + return NumPyEigensolver.supports_aux_operators() + + def compute_minimum_eigenvalue( + self, + operator: BaseOperator | PauliSumOp, + aux_operators: ListOrDict[BaseOperator | PauliSumOp] | None = None, + ) -> NumPyMinimumEigensolverResult: + super().compute_minimum_eigenvalue(operator, aux_operators) + eigensolver_result = self._eigensolver.compute_eigenvalues(operator, aux_operators) + result = NumPyMinimumEigensolverResult() + if eigensolver_result.eigenvalues is not None and len(eigensolver_result.eigenvalues) > 0: + result.eigenvalue = eigensolver_result.eigenvalues[0] + result.eigenstate = eigensolver_result.eigenstates[0] + if eigensolver_result.aux_operator_eigenvalues: + result.aux_operators_evaluated = eigensolver_result.aux_operator_eigenvalues[0] + + logger.debug("NumPy minimum eigensolver result: %s", result) + + return result + + +class NumPyMinimumEigensolverResult(MinimumEigensolverResult): + """NumPy minimum eigensolver result.""" + + def __init__(self) -> None: + super().__init__() + self._eigenstate = None + + @property + def eigenstate(self) -> np.ndarray | None: + """Returns the eigenstate corresponding to the computed minimum eigenvalue.""" + return self._eigenstate + + @eigenstate.setter + def eigenstate(self, value: np.ndarray) -> None: + self._eigenstate = value diff --git a/qiskit/algorithms/minimum_eigensolvers/vqe.py b/qiskit/algorithms/minimum_eigensolvers/vqe.py new file mode 100644 index 0000000000..7d4c76d981 --- /dev/null +++ b/qiskit/algorithms/minimum_eigensolvers/vqe.py @@ -0,0 +1,341 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. + +"""The variational quantum eigensolver algorithm.""" + +from __future__ import annotations + +import logging +from time import time +from collections.abc import Callable, Sequence +from typing import Any + +import numpy as np + +from qiskit.algorithms.gradients import BaseEstimatorGradient +from qiskit.circuit import QuantumCircuit +from qiskit.opflow import PauliSumOp +from qiskit.primitives import BaseEstimator +from qiskit.quantum_info.operators.base_operator import BaseOperator + +from ..exceptions import AlgorithmError +from ..list_or_dict import ListOrDict +from ..optimizers import Optimizer, Minimizer, OptimizerResult +from ..variational_algorithm import VariationalAlgorithm, VariationalResult +from .minimum_eigensolver import MinimumEigensolver, MinimumEigensolverResult +from ..observables_evaluator import estimate_observables +from ..utils import validate_initial_point, validate_bounds + +logger = logging.getLogger(__name__) + + +class VQE(VariationalAlgorithm, MinimumEigensolver): + r"""The variational quantum eigensolver (VQE) algorithm. + + VQE is a hybrid quantum-classical algorithm that uses a variational technique to find the + minimum eigenvalue of a given Hamiltonian operator :math:`H`. + + The ``VQE`` algorithm is executed using an :attr:`estimator` primitive, which computes + expectation values of operators (observables). + + An instance of ``VQE`` also requires an :attr:`ansatz`, a parameterized + :class:`.QuantumCircuit`, to prepare the trial state :math:`|\psi(\vec\theta)\rangle`. It also + needs a classical :attr:`optimizer` which varies the circuit parameters :math:`\vec\theta` such + that the expectation value of the operator on the corresponding state approaches a minimum, + + .. math:: + + \min_{\vec\theta} \langle\psi(\vec\theta)|H|\psi(\vec\theta)\rangle. + + The :attr:`estimator` is used to compute this expectation value for every optimization step. + + The optimizer can either be one of Qiskit's optimizers, such as + :class:`~qiskit.algorithms.optimizers.SPSA` or a callable with the following signature: + + .. code-block:: python + + from qiskit.algorithms.optimizers import OptimizerResult + + def my_minimizer(fun, x0, jac=None, bounds=None) -> OptimizerResult: + # Note that the callable *must* have these argument names! + # Args: + # fun (callable): the function to minimize + # x0 (np.ndarray): the initial point for the optimization + # jac (callable, optional): the gradient of the objective function + # bounds (list, optional): a list of tuples specifying the parameter bounds + + result = OptimizerResult() + result.x = # optimal parameters + result.fun = # optimal function value + return result + + The above signature also allows one to use any SciPy minimizer, for instance as + + .. code-block:: python + + from functools import partial + from scipy.optimize import minimize + + optimizer = partial(minimize, method="L-BFGS-B") + + The following attributes can be set via the initializer but can also be read and updated once + the VQE object has been constructed. + + Attributes: + estimator (BaseEstimator): The estimator primitive to compute the expectation value of the + Hamiltonian operator. + ansatz (QuantumCircuit): A parameterized quantum circuit to prepare the trial state. + optimizer (Optimizer | Minimizer): A classical optimizer to find the minimum energy. This + can either be a Qiskit :class:`.Optimizer` or a callable implementing the + :class:`.Minimizer` protocol. + gradient (BaseEstimatorGradient | None): An optional estimator gradient to be used with the + optimizer. + callback (Callable[[int, np.ndarray, float, dict[str, Any]], None] | None): A callback that + can access the intermediate data at each optimization step. These data are: the + evaluation count, the optimizer parameters for the ansatz, the evaluated mean, and the + metadata dictionary. + + References: + [1] Peruzzo et al, "A variational eigenvalue solver on a quantum processor" + `arXiv:1304.3061 `__ + """ + + def __init__( + self, + estimator: BaseEstimator, + ansatz: QuantumCircuit, + optimizer: Optimizer | Minimizer, + *, + gradient: BaseEstimatorGradient | None = None, + initial_point: Sequence[float] | None = None, + callback: Callable[[int, np.ndarray, float, dict[str, Any]], None] | None = None, + ) -> None: + r""" + Args: + estimator: The estimator primitive to compute the expectation value of the + Hamiltonian operator. + ansatz: A parameterized quantum circuit to prepare the trial state. + optimizer: A classical optimizer to find the minimum energy. This + can either be a Qiskit :class:`.Optimizer` or a callable implementing the + :class:`.Minimizer` protocol. + gradient: An optional estimator gradient to be used with the optimizer. + initial_point: An optional initial point (i.e. initial parameter values) for the + optimizer. The length of the initial point must match the number of :attr:`ansatz` + parameters. If ``None``, a random point will be generated within certain parameter + bounds. ``VQE`` will look to the ansatz for these bounds. If the ansatz does not + specify bounds, bounds of :math:`-2\pi`, :math:`2\pi` will be used. + callback: A callback that can access the intermediate data at each optimization step. + These data are: the evaluation count, the optimizer parameters for the ansatz, the + estimated value, and the metadata dictionary. + """ + super().__init__() + + self.estimator = estimator + self.ansatz = ansatz + self.optimizer = optimizer + self.gradient = gradient + # this has to go via getters and setters due to the VariationalAlgorithm interface + self.initial_point = initial_point + self.callback = callback + + @property + def initial_point(self) -> Sequence[float] | None: + return self._initial_point + + @initial_point.setter + def initial_point(self, value: Sequence[float] | None) -> None: + self._initial_point = value + + def compute_minimum_eigenvalue( + self, + operator: BaseOperator | PauliSumOp, + aux_operators: ListOrDict[BaseOperator | PauliSumOp] | None = None, + ) -> VQEResult: + self._check_operator_ansatz(operator) + + initial_point = validate_initial_point(self.initial_point, self.ansatz) + + bounds = validate_bounds(self.ansatz) + + start_time = time() + + evaluate_energy = self._get_evaluate_energy(self.ansatz, operator) + + if self.gradient is not None: + evaluate_gradient = self._get_evaluate_gradient(self.ansatz, operator) + else: + evaluate_gradient = None + + # perform optimization + if callable(self.optimizer): + optimizer_result = self.optimizer( + fun=evaluate_energy, x0=initial_point, jac=evaluate_gradient, bounds=bounds + ) + else: + optimizer_result = self.optimizer.minimize( + fun=evaluate_energy, x0=initial_point, jac=evaluate_gradient, bounds=bounds + ) + + optimizer_time = time() - start_time + + logger.info( + "Optimization complete in %s seconds.\nFound optimal point %s", + optimizer_time, + optimizer_result.x, + ) + + if aux_operators is not None: + aux_operators_evaluated = estimate_observables( + self.estimator, self.ansatz, aux_operators, optimizer_result.x + ) + else: + aux_operators_evaluated = None + + return self._build_vqe_result(optimizer_result, aux_operators_evaluated, optimizer_time) + + @classmethod + def supports_aux_operators(cls) -> bool: + return True + + def _get_evaluate_energy( + self, + ansatz: QuantumCircuit, + operator: BaseOperator | PauliSumOp, + ) -> Callable[[np.ndarray], np.ndarray | float]: + """Returns a function handle to evaluate the energy at given parameters for the ansatz. + This is the objective function to be passed to the optimizer that is used for evaluation. + + Args: + ansatz: The ansatz preparing the quantum state. + operator: The operator whose energy to evaluate. + + Returns: + A callable that computes and returns the energy of the hamiltonian of each parameter. + + Raises: + AlgorithmError: If the primitive job to evaluate the energy fails. + """ + num_parameters = ansatz.num_parameters + + # avoid creating an instance variable to remain stateless regarding results + eval_count = 0 + + def evaluate_energy(parameters: np.ndarray) -> np.ndarray | float: + nonlocal eval_count + + # handle broadcasting: ensure parameters is of shape [array, array, ...] + parameters = np.reshape(parameters, (-1, num_parameters)).tolist() + batch_size = len(parameters) + + try: + job = self.estimator.run(batch_size * [ansatz], batch_size * [operator], parameters) + estimator_result = job.result() + except Exception as exc: + raise AlgorithmError("The primitive job to evaluate the energy failed!") from exc + + values = estimator_result.values + + if self.callback is not None: + metadata = estimator_result.metadata + for params, value, meta in zip(parameters, values, metadata): + eval_count += 1 + self.callback(eval_count, params, value, meta) + + energy = values[0] if len(values) == 1 else values + + return energy + + return evaluate_energy + + def _get_evaluate_gradient( + self, + ansatz: QuantumCircuit, + operator: BaseOperator | PauliSumOp, + ) -> Callable[[np.ndarray], np.ndarray]: + """Get a function handle to evaluate the gradient at given parameters for the ansatz. + + Args: + ansatz: The ansatz preparing the quantum state. + operator: The operator whose energy to evaluate. + + Returns: + A function handle to evaluate the gradient at given parameters for the ansatz. + + Raises: + AlgorithmError: If the primitive job to evaluate the gradient fails. + """ + + def evaluate_gradient(parameters: np.ndarray) -> np.ndarray: + # broadcasting not required for the estimator gradients + try: + job = self.gradient.run([ansatz], [operator], [parameters]) + gradients = job.result().gradients + except Exception as exc: + raise AlgorithmError("The primitive job to evaluate the gradient failed!") from exc + + return gradients[0] + + return evaluate_gradient + + def _check_operator_ansatz(self, operator: BaseOperator | PauliSumOp): + """Check that the number of qubits of operator and ansatz match and that the ansatz is + parameterized. + """ + if operator.num_qubits != self.ansatz.num_qubits: + try: + logger.info( + "Trying to resize ansatz to match operator on %s qubits.", operator.num_qubits + ) + self.ansatz.num_qubits = operator.num_qubits + except AttributeError as error: + raise AlgorithmError( + "The number of qubits of the ansatz does not match the " + "operator, and the ansatz does not allow setting the " + "number of qubits using `num_qubits`." + ) from error + + if self.ansatz.num_parameters == 0: + raise AlgorithmError("The ansatz must be parameterized, but has no free parameters.") + + def _build_vqe_result( + self, + optimizer_result: OptimizerResult, + aux_operators_evaluated: ListOrDict[tuple[complex, tuple[complex, int]]], + optimizer_time: float, + ) -> VQEResult: + result = VQEResult() + result.eigenvalue = optimizer_result.fun + result.cost_function_evals = optimizer_result.nfev + result.optimal_point = optimizer_result.x + result.optimal_parameters = dict(zip(self.ansatz.parameters, optimizer_result.x)) + result.optimal_value = optimizer_result.fun + result.optimizer_time = optimizer_time + result.aux_operators_evaluated = aux_operators_evaluated + result.optimizer_result = optimizer_result + return result + + +class VQEResult(VariationalResult, MinimumEigensolverResult): + """Variational quantum eigensolver result.""" + + def __init__(self) -> None: + super().__init__() + self._cost_function_evals = None + + @property + def cost_function_evals(self) -> int | None: + """The number of cost optimizer evaluations.""" + return self._cost_function_evals + + @cost_function_evals.setter + def cost_function_evals(self, value: int) -> None: + self._cost_function_evals = value diff --git a/qiskit/algorithms/observables_evaluator.py b/qiskit/algorithms/observables_evaluator.py index f6a6c3d809..f1ca00c4b5 100644 --- a/qiskit/algorithms/observables_evaluator.py +++ b/qiskit/algorithms/observables_evaluator.py @@ -9,8 +9,12 @@ # Any modifications or derivative works of this code must retain this # copyright notice, and modified files need to carry a notice indicating # that they have been altered from the originals. + """Evaluator of observables for algorithms.""" + from __future__ import annotations +from collections.abc import Sequence +from typing import Any import numpy as np @@ -18,7 +22,7 @@ from qiskit import QuantumCircuit from qiskit.opflow import PauliSumOp from .exceptions import AlgorithmError from .list_or_dict import ListOrDict -from ..primitives import EstimatorResult, BaseEstimator +from ..primitives import BaseEstimator from ..quantum_info.operators.base_operator import BaseOperator @@ -26,38 +30,31 @@ def estimate_observables( estimator: BaseEstimator, quantum_state: QuantumCircuit, observables: ListOrDict[BaseOperator | PauliSumOp], + parameter_values: Sequence[float] | None = None, threshold: float = 1e-12, -) -> ListOrDict[tuple[complex, tuple[complex, int]]]: +) -> ListOrDict[tuple[complex, dict[str, Any]]]: """ Accepts a sequence of operators and calculates their expectation values - means - and standard deviations. They are calculated with respect to a quantum state provided. A user + and metadata. They are calculated with respect to a quantum state provided. A user can optionally provide a threshold value which filters mean values falling below the threshold. Args: estimator: An estimator primitive used for calculations. - quantum_state: An unparametrized quantum circuit representing a quantum state that - expectation values are computed against. + quantum_state: A (parameterized) quantum circuit preparing a quantum state that expectation + values are computed against. observables: A list or a dictionary of operators whose expectation values are to be calculated. + parameter_values: Optional list of parameters values to evaluate the quantum circuit on. threshold: A threshold value that defines which mean values should be neglected (helpful for ignoring numerical instabilities close to 0). Returns: - A list or a dictionary of tuples (mean, (variance, shots)). + A list or a dictionary of tuples (mean, metadata). Raises: - ValueError: If a ``quantum_state`` with free parameters is provided. AlgorithmError: If a primitive job is not successful. """ - if ( - isinstance(quantum_state, QuantumCircuit) # State cannot be parametrized - and len(quantum_state.parameters) > 0 - ): - raise ValueError( - "A parametrized representation of a quantum_state was provided. It is not " - "allowed - it cannot have free parameters." - ) if isinstance(observables, dict): observables_list = list(observables.values()) else: @@ -65,18 +62,19 @@ def estimate_observables( observables_list = _handle_zero_ops(observables_list) quantum_state = [quantum_state] * len(observables) + if parameter_values is not None: + parameter_values = [parameter_values] * len(observables) try: - estimator_job = estimator.run(quantum_state, observables_list) + estimator_job = estimator.run(quantum_state, observables_list, parameter_values) expectation_values = estimator_job.result().values except Exception as exc: raise AlgorithmError("The primitive job failed!") from exc - variance_and_shots = _prep_variance_and_shots(estimator_job, len(expectation_values)) - + metadata = estimator_job.result().metadata # Discard values below threshold observables_means = expectation_values * (np.abs(expectation_values) > threshold) - # zip means and standard deviations into tuples - observables_results = list(zip(observables_means, variance_and_shots)) + # zip means and metadata into tuples + observables_results = list(zip(observables_means, metadata)) return _prepare_result(observables_results, observables) @@ -95,20 +93,20 @@ def _handle_zero_ops( def _prepare_result( - observables_results: list[tuple[complex, tuple[complex, int]]], + observables_results: list[tuple[complex, dict]], observables: ListOrDict[BaseOperator | PauliSumOp], -) -> ListOrDict[tuple[complex, tuple[complex, int]]]: +) -> ListOrDict[tuple[complex, dict[str, Any]]]: """ - Prepares a list of tuples of eigenvalues and (variance, shots) tuples from + Prepares a list of tuples of eigenvalues and metadata tuples from ``observables_results`` and ``observables``. Args: - observables_results: A list of tuples (mean, (variance, shots)). + observables_results: A list of tuples (mean, metadata). observables: A list or a dictionary of operators whose expectation values are to be calculated. Returns: - A list or a dictionary of tuples (mean, (variance, shots)). + A list or a dictionary of tuples (mean, metadata). """ if isinstance(observables, list): @@ -122,36 +120,3 @@ def _prepare_result( for key, value in key_value_iterator: observables_eigenvalues[key] = value return observables_eigenvalues - - -def _prep_variance_and_shots( - estimator_result: EstimatorResult, - results_length: int, -) -> list[tuple[complex, int]]: - """ - Prepares a list of tuples with variances and shots from results provided by expectation values - calculations. If there is no variance or shots data available from a primitive, the values will - be set to ``0``. - - Args: - estimator_result: An estimator result. - results_length: Number of expectation values calculated. - - Returns: - A list of tuples of the form (variance, shots). - """ - if not estimator_result.metadata: - return [(0, 0)] * results_length - - results = [] - for metadata in estimator_result.metadata: - variance, shots = 0.0, 0 - if metadata: - if "variance" in metadata.keys(): - variance = metadata["variance"] - if "shots" in metadata.keys(): - shots = metadata["shots"] - - results.append((variance, shots)) - - return results diff --git a/qiskit/algorithms/optimizers/qnspsa.py b/qiskit/algorithms/optimizers/qnspsa.py index f74f91e51f..42296f9ddf 100644 --- a/qiskit/algorithms/optimizers/qnspsa.py +++ b/qiskit/algorithms/optimizers/qnspsa.py @@ -187,7 +187,9 @@ class QNSPSA(SPSA): self._nfev += 6 loss_values = _batch_evaluate(loss, loss_points, self._max_evals_grouped) - fidelity_values = _batch_evaluate(self.fidelity, fidelity_points, self._max_evals_grouped) + fidelity_values = _batch_evaluate( + self.fidelity, fidelity_points, self._max_evals_grouped, unpack_points=True + ) # compute the gradient approximation and additionally return the loss function evaluations gradient_estimate = (loss_values[0] - loss_values[1]) / (2 * eps) * delta1 @@ -268,11 +270,13 @@ class QNSPSA(SPSA): sampler = CircuitSampler(backend) def fidelity(values_x, values_y=None): - if values_y is not None: # no batches + # no batches + if isinstance(values_x, np.ndarray) and isinstance(values_y, np.ndarray): value_dict = dict( zip(params_x[:] + params_y[:], values_x.tolist() + values_y.tolist()) ) - else: + # legacy batching -- remove once QNSPSA.get_fidelity is only supported with sampler + elif values_y is None: value_dict = {p: [] for p in params_x[:] + params_y[:]} for values_xy in values_x: for value_x, param_x in zip(values_xy[0, :], params_x): @@ -280,6 +284,14 @@ class QNSPSA(SPSA): for value_y, param_y in zip(values_xy[1, :], params_y): value_dict[param_y].append(value_y) + else: + value_dict = {p: [] for p in params_x[:] + params_y[:]} + for values_i_x, values_i_y in zip(values_x, values_y): + for value_x, param_x in zip(values_i_x, params_x): + value_dict[param_x].append(value_x) + + for value_y, param_y in zip(values_i_y, params_y): + value_dict[param_y].append(value_y) return np.abs(sampler.convert(expression, params=value_dict).eval()) ** 2 diff --git a/qiskit/algorithms/optimizers/spsa.py b/qiskit/algorithms/optimizers/spsa.py index 10d81dfdd4..f12945762c 100644 --- a/qiskit/algorithms/optimizers/spsa.py +++ b/qiskit/algorithms/optimizers/spsa.py @@ -711,7 +711,7 @@ def constant(eta=0.01): yield eta -def _batch_evaluate(function, points, max_evals_grouped): +def _batch_evaluate(function, points, max_evals_grouped, unpack_points=False): # if the function cannot handle lists of points as input, cover this case immediately if max_evals_grouped == 1: # support functions with multiple arguments where the points are given in a tuple @@ -731,11 +731,32 @@ def _batch_evaluate(function, points, max_evals_grouped): results = [] for batch in batched_points: - results += function(batch).tolist() + if unpack_points: + batch = _repack_points(batch) + results += _as_list(function(*batch)) + else: + results += _as_list(function(batch)) return results +def _as_list(obj): + """Convert a list or numpy array into a list.""" + return obj.tolist() if isinstance(obj, np.ndarray) else obj + + +def _repack_points(points): + """Turn a list of tuples of points into a tuple of lists of points. + E.g. turns + [(a1, a2, a3), (b1, b2, b3)] + into + ([a1, b1], [a2, b2], [a3, b3]) + where all elements are np.ndarray. + """ + num_sets = len(points[0]) # length of (a1, a2, a3) + return ([x[i] for x in points] for i in range(num_sets)) + + def _make_spd(matrix, bias=0.01): identity = np.identity(matrix.shape[0]) psd = scipy.linalg.sqrtm(matrix.dot(matrix)) diff --git a/qiskit/algorithms/time_evolvers/trotterization/trotter_qrte.py b/qiskit/algorithms/time_evolvers/trotterization/trotter_qrte.py index d244f718ce..8ee855d60c 100644 --- a/qiskit/algorithms/time_evolvers/trotterization/trotter_qrte.py +++ b/qiskit/algorithms/time_evolvers/trotterization/trotter_qrte.py @@ -168,6 +168,7 @@ class TrotterQRTE(RealTimeEvolver): self.estimator, evolved_state, evolution_problem.aux_operators, + None, evolution_problem.truncation_threshold, ) diff --git a/qiskit/algorithms/utils/__init__.py b/qiskit/algorithms/utils/__init__.py new file mode 100644 index 0000000000..2b49396270 --- /dev/null +++ b/qiskit/algorithms/utils/__init__.py @@ -0,0 +1,21 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. + +"""Common Qiskit algorithms utility functions.""" + +from .validate_initial_point import validate_initial_point +from .validate_bounds import validate_bounds + +__all__ = [ + "validate_initial_point", + "validate_bounds", +] diff --git a/qiskit/algorithms/utils/validate_bounds.py b/qiskit/algorithms/utils/validate_bounds.py new file mode 100644 index 0000000000..65f00cfcbd --- /dev/null +++ b/qiskit/algorithms/utils/validate_bounds.py @@ -0,0 +1,44 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. + +"""Validate parameter bounds.""" + +from __future__ import annotations + +from qiskit.circuit import QuantumCircuit + + +def validate_bounds(circuit: QuantumCircuit) -> list[tuple(float | None, float | None)]: + """ + Validate the bounds provided by a quantum circuit against its number of parameters. + If no bounds are obtained, return ``None`` for all lower and upper bounds. + + Args: + circuit: A parameterized quantum circuit. + + Returns: + A list of tuples (lower_bound, upper_bound)). + + Raises: + ValueError: If the number of bounds does not the match the number of circuit parameters. + """ + if hasattr(circuit, "parameter_bounds") and circuit.parameter_bounds is not None: + bounds = circuit.parameter_bounds + if len(bounds) != circuit.num_parameters: + raise ValueError( + f"The number of bounds ({len(bounds)}) does not match the number of " + f"parameters in the circuit ({circuit.num_parameters})." + ) + else: + bounds = [(None, None)] * circuit.num_parameters + + return bounds diff --git a/qiskit/algorithms/utils/validate_initial_point.py b/qiskit/algorithms/utils/validate_initial_point.py new file mode 100644 index 0000000000..15e871cdd4 --- /dev/null +++ b/qiskit/algorithms/utils/validate_initial_point.py @@ -0,0 +1,69 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. + +"""Validate an initial point.""" + +from __future__ import annotations + +from collections.abc import Sequence + +import numpy as np + +from qiskit.circuit import QuantumCircuit +from qiskit.utils import algorithm_globals + + +def validate_initial_point( + point: Sequence[float] | None, circuit: QuantumCircuit +) -> Sequence[float]: + r""" + Validate a choice of initial point against a choice of circuit. If no point is provided, a + random point will be generated within certain parameter bounds. It will first look to the + circuit for these bounds. If the circuit does not specify bounds, bounds of :math:`-2\pi`, + :math:`2\pi` will be used. + + Args: + point: An initial point. + circuit: A parameterized quantum circuit. + + Returns: + A validated initial point. + + Raises: + ValueError: If the dimension of the initial point does not match the number of circuit + parameters. + """ + expected_size = circuit.num_parameters + + if point is None: + # get bounds if circuit has them set, otherwise use [-2pi, 2pi] for each parameter + bounds = getattr(circuit, "parameter_bounds", None) + if bounds is None: + bounds = [(-2 * np.pi, 2 * np.pi)] * expected_size + + # replace all Nones by [-2pi, 2pi] + lower_bounds = [] + upper_bounds = [] + for lower, upper in bounds: + lower_bounds.append(lower if lower is not None else -2 * np.pi) + upper_bounds.append(upper if upper is not None else 2 * np.pi) + + # sample from within bounds + point = algorithm_globals.random.uniform(lower_bounds, upper_bounds) + + elif len(point) != expected_size: + raise ValueError( + f"The dimension of the initial point ({len(point)}) does not match the " + f"number of parameters in the circuit ({expected_size})." + ) + + return point diff --git a/qiskit/algorithms/variational_algorithm.py b/qiskit/algorithms/variational_algorithm.py index 92d8976ffc..0be89abf6b 100644 --- a/qiskit/algorithms/variational_algorithm.py +++ b/qiskit/algorithms/variational_algorithm.py @@ -31,6 +31,7 @@ from abc import ABC, abstractmethod import numpy as np from .algorithm_result import AlgorithmResult +from .optimizers import OptimizerResult class VariationalAlgorithm(ABC): @@ -109,3 +110,13 @@ class VariationalResult(AlgorithmResult): def optimal_parameters(self, value: Dict) -> None: """Sets optimal parameters""" self._optimal_parameters = value + + @property + def optimizer_result(self) -> Optional[OptimizerResult]: + """Returns the optimizer result""" + return self._optimizer_result + + @optimizer_result.setter + def optimizer_result(self, value: OptimizerResult) -> None: + """Sets optimizer result""" + self._optimizer_result = value diff --git a/releasenotes/notes/0.19/add-getters-and-setters-for-vqe-edc753591b368980.yaml b/releasenotes/notes/0.19/add-getters-and-setters-for-vqe-edc753591b368980.yaml index a1f20e9220..1d6ca9b826 100644 --- a/releasenotes/notes/0.19/add-getters-and-setters-for-vqe-edc753591b368980.yaml +++ b/releasenotes/notes/0.19/add-getters-and-setters-for-vqe-edc753591b368980.yaml @@ -3,8 +3,10 @@ features: - | Every attribute of the :class:`~qiskit.algorithms.VQE` class that is set at the initialization is now accessible with getters and setters. Further, the - default values of the VQE attributes :attr:`~.VQE.ansatz` and - :attr:`~.VQE.optimizer` can be reset by assigning ``None`` to them:: + default values of the VQE attributes + :attr:`~qiskit.algorithms.minimimum_eigen_solvers.VQE.ansatz` and + :attr:`~qiskit.algorithms.minimimum_eigen_solvers.VQE.optimizer` can be + reset by assigning ``None`` to them:: vqe = VQE(my_ansatz, my_optimizer) vqe.ansatz = None # reset to default: RealAmplitudes ansatz diff --git a/releasenotes/notes/0.19/support-dict-for-aux-operators-c3c9ad380c208afd.yaml b/releasenotes/notes/0.19/support-dict-for-aux-operators-c3c9ad380c208afd.yaml index 5fef8f5d64..63af5ca0dc 100644 --- a/releasenotes/notes/0.19/support-dict-for-aux-operators-c3c9ad380c208afd.yaml +++ b/releasenotes/notes/0.19/support-dict-for-aux-operators-c3c9ad380c208afd.yaml @@ -1,9 +1,11 @@ --- features: - | - The :obj:`.Eigensolver` and :obj:`.MinimumEigensolver` interfaces now support the type + The :obj:`.Eigensolver` and :obj:`~qiskit.algorithms.minimimum_eigen_solvers.MinimumEigensolver` + interfaces now support the type ``Dict[str, Optional[OperatorBase]]`` for the ``aux_operators`` parameter in their respective - :meth:`~.Eigensolver.compute_eigenvalues` and :meth:`~.MinimumEigensolver.compute_minimum_eigenvalue` methods. + :meth:`~.Eigensolver.compute_eigenvalues` and + :meth:`~qiskit.algorithms.minimimum_eigen_solvers.MinimumEigensolver.compute_minimum_eigenvalue` methods. In this case, the auxiliary eigenvalues are also stored in a dictionary under the same keys provided by the ``aux_operators`` dictionary. Keys that correspond to an operator that does not commute with the main operator are dropped. diff --git a/releasenotes/notes/0.19/taper_empty_operator_fix-53ce20e5d2b68fd6.yaml b/releasenotes/notes/0.19/taper_empty_operator_fix-53ce20e5d2b68fd6.yaml index b02a9dde71..a719f78310 100644 --- a/releasenotes/notes/0.19/taper_empty_operator_fix-53ce20e5d2b68fd6.yaml +++ b/releasenotes/notes/0.19/taper_empty_operator_fix-53ce20e5d2b68fd6.yaml @@ -1,11 +1,11 @@ --- fixes: - | - When tapering an empty zero operator in :mod:`qiskit.opflow`, the code, on detecting it was zero, logged a - warning and returned the original operator. Such operators are commonly found in - the auxiliary operators, when using Qiskit Nature, and the above behavior caused :obj:`.VQE` - to throw an exception as tapered non-zero operators were a different number of qubits - from the tapered zero operators (since taper has returned the input operator unchanged). - The code will now correctly taper a zero operator such that the number of qubits is - reduced as expected and matches to tapered non-zero operators e.g ```0*"IIII"``` when we are - tapering by 3 qubits will become ``0*"I"``. + When tapering an empty zero operator in :mod:`qiskit.opflow`, the code, on detecting it was + zero, logged a warning and returned the original operator. Such operators are commonly found in + the auxiliary operators, when using Qiskit Nature, and the above behavior caused + :obj:`~qiskit.algorithms.minimimum_eigen_solvers.VQE` to throw an exception as tapered non-zero + operators were a different number of qubits from the tapered zero operators (since taper has + returned the input operator unchanged). The code will now correctly taper a zero operator such + that the number of qubits is reduced as expected and matches to tapered non-zero operators e.g + ```0*"IIII"``` when we are tapering by 3 qubits will become ``0*"I"``. diff --git a/releasenotes/notes/vqe-with-estimator-primitive-7cbcc462ad4dc593.yaml b/releasenotes/notes/vqe-with-estimator-primitive-7cbcc462ad4dc593.yaml new file mode 100644 index 0000000000..b500ac2055 --- /dev/null +++ b/releasenotes/notes/vqe-with-estimator-primitive-7cbcc462ad4dc593.yaml @@ -0,0 +1,34 @@ +--- +features: + - | + Added the :class:`qiskit.algorithms.minimum_eigensolvers` package to include interfaces for + primitive-enabled algorithms. :class:`~qiskit.algorithms.minimum_eigensolvers.VQE` has been + refactored in this implementation to leverage primitives. + + To use the new implementation with a reference primitive, one can do, for example: + + .. code-block:: python + + from qiskit.algorithms.minimum_eigensolvers import VQE + from qiskit.algorithms.optimizers import SLSQP + from qiskit.circuit.library import TwoLocal + from qiskit.primitives import Estimator + + h2_op = SparsePauliOp( + ["II", "IZ", "ZI", "ZZ", "XX"], + coeffs=[ + -1.052373245772859, + 0.39793742484318045, + -0.39793742484318045, + -0.01128010425623538, + 0.18093119978423156, + ], + ) + + estimator = Estimator() + ansatz = TwoLocal(rotation_blocks=["ry", "rz"], entanglement_blocks="cz") + optimizer = SLSQP() + + vqe = VQE(estimator, ansatz, optimizer) + result = vqe.compute_minimum_eigenvalue(h2_op) + eigenvalue = result.eigenvalue diff --git a/test/python/algorithms/minimum_eigensolvers/__init__.py b/test/python/algorithms/minimum_eigensolvers/__init__.py new file mode 100644 index 0000000000..fdb172d367 --- /dev/null +++ b/test/python/algorithms/minimum_eigensolvers/__init__.py @@ -0,0 +1,11 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. diff --git a/test/python/algorithms/minimum_eigensolvers/test_numpy_minimum_eigensolver.py b/test/python/algorithms/minimum_eigensolvers/test_numpy_minimum_eigensolver.py new file mode 100644 index 0000000000..8d2212e1d7 --- /dev/null +++ b/test/python/algorithms/minimum_eigensolvers/test_numpy_minimum_eigensolver.py @@ -0,0 +1,241 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. + +""" Test NumPy minimum eigensolver """ + +import unittest +from test.python.algorithms import QiskitAlgorithmsTestCase + +import numpy as np +from ddt import ddt, data + +from qiskit.algorithms.minimum_eigensolvers import NumPyMinimumEigensolver +from qiskit.opflow import PauliSumOp + + +@ddt +class TestNumPyMinimumEigensolver(QiskitAlgorithmsTestCase): + """Test NumPy minimum eigensolver""" + + def setUp(self): + super().setUp() + self.qubit_op = PauliSumOp.from_list( + [ + ("II", -1.052373245772859), + ("ZI", 0.39793742484318045), + ("IZ", -0.39793742484318045), + ("ZZ", -0.01128010425623538), + ("XX", 0.18093119978423156), + ] + ) + + aux_op1 = PauliSumOp.from_list([("II", 2.0)]) + aux_op2 = PauliSumOp.from_list([("II", 0.5), ("ZZ", 0.5), ("YY", 0.5), ("XX", -0.5)]) + self.aux_ops_list = [aux_op1, aux_op2] + self.aux_ops_dict = {"aux_op1": aux_op1, "aux_op2": aux_op2} + + def test_cme(self): + """Basic test""" + algo = NumPyMinimumEigensolver() + result = algo.compute_minimum_eigenvalue( + operator=self.qubit_op, aux_operators=self.aux_ops_list + ) + self.assertAlmostEqual(result.eigenvalue, -1.85727503 + 0j) + self.assertEqual(len(result.aux_operators_evaluated), 2) + np.testing.assert_array_almost_equal(result.aux_operators_evaluated[0], [2, 0]) + np.testing.assert_array_almost_equal(result.aux_operators_evaluated[1], [0, 0]) + + def test_cme_reuse(self): + """Test reuse""" + algo = NumPyMinimumEigensolver() + + with self.subTest("Test with no operator or aux_operators, give via compute method"): + result = algo.compute_minimum_eigenvalue(operator=self.qubit_op) + self.assertEqual(result.eigenvalue.dtype, np.float64) + self.assertAlmostEqual(result.eigenvalue, -1.85727503) + self.assertIsNone(result.aux_operators_evaluated) + + with self.subTest("Test with added aux_operators"): + result = algo.compute_minimum_eigenvalue( + operator=self.qubit_op, aux_operators=self.aux_ops_list + ) + self.assertAlmostEqual(result.eigenvalue, -1.85727503 + 0j) + self.assertEqual(len(result.aux_operators_evaluated), 2) + np.testing.assert_array_almost_equal(result.aux_operators_evaluated[0], [2, 0]) + np.testing.assert_array_almost_equal(result.aux_operators_evaluated[1], [0, 0]) + + with self.subTest("Test with aux_operators removed"): + result = algo.compute_minimum_eigenvalue(operator=self.qubit_op, aux_operators=[]) + self.assertEqual(result.eigenvalue.dtype, np.float64) + self.assertAlmostEqual(result.eigenvalue, -1.85727503) + self.assertIsNone(result.aux_operators_evaluated) + + with self.subTest("Test with aux_operators set again"): + result = algo.compute_minimum_eigenvalue( + operator=self.qubit_op, aux_operators=self.aux_ops_list + ) + self.assertAlmostEqual(result.eigenvalue, -1.85727503 + 0j) + self.assertEqual(len(result.aux_operators_evaluated), 2) + np.testing.assert_array_almost_equal(result.aux_operators_evaluated[0], [2, 0]) + np.testing.assert_array_almost_equal(result.aux_operators_evaluated[1], [0, 0]) + + with self.subTest("Test after setting first aux_operators as main operator"): + result = algo.compute_minimum_eigenvalue( + operator=self.aux_ops_list[0], aux_operators=[] + ) + self.assertAlmostEqual(result.eigenvalue, 2 + 0j) + self.assertIsNone(result.aux_operators_evaluated) + + def test_cme_filter(self): + """Basic test""" + + # define filter criterion + # pylint: disable=unused-argument + def criterion(x, v, a_v): + return v >= -0.5 + + algo = NumPyMinimumEigensolver(filter_criterion=criterion) + result = algo.compute_minimum_eigenvalue( + operator=self.qubit_op, aux_operators=self.aux_ops_list + ) + self.assertAlmostEqual(result.eigenvalue, -0.22491125 + 0j) + self.assertEqual(len(result.aux_operators_evaluated), 2) + np.testing.assert_array_almost_equal(result.aux_operators_evaluated[0], [2, 0]) + np.testing.assert_array_almost_equal(result.aux_operators_evaluated[1], [0, 0]) + + def test_cme_filter_empty(self): + """Test with filter always returning False""" + + # define filter criterion + # pylint: disable=unused-argument + def criterion(x, v, a_v): + return False + + algo = NumPyMinimumEigensolver(filter_criterion=criterion) + result = algo.compute_minimum_eigenvalue( + operator=self.qubit_op, aux_operators=self.aux_ops_list + ) + self.assertEqual(result.eigenvalue, None) + self.assertEqual(result.eigenstate, None) + self.assertEqual(result.aux_operators_evaluated, None) + + @data("X", "Y", "Z") + def test_cme_1q(self, op): + """Test for 1 qubit operator""" + algo = NumPyMinimumEigensolver() + operator = PauliSumOp.from_list([(op, 1.0)]) + result = algo.compute_minimum_eigenvalue(operator=operator) + self.assertAlmostEqual(result.eigenvalue, -1) + + def test_cme_aux_ops_dict(self): + """Test dictionary compatibility of aux_operators""" + # Start with an empty dictionary + algo = NumPyMinimumEigensolver() + + with self.subTest("Test with an empty dictionary."): + result = algo.compute_minimum_eigenvalue(operator=self.qubit_op, aux_operators={}) + self.assertAlmostEqual(result.eigenvalue, -1.85727503 + 0j) + self.assertIsNone(result.aux_operators_evaluated) + + with self.subTest("Test with two auxiliary operators."): + result = algo.compute_minimum_eigenvalue( + operator=self.qubit_op, aux_operators=self.aux_ops_dict + ) + self.assertAlmostEqual(result.eigenvalue, -1.85727503 + 0j) + self.assertEqual(len(result.aux_operators_evaluated), 2) + np.testing.assert_array_almost_equal(result.aux_operators_evaluated["aux_op1"], [2, 0]) + np.testing.assert_array_almost_equal(result.aux_operators_evaluated["aux_op2"], [0, 0]) + + with self.subTest("Test with additional zero and None operators."): + extra_ops = {"None_op": None, "zero_op": 0, **self.aux_ops_dict} + result = algo.compute_minimum_eigenvalue( + operator=self.qubit_op, aux_operators=extra_ops + ) + self.assertAlmostEqual(result.eigenvalue, -1.85727503 + 0j) + self.assertEqual(len(result.aux_operators_evaluated), 3) + np.testing.assert_array_almost_equal(result.aux_operators_evaluated["aux_op1"], [2, 0]) + np.testing.assert_array_almost_equal(result.aux_operators_evaluated["aux_op2"], [0, 0]) + self.assertEqual(result.aux_operators_evaluated["zero_op"], (0.0, 0)) + + def test_aux_operators_list(self): + """Test list-based aux_operators.""" + aux_op1 = PauliSumOp.from_list([("II", 2.0)]) + aux_op2 = PauliSumOp.from_list([("II", 0.5), ("ZZ", 0.5), ("YY", 0.5), ("XX", -0.5)]) + aux_ops = [aux_op1, aux_op2] + algo = NumPyMinimumEigensolver() + + with self.subTest("Test with two auxiliary operators."): + result = algo.compute_minimum_eigenvalue(operator=self.qubit_op, aux_operators=aux_ops) + self.assertAlmostEqual(result.eigenvalue, -1.85727503 + 0j) + self.assertEqual(len(result.aux_operators_evaluated), 2) + # expectation values + self.assertAlmostEqual(result.aux_operators_evaluated[0][0], 2, places=6) + self.assertAlmostEqual(result.aux_operators_evaluated[1][0], 0, places=6) + # standard deviations + self.assertAlmostEqual(result.aux_operators_evaluated[0][1], 0.0) + self.assertAlmostEqual(result.aux_operators_evaluated[1][1], 0.0) + + with self.subTest("Test with additional zero and None operators."): + extra_ops = [*aux_ops, None, 0] + result = algo.compute_minimum_eigenvalue( + operator=self.qubit_op, aux_operators=extra_ops + ) + self.assertAlmostEqual(result.eigenvalue, -1.85727503 + 0j) + self.assertEqual(len(result.aux_operators_evaluated), 4) + # expectation values + self.assertAlmostEqual(result.aux_operators_evaluated[0][0], 2, places=6) + self.assertAlmostEqual(result.aux_operators_evaluated[1][0], 0, places=6) + self.assertIsNone(result.aux_operators_evaluated[2], None) + self.assertEqual(result.aux_operators_evaluated[3][0], 0.0) + # standard deviations + self.assertAlmostEqual(result.aux_operators_evaluated[0][1], 0.0) + self.assertAlmostEqual(result.aux_operators_evaluated[1][1], 0.0) + self.assertEqual(result.aux_operators_evaluated[3][1], 0.0) + + def test_aux_operators_dict(self): + """Test dict-based aux_operators.""" + aux_op1 = PauliSumOp.from_list([("II", 2.0)]) + aux_op2 = PauliSumOp.from_list([("II", 0.5), ("ZZ", 0.5), ("YY", 0.5), ("XX", -0.5)]) + aux_ops = {"aux_op1": aux_op1, "aux_op2": aux_op2} + algo = NumPyMinimumEigensolver() + + with self.subTest("Test with two auxiliary operators."): + result = algo.compute_minimum_eigenvalue(operator=self.qubit_op, aux_operators=aux_ops) + self.assertAlmostEqual(result.eigenvalue, -1.85727503 + 0j) + self.assertEqual(len(result.aux_operators_evaluated), 2) + # expectation values + self.assertAlmostEqual(result.aux_operators_evaluated["aux_op1"][0], 2, places=6) + self.assertAlmostEqual(result.aux_operators_evaluated["aux_op2"][0], 0, places=6) + # standard deviations + self.assertAlmostEqual(result.aux_operators_evaluated["aux_op1"][1], 0.0) + self.assertAlmostEqual(result.aux_operators_evaluated["aux_op2"][1], 0.0) + + with self.subTest("Test with additional zero and None operators."): + extra_ops = {**aux_ops, "None_operator": None, "zero_operator": 0} + result = algo.compute_minimum_eigenvalue( + operator=self.qubit_op, aux_operators=extra_ops + ) + self.assertAlmostEqual(result.eigenvalue, -1.85727503 + 0j) + self.assertEqual(len(result.aux_operators_evaluated), 3) + # expectation values + self.assertAlmostEqual(result.aux_operators_evaluated["aux_op1"][0], 2, places=6) + self.assertAlmostEqual(result.aux_operators_evaluated["aux_op2"][0], 0, places=6) + self.assertEqual(result.aux_operators_evaluated["zero_operator"][0], 0.0) + self.assertTrue("None_operator" not in result.aux_operators_evaluated.keys()) + # standard deviations + self.assertAlmostEqual(result.aux_operators_evaluated["aux_op1"][1], 0.0) + self.assertAlmostEqual(result.aux_operators_evaluated["aux_op2"][1], 0.0) + self.assertAlmostEqual(result.aux_operators_evaluated["zero_operator"][1], 0.0) + + +if __name__ == "__main__": + unittest.main() diff --git a/test/python/algorithms/minimum_eigensolvers/test_vqe.py b/test/python/algorithms/minimum_eigensolvers/test_vqe.py new file mode 100644 index 0000000000..1ecc870237 --- /dev/null +++ b/test/python/algorithms/minimum_eigensolvers/test_vqe.py @@ -0,0 +1,444 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. + +"""Test the variational quantum eigensolver algorithm.""" + +import unittest +from test.python.algorithms import QiskitAlgorithmsTestCase + +from functools import partial +import numpy as np +from scipy.optimize import minimize as scipy_minimize +from ddt import data, ddt + +from qiskit import QuantumCircuit +from qiskit.algorithms import AlgorithmError +from qiskit.algorithms.gradients import ParamShiftEstimatorGradient +from qiskit.algorithms.minimum_eigensolvers import VQE +from qiskit.algorithms.optimizers import ( + CG, + COBYLA, + GradientDescent, + L_BFGS_B, + OptimizerResult, + P_BFGS, + QNSPSA, + SLSQP, + SPSA, + TNC, +) +from qiskit.algorithms.state_fidelities import ComputeUncompute +from qiskit.circuit.library import RealAmplitudes, TwoLocal +from qiskit.opflow import PauliSumOp, TwoQubitReduction +from qiskit.quantum_info import SparsePauliOp, Operator, Pauli +from qiskit.primitives import Estimator, Sampler +from qiskit.utils import algorithm_globals + +# pylint: disable=invalid-name +def _mock_optimizer(fun, x0, jac=None, bounds=None, inputs=None) -> OptimizerResult: + """A mock of a callable that can be used as minimizer in the VQE.""" + result = OptimizerResult() + result.x = np.zeros_like(x0) + result.fun = fun(result.x) + result.nit = 0 + + if inputs is not None: + inputs.update({"fun": fun, "x0": x0, "jac": jac, "bounds": bounds}) + return result + + +@ddt +class TestVQE(QiskitAlgorithmsTestCase): + """Test VQE""" + + def setUp(self): + super().setUp() + self.seed = 50 + algorithm_globals.random_seed = self.seed + self.h2_op = SparsePauliOp( + ["II", "IZ", "ZI", "ZZ", "XX"], + coeffs=[ + -1.052373245772859, + 0.39793742484318045, + -0.39793742484318045, + -0.01128010425623538, + 0.18093119978423156, + ], + ) + self.h2_energy = -1.85727503 + + self.ryrz_wavefunction = TwoLocal(rotation_blocks=["ry", "rz"], entanglement_blocks="cz") + self.ry_wavefunction = TwoLocal(rotation_blocks="ry", entanglement_blocks="cz") + + @data(L_BFGS_B(), COBYLA()) + def test_basic_aer_statevector(self, estimator): + """Test VQE using reference Estimator.""" + vqe = VQE(Estimator(), self.ryrz_wavefunction, estimator) + + result = vqe.compute_minimum_eigenvalue(operator=self.h2_op) + + with self.subTest(msg="test eigenvalue"): + self.assertAlmostEqual(result.eigenvalue.real, self.h2_energy, places=5) + + with self.subTest(msg="test optimal_value"): + self.assertAlmostEqual(result.optimal_value, self.h2_energy) + + with self.subTest(msg="test dimension of optimal point"): + self.assertEqual(len(result.optimal_point), 16) + + with self.subTest(msg="assert cost_function_evals is set"): + self.assertIsNotNone(result.cost_function_evals) + + with self.subTest(msg="assert optimizer_time is set"): + self.assertIsNotNone(result.optimizer_time) + + with self.subTest(msg="assert optimizer_result is set"): + self.assertIsNotNone(result.optimizer_result) + + with self.subTest(msg="assert optimizer_result."): + self.assertAlmostEqual(result.optimizer_result.fun, self.h2_energy, places=5) + + def test_invalid_initial_point(self): + """Test the proper error is raised when the initial point has the wrong size.""" + ansatz = self.ryrz_wavefunction + initial_point = np.array([1]) + + vqe = VQE( + Estimator(), + ansatz, + SLSQP(), + initial_point=initial_point, + ) + + with self.assertRaises(ValueError): + _ = vqe.compute_minimum_eigenvalue(operator=self.h2_op) + + def test_ansatz_resize(self): + """Test the ansatz is properly resized if it's a blueprint circuit.""" + ansatz = RealAmplitudes(1, reps=1) + vqe = VQE(Estimator(), ansatz, SLSQP()) + result = vqe.compute_minimum_eigenvalue(self.h2_op) + self.assertAlmostEqual(result.eigenvalue.real, self.h2_energy, places=5) + + def test_invalid_ansatz_size(self): + """Test an error is raised if the ansatz has the wrong number of qubits.""" + ansatz = QuantumCircuit(1) + ansatz.compose(RealAmplitudes(1, reps=2)) + vqe = VQE(Estimator(), ansatz, SLSQP()) + + with self.assertRaises(AlgorithmError): + _ = vqe.compute_minimum_eigenvalue(operator=self.h2_op) + + def test_missing_ansatz_params(self): + """Test specifying an ansatz with no parameters raises an error.""" + ansatz = QuantumCircuit(self.h2_op.num_qubits) + vqe = VQE(Estimator(), ansatz, SLSQP()) + with self.assertRaises(AlgorithmError): + vqe.compute_minimum_eigenvalue(operator=self.h2_op) + + def test_max_evals_grouped(self): + """Test with SLSQP with max_evals_grouped.""" + optimizer = SLSQP(maxiter=50, max_evals_grouped=5) + vqe = VQE( + Estimator(), + self.ryrz_wavefunction, + optimizer, + ) + result = vqe.compute_minimum_eigenvalue(operator=self.h2_op) + self.assertAlmostEqual(result.eigenvalue.real, self.h2_energy, places=5) + + @data( + CG(), + L_BFGS_B(), + P_BFGS(), + SLSQP(), + TNC(), + ) + def test_with_gradient(self, optimizer): + """Test VQE using gradient primitive.""" + estimator = Estimator() + vqe = VQE( + estimator, + self.ry_wavefunction, + optimizer, + gradient=ParamShiftEstimatorGradient(estimator), + ) + result = vqe.compute_minimum_eigenvalue(operator=self.h2_op) + self.assertAlmostEqual(result.eigenvalue.real, self.h2_energy, places=5) + + def test_gradient_passed(self): + """Test the gradient is properly passed into the optimizer.""" + inputs = {} + estimator = Estimator() + vqe = VQE( + estimator, + RealAmplitudes(), + partial(_mock_optimizer, inputs=inputs), + gradient=ParamShiftEstimatorGradient(estimator), + ) + _ = vqe.compute_minimum_eigenvalue(operator=self.h2_op) + + self.assertIsNotNone(inputs["jac"]) + + def test_gradient_run(self): + """Test using the gradient to calculate the minimum.""" + estimator = Estimator() + vqe = VQE( + estimator, + RealAmplitudes(), + GradientDescent(maxiter=200, learning_rate=0.1), + gradient=ParamShiftEstimatorGradient(estimator), + ) + result = vqe.compute_minimum_eigenvalue(operator=self.h2_op) + self.assertAlmostEqual(result.eigenvalue.real, self.h2_energy, places=5) + + def test_with_two_qubit_reduction(self): + """Test the VQE using TwoQubitReduction.""" + qubit_op = PauliSumOp.from_list( + [ + ("IIII", -0.8105479805373266), + ("IIIZ", 0.17218393261915552), + ("IIZZ", -0.22575349222402472), + ("IZZI", 0.1721839326191556), + ("ZZII", -0.22575349222402466), + ("IIZI", 0.1209126326177663), + ("IZZZ", 0.16892753870087912), + ("IXZX", -0.045232799946057854), + ("ZXIX", 0.045232799946057854), + ("IXIX", 0.045232799946057854), + ("ZXZX", -0.045232799946057854), + ("ZZIZ", 0.16614543256382414), + ("IZIZ", 0.16614543256382414), + ("ZZZZ", 0.17464343068300453), + ("ZIZI", 0.1209126326177663), + ] + ) + tapered_qubit_op = TwoQubitReduction(num_particles=2).convert(qubit_op) + vqe = VQE( + Estimator(), + self.ry_wavefunction, + SPSA(maxiter=300, last_avg=5), + ) + result = vqe.compute_minimum_eigenvalue(tapered_qubit_op) + self.assertAlmostEqual(result.eigenvalue.real, self.h2_energy, places=2) + + def test_callback(self): + """Test the callback on VQE.""" + history = {"eval_count": [], "parameters": [], "mean": [], "metadata": []} + # expected_shots = options["shots"] if "shots" in options else 0 + + def store_intermediate_result(eval_count, parameters, mean, metadata): + history["eval_count"].append(eval_count) + history["parameters"].append(parameters) + history["mean"].append(mean) + history["metadata"].append(metadata) + + optimizer = COBYLA(maxiter=3) + wavefunction = self.ry_wavefunction + + estimator = Estimator() + + vqe = VQE( + estimator, + wavefunction, + optimizer, + callback=store_intermediate_result, + ) + vqe.compute_minimum_eigenvalue(operator=self.h2_op) + + self.assertTrue(all(isinstance(count, int) for count in history["eval_count"])) + self.assertTrue(all(isinstance(mean, float) for mean in history["mean"])) + self.assertTrue(all(isinstance(metadata, dict) for metadata in history["metadata"])) + for params in history["parameters"]: + self.assertTrue(all(isinstance(param, float) for param in params)) + + def test_reuse(self): + """Test re-using a VQE algorithm instance.""" + ansatz = TwoLocal(rotation_blocks=["ry", "rz"], entanglement_blocks="cz") + vqe = VQE(Estimator(), ansatz, SLSQP(maxiter=300)) + with self.subTest(msg="assert VQE works once all info is available"): + result = vqe.compute_minimum_eigenvalue(operator=self.h2_op) + self.assertAlmostEqual(result.eigenvalue.real, self.h2_energy, places=5) + + operator = Operator(np.array([[1, 0, 0, 0], [0, -1, 0, 0], [0, 0, 2, 0], [0, 0, 0, 3]])) + + with self.subTest(msg="assert vqe works on re-use."): + result = vqe.compute_minimum_eigenvalue(operator=operator) + self.assertAlmostEqual(result.eigenvalue.real, -1.0, places=5) + + def test_vqe_optimizer_reuse(self): + """Test running same VQE twice to re-use optimizer, then switch optimizer""" + vqe = VQE( + Estimator(), + self.ryrz_wavefunction, + SLSQP(), + ) + + def run_check(): + result = vqe.compute_minimum_eigenvalue(operator=self.h2_op) + self.assertAlmostEqual(result.eigenvalue.real, self.h2_energy, places=5) + + run_check() + + with self.subTest("Optimizer re-use."): + run_check() + + with self.subTest("Optimizer replace."): + vqe.optimizer = L_BFGS_B() + run_check() + + def test_batch_evaluate_with_qnspsa(self): + """Test batch evaluating with QNSPSA works.""" + ansatz = TwoLocal(2, rotation_blocks=["ry", "rz"], entanglement_blocks="cz") + + wrapped_sampler = Sampler() + inner_sampler = Sampler() + + wrapped_estimator = Estimator() + inner_estimator = Estimator() + + callcount = {"sampler": 0, "estimator": 0} + + def wrapped_estimator_run(*args, **kwargs): + kwargs["callcount"]["estimator"] += 1 + return inner_estimator.run(*args, **kwargs) + + def wrapped_sampler_run(*args, **kwargs): + kwargs["callcount"]["sampler"] += 1 + return inner_sampler.run(*args, **kwargs) + + wrapped_estimator.run = partial(wrapped_estimator_run, callcount=callcount) + wrapped_sampler.run = partial(wrapped_sampler_run, callcount=callcount) + + fidelity = ComputeUncompute(wrapped_sampler) + + def fidelity_callable(left, right): + batchsize = np.asarray(left).shape[0] + job = fidelity.run(batchsize * [ansatz], batchsize * [ansatz], left, right) + return job.result().fidelities + + qnspsa = QNSPSA(fidelity_callable, maxiter=5) + qnspsa.set_max_evals_grouped(100) + + vqe = VQE( + wrapped_estimator, + ansatz, + qnspsa, + ) + _ = vqe.compute_minimum_eigenvalue(Pauli("ZZ")) + + # 5 (fidelity) + expected_sampler_runs = 5 + # 1 calibration + 1 stddev estimation + 1 initial blocking + # + 5 (1 loss + 1 blocking) + 1 return loss + expected_estimator_runs = 1 + 1 + 1 + 5 * 2 + 1 + + self.assertEqual(callcount["sampler"], expected_sampler_runs) + self.assertEqual(callcount["estimator"], expected_estimator_runs) + + def test_optimizer_scipy_callable(self): + """Test passing a SciPy optimizer directly as callable.""" + vqe = VQE( + Estimator(), + self.ryrz_wavefunction, + partial(scipy_minimize, method="L-BFGS-B", options={"maxiter": 10}), + ) + result = vqe.compute_minimum_eigenvalue(self.h2_op) + self.assertAlmostEqual(result.eigenvalue.real, self.h2_energy, places=2) + + def test_optimizer_callable(self): + """Test passing a optimizer directly as callable.""" + ansatz = RealAmplitudes(1, reps=1) + vqe = VQE(Estimator(), ansatz, _mock_optimizer) + result = vqe.compute_minimum_eigenvalue(SparsePauliOp("Z")) + self.assertTrue(np.all(result.optimal_point == np.zeros(ansatz.num_parameters))) + + def test_aux_operators_list(self): + """Test list-based aux_operators.""" + vqe = VQE(Estimator(), self.ry_wavefunction, SLSQP(maxiter=300)) + + with self.subTest("Test with an empty list."): + result = vqe.compute_minimum_eigenvalue(self.h2_op, aux_operators=[]) + self.assertAlmostEqual(result.eigenvalue.real, self.h2_energy, places=6) + self.assertIsInstance(result.aux_operators_evaluated, list) + self.assertEqual(len(result.aux_operators_evaluated), 0) + + with self.subTest("Test with two auxiliary operators."): + aux_op1 = PauliSumOp.from_list([("II", 2.0)]) + aux_op2 = PauliSumOp.from_list([("II", 0.5), ("ZZ", 0.5), ("YY", 0.5), ("XX", -0.5)]) + aux_ops = [aux_op1, aux_op2] + result = vqe.compute_minimum_eigenvalue(self.h2_op, aux_operators=aux_ops) + self.assertAlmostEqual(result.eigenvalue.real, self.h2_energy, places=5) + self.assertEqual(len(result.aux_operators_evaluated), 2) + # expectation values + self.assertAlmostEqual(result.aux_operators_evaluated[0][0], 2.0, places=6) + self.assertAlmostEqual(result.aux_operators_evaluated[1][0], 0.0, places=6) + # metadata + self.assertIsInstance(result.aux_operators_evaluated[0][1], dict) + self.assertIsInstance(result.aux_operators_evaluated[1][1], dict) + + with self.subTest("Test with additional zero operator."): + extra_ops = [*aux_ops, 0] + result = vqe.compute_minimum_eigenvalue(self.h2_op, aux_operators=extra_ops) + self.assertAlmostEqual(result.eigenvalue.real, self.h2_energy, places=5) + self.assertEqual(len(result.aux_operators_evaluated), 3) + # expectation values + self.assertAlmostEqual(result.aux_operators_evaluated[0][0], 2.0, places=6) + self.assertAlmostEqual(result.aux_operators_evaluated[1][0], 0.0, places=6) + self.assertAlmostEqual(result.aux_operators_evaluated[2][0], 0.0) + # metadata + self.assertIsInstance(result.aux_operators_evaluated[0][1], dict) + self.assertIsInstance(result.aux_operators_evaluated[1][1], dict) + self.assertIsInstance(result.aux_operators_evaluated[2][1], dict) + + def test_aux_operators_dict(self): + """Test dictionary compatibility of aux_operators""" + vqe = VQE(Estimator(), self.ry_wavefunction, SLSQP(maxiter=300)) + + with self.subTest("Test with an empty dictionary."): + result = vqe.compute_minimum_eigenvalue(self.h2_op, aux_operators={}) + self.assertAlmostEqual(result.eigenvalue.real, self.h2_energy, places=6) + self.assertIsInstance(result.aux_operators_evaluated, dict) + self.assertEqual(len(result.aux_operators_evaluated), 0) + + with self.subTest("Test with two auxiliary operators."): + aux_op1 = PauliSumOp.from_list([("II", 2.0)]) + aux_op2 = PauliSumOp.from_list([("II", 0.5), ("ZZ", 0.5), ("YY", 0.5), ("XX", -0.5)]) + aux_ops = {"aux_op1": aux_op1, "aux_op2": aux_op2} + result = vqe.compute_minimum_eigenvalue(self.h2_op, aux_operators=aux_ops) + self.assertAlmostEqual(result.eigenvalue.real, self.h2_energy, places=6) + self.assertEqual(len(result.aux_operators_evaluated), 2) + + # expectation values + self.assertAlmostEqual(result.aux_operators_evaluated["aux_op1"][0], 2.0, places=5) + self.assertAlmostEqual(result.aux_operators_evaluated["aux_op2"][0], 0.0, places=5) + # metadata + self.assertIsInstance(result.aux_operators_evaluated["aux_op1"][1], dict) + self.assertIsInstance(result.aux_operators_evaluated["aux_op2"][1], dict) + + with self.subTest("Test with additional zero operator."): + extra_ops = {**aux_ops, "zero_operator": 0} + result = vqe.compute_minimum_eigenvalue(self.h2_op, aux_operators=extra_ops) + self.assertAlmostEqual(result.eigenvalue.real, self.h2_energy, places=6) + self.assertEqual(len(result.aux_operators_evaluated), 3) + # expectation values + self.assertAlmostEqual(result.aux_operators_evaluated["aux_op1"][0], 2.0, places=5) + self.assertAlmostEqual(result.aux_operators_evaluated["aux_op2"][0], 0.0, places=5) + self.assertAlmostEqual(result.aux_operators_evaluated["zero_operator"][0], 0.0) + # metadata + self.assertIsInstance(result.aux_operators_evaluated["aux_op1"][1], dict) + self.assertIsInstance(result.aux_operators_evaluated["aux_op2"][1], dict) + self.assertIsInstance(result.aux_operators_evaluated["zero_operator"][1], dict) + + +if __name__ == "__main__": + unittest.main() diff --git a/test/python/algorithms/test_observables_evaluator.py b/test/python/algorithms/test_observables_evaluator.py index bcff77cdd5..9db5f02b31 100644 --- a/test/python/algorithms/test_observables_evaluator.py +++ b/test/python/algorithms/test_observables_evaluator.py @@ -9,7 +9,9 @@ # Any modifications or derivative works of this code must retain this # copyright notice, and modified files need to carry a notice indicating # that they have been altered from the originals. + """Tests evaluator of auxiliary operators for algorithms.""" + from __future__ import annotations import unittest from typing import Tuple @@ -53,7 +55,7 @@ class TestObservablesEvaluator(QiskitAlgorithmsTestCase): # the exact value is a list of (mean, (variance, shots)) where we expect 0 variance and # 0 shots exact = [ - (Statevector(ansatz).expectation_value(observable), (0, 0)) + (Statevector(ansatz).expectation_value(observable), {}) for observable in observables_list ] @@ -70,7 +72,7 @@ class TestObservablesEvaluator(QiskitAlgorithmsTestCase): observables: ListOrDict[BaseOperator | PauliSumOp], estimator: Estimator, ): - result = estimate_observables(estimator, quantum_state, observables, self.threshold) + result = estimate_observables(estimator, quantum_state, observables, None, self.threshold) if isinstance(observables, dict): np.testing.assert_equal(list(result.keys()), list(expected_result.keys())) @@ -142,8 +144,8 @@ class TestObservablesEvaluator(QiskitAlgorithmsTestCase): state = bound_ansatz estimator = Estimator() observables = [SparsePauliOp(["XX", "YY"]), 0] - result = estimate_observables(estimator, state, observables, self.threshold) - expected_result = [(0.015607318055509564, (0, 0)), (0.0, (0, 0))] + result = estimate_observables(estimator, state, observables, None, self.threshold) + expected_result = [(0.015607318055509564, {}), (0.0, {})] means = [element[0] for element in result] expected_means = [element[0] for element in expected_result] np.testing.assert_array_almost_equal(means, expected_means, decimal=0.01) @@ -152,6 +154,32 @@ class TestObservablesEvaluator(QiskitAlgorithmsTestCase): expected_vars_and_shots = [element[1] for element in expected_result] np.testing.assert_array_equal(vars_and_shots, expected_vars_and_shots) + def test_estimate_observables_shots(self): + """Tests that variances and shots are returned properly.""" + ansatz = EfficientSU2(2) + parameters = np.array( + [1.2, 4.2, 1.4, 2.0, 1.2, 4.2, 1.4, 2.0, 1.2, 4.2, 1.4, 2.0, 1.2, 4.2, 1.4, 2.0], + dtype=float, + ) + + bound_ansatz = ansatz.bind_parameters(parameters) + state = bound_ansatz + estimator = Estimator(options={"shots": 2048}) + observables = [PauliSumOp.from_list([("ZZ", 2.0)])] + result = estimate_observables(estimator, state, observables, None, self.threshold) + exact_result = self.get_exact_expectation(bound_ansatz, observables) + expected_result = [(exact_result[0][0], {"variance": 1.0898, "shots": 2048})] + + means = [element[0] for element in result] + expected_means = [element[0] for element in expected_result] + np.testing.assert_array_almost_equal(means, expected_means, decimal=0.01) + + vars_and_shots = [element[1] for element in result] + expected_vars_and_shots = [element[1] for element in expected_result] + for computed, expected in zip(vars_and_shots, expected_vars_and_shots): + self.assertAlmostEqual(computed.pop("variance"), expected.pop("variance"), 2) + self.assertEqual(computed.pop("shots"), expected.pop("shots")) + if __name__ == "__main__": unittest.main() diff --git a/test/python/algorithms/time_evolvers/test_trotter_qrte.py b/test/python/algorithms/time_evolvers/test_trotter_qrte.py index fbdbb01615..16efa06549 100644 --- a/test/python/algorithms/time_evolvers/test_trotter_qrte.py +++ b/test/python/algorithms/time_evolvers/test_trotter_qrte.py @@ -88,7 +88,10 @@ class TestTrotterQRTE(QiskitAlgorithmsTestCase): ) aux_ops_result = evolution_result.aux_ops_evaluated - expected_aux_ops_result = [(0.078073, (0.0, 0.0)), (0.268286, (0.0, 0.0))] + expected_aux_ops_result = [ + (0.078073, {"variance": 0, "shots": 0}), + (0.268286, {"variance": 0, "shots": 0}), + ] means = [element[0] for element in aux_ops_result] expected_means = [element[0] for element in expected_aux_ops_result] @@ -96,7 +99,10 @@ class TestTrotterQRTE(QiskitAlgorithmsTestCase): vars_and_shots = [element[1] for element in aux_ops_result] expected_vars_and_shots = [element[1] for element in expected_aux_ops_result] - np.testing.assert_array_equal(vars_and_shots, expected_vars_and_shots) + + for computed, expected in zip(vars_and_shots, expected_vars_and_shots): + self.assertAlmostEqual(computed.pop("variance", 0), expected["variance"], 2) + self.assertEqual(computed.pop("shots", 0), expected["shots"]) @data( ( diff --git a/test/python/algorithms/utils/__init__.py b/test/python/algorithms/utils/__init__.py new file mode 100644 index 0000000000..fdb172d367 --- /dev/null +++ b/test/python/algorithms/utils/__init__.py @@ -0,0 +1,11 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. diff --git a/test/python/algorithms/utils/test_validate_bounds.py b/test/python/algorithms/utils/test_validate_bounds.py new file mode 100644 index 0000000000..80fe54ee07 --- /dev/null +++ b/test/python/algorithms/utils/test_validate_bounds.py @@ -0,0 +1,53 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. + +"""Test validate bounds.""" + +from test.python.algorithms import QiskitAlgorithmsTestCase + +from unittest.mock import Mock + +import numpy as np + +from qiskit.algorithms.utils import validate_bounds +from qiskit.utils import algorithm_globals + + +class TestValidateBounds(QiskitAlgorithmsTestCase): + """Test the ``validate_bounds`` utility function.""" + + def setUp(self): + super().setUp() + algorithm_globals.random_seed = 0 + self.bounds = [(-np.pi / 2, np.pi / 2)] + self.ansatz = Mock() + + def test_with_no_ansatz_bounds(self): + """Test with no ansatz bounds.""" + self.ansatz.num_parameters = 1 + self.ansatz.parameter_bounds = None + bounds = validate_bounds(self.ansatz) + self.assertEqual(bounds, [(None, None)]) + + def test_with_ansatz_bounds(self): + """Test with ansatz bounds.""" + self.ansatz.num_parameters = 1 + self.ansatz.parameter_bounds = self.bounds + bounds = validate_bounds(self.ansatz) + self.assertEqual(bounds, self.bounds) + + def test_with_mismatched_num_params(self): + """Test with a mismatched number of parameters and bounds""" + self.ansatz.num_parameters = 2 + self.ansatz.parameter_bounds = self.bounds + with self.assertRaises(ValueError): + _ = validate_bounds(self.ansatz) diff --git a/test/python/algorithms/utils/test_validate_initial_point.py b/test/python/algorithms/utils/test_validate_initial_point.py new file mode 100644 index 0000000000..32dd48cf2c --- /dev/null +++ b/test/python/algorithms/utils/test_validate_initial_point.py @@ -0,0 +1,50 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2022. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. + +"""Test validate initial point.""" + +from test.python.algorithms import QiskitAlgorithmsTestCase + +from unittest.mock import Mock + +import numpy as np + +from qiskit.algorithms.utils import validate_initial_point +from qiskit.utils import algorithm_globals + + +class TestValidateInitialPoint(QiskitAlgorithmsTestCase): + """Test the ``validate_initial_point`` utility function.""" + + def setUp(self): + super().setUp() + algorithm_globals.random_seed = 0 + self.ansatz = Mock() + self.ansatz.num_parameters = 1 + + def test_with_no_initial_point_or_bounds(self): + """Test with no user-defined initial point and no ansatz bounds.""" + self.ansatz.parameter_bounds = None + initial_point = validate_initial_point(None, self.ansatz) + np.testing.assert_array_almost_equal(initial_point, [1.721111]) + + def test_with_no_initial_point(self): + """Test with no user-defined initial point with ansatz bounds.""" + self.ansatz.parameter_bounds = [(-np.pi / 2, np.pi / 2)] + initial_point = validate_initial_point(None, self.ansatz) + np.testing.assert_array_almost_equal(initial_point, [0.430278]) + + def test_with_mismatched_params(self): + """Test with mistmatched parameters and bounds..""" + self.ansatz.parameter_bounds = None + with self.assertRaises(ValueError): + _ = validate_initial_point([1.0, 2.0], self.ansatz)