Merge pull request #2 from Cryoris/qpe

Update to current algorithm paradigm + small fixes

qiskit/algorithms/phase_estimators/hamiltonian_phase_estimation.py

@@ -49,7 +49,6 @@ class HamiltonianPhaseEstimation:
          `arXiv:1809.09697 <https://arxiv.org/abs/1809.09697>`_
     """
 
-
     def __init__(self,
                  num_evaluation_qubits: int,
                  quantum_instance: Optional[Union[QuantumInstance, BaseBackend]] = None) -> None:
@@ -58,31 +57,25 @@ def __init__(self,
                 be estimated as a binary string with this many bits.
             quantum_instance: The quantum instance on which the circuit will be run.
         """
-        self._pe_scale = None
-        self._bound = None
-        self._evolution = None
-        self._hamiltonian = None
-        self._id_coefficient = None
         self._phase_estimation = PhaseEstimation(
             num_evaluation_qubits=num_evaluation_qubits,
             quantum_instance=quantum_instance)
 
-    def _set_scale(self) -> None:
-        if self._bound is None:
-            pe_scale = PhaseEstimationScale.from_pauli_sum(self._hamiltonian)
-            self._pe_scale = pe_scale
-        else:
-            self._pe_scale = PhaseEstimationScale(self._bound)
+    def _get_scale(self, hamiltonian, bound=None) -> None:
+        if bound is None:
+            return PhaseEstimationScale.from_pauli_sum(hamiltonian)
+
+        return PhaseEstimationScale(bound)
 
-    def _get_unitary(self) -> QuantumCircuit:
+    def _get_unitary(self, hamiltonian, pe_scale, evolution) -> QuantumCircuit:
         """Evolve the Hamiltonian to obtain a unitary.
 
         Apply the scaling to the Hamiltonian that has been computed from an eigenvalue bound
         and compute the unitary by applying the evolution object.
         """
         # scale so that phase does not wrap.
-        scaled_hamiltonian = -self._pe_scale.scale * self._hamiltonian
-        unitary = self._evolution.convert(scaled_hamiltonian.exp_i())
+        scaled_hamiltonian = -pe_scale.scale * hamiltonian
+        unitary = evolution.convert(scaled_hamiltonian.exp_i())
         if not isinstance(unitary, QuantumCircuit):
             unitary_circuit = unitary.to_circuit()
         else:
@@ -110,16 +103,13 @@ def estimate(self, hamiltonian: OperatorBase,
                 then a bound will be computed.
 
         Returns:
-               HamiltonianPhaseEstimationResult instance containing the result of the estimation
+            HamiltonianPhaseEstimationResult instance containing the result of the estimation
             and diagnostic information.
 
         Raises:
             ValueError: if `bound` is `None` and `hamiltonian` is not a Pauli sum (i.e. a
             `PauliSumOp` or a `SummedOp` whose terms are `PauliOp`s.)
         """
-
-        self._evolution = evolution
-        self._bound = bound
         # The term propto the identity is removed from hamiltonian.
         # This is done for three reasons:
         # 1. Work around an unknown bug that otherwise causes the energies to be wrong in some
@@ -131,17 +121,24 @@ def estimate(self, hamiltonian: OperatorBase,
         if isinstance(hamiltonian, PauliSumOp):
             hamiltonian = hamiltonian.to_pauli_op()
 
+        # remove identitiy terms
         id_coefficient, hamiltonian_no_id = _remove_identity(hamiltonian)
-        self._hamiltonian = hamiltonian_no_id
-        self._id_coefficient = id_coefficient
-        self._set_scale()
-        unitary = self._get_unitary()
+
+        # get the rescaling object
+        pe_scale = self._get_scale(hamiltonian_no_id, bound)
+
+        # get the unitary
+        unitary = self._get_unitary(hamiltonian_no_id, pe_scale, evolution)
+
+        # run phase estimation
         phase_estimation_result = self._phase_estimation.estimate(
             unitary=unitary, state_preparation=state_preparation)
+
         return HamiltonianPhaseEstimationResult(
             phase_estimation_result=phase_estimation_result,
-            id_coefficient=self._id_coefficient,
-            phase_estimation_scale=self._pe_scale)
+            id_coefficient=id_coefficient,
+            phase_estimation_scale=pe_scale)
+
 
 def _remove_identity(pauli_sum):
     """Remove any identity operators from `pauli_sum`. Return

qiskit/algorithms/phase_estimators/hamiltonian_phase_estimation_result.py

@@ -14,13 +14,12 @@
 
 
 from typing import Dict, Union, cast
-import numpy
-from qiskit.result import Result
+from qiskit.algorithms import AlgorithmResult
 from .phase_estimation_result import PhaseEstimationResult
 from .phase_estimation_scale import PhaseEstimationScale
 
 
-class HamiltonianPhaseEstimationResult:
+class HamiltonianPhaseEstimationResult(AlgorithmResult):
     """Store and manipulate results from running `HamiltonianPhaseEstimation`.
 
     This API of this class is nearly the same as `PhaseEstimatorResult`, differing only in
@@ -31,6 +30,7 @@ class HamiltonianPhaseEstimationResult:
 
     This class is meant to be instantiated via `HamiltonianPhaseEstimation.estimate`.
     """
+
     def __init__(self,
                  phase_estimation_result: PhaseEstimationResult,
                  phase_estimation_scale: PhaseEstimationScale,
@@ -45,6 +45,7 @@ def __init__(self,
                             coefficient must added to give correct eigenvalues.
                             This is done automatically when retrieving eigenvalues.
         """
+        super().__init__()
         self._phase_estimation_scale = phase_estimation_scale
         self._id_coefficient = id_coefficient
         self._phase_estimation_result = phase_estimation_result

qiskit/algorithms/phase_estimators/phase_estimation.py

@@ -60,7 +60,6 @@ def __init__(self,
                  num_evaluation_qubits: int,
                  quantum_instance: Optional[Union[QuantumInstance,
                                                   BaseBackend, Backend]] = None) -> None:
-
         """Args:
             num_evaluation_qubits: The number of qubits used in estimating the phase. The phase will
                 be estimated as a binary string with this many bits.
@@ -79,29 +78,46 @@ def __init__(self,
             quantum_instance = QuantumInstance(quantum_instance)
         self._quantum_instance = quantum_instance
 
-    def _add_classical_register(self) -> None:
-        """Add measurement instructions only if we are not using a state vector simulator."""
-        if not self._quantum_instance.is_statevector and not self._measurements_added:
-            # Measure only the evaluation qubits.
-            regname = 'meas'
-            circ = self._pe_circuit
-            creg = ClassicalRegister(self._num_evaluation_qubits, regname)
-            circ.add_register(creg)
-            circ.barrier()
-            circ.measure(range(self._num_evaluation_qubits), range(self._num_evaluation_qubits))
-            self._measurements_added = True
-
-    def construct_circuit(self) -> QuantumCircuit:
+    def construct_circuit(self,
+                          unitary: QuantumCircuit,
+                          state_preparation: Optional[QuantumCircuit] = None) -> QuantumCircuit:
         """Return the circuit to be executed to estimate phases.
 
         This circuit includes as sub-circuits the core phase estimation circuit,
         with the addition of the state-preparation circuit and possibly measurement instructions.
         """
-        self._add_classical_register()
-        return self._pe_circuit
+        num_evaluation_qubits = self._num_evaluation_qubits
+        num_unitary_qubits = unitary.num_qubits
+
+        pe_circuit = circuit.library.PhaseEstimation(num_evaluation_qubits, unitary)
+
+        if state_preparation is not None:
+            pe_circuit.compose(
+                state_preparation,
+                qubits=range(num_evaluation_qubits,
+                             num_evaluation_qubits + num_unitary_qubits),
+                inplace=True,
+                front=True)
 
-    def _compute_phases(self, circuit_result: Result) -> Union[numpy.ndarray,
-                                                               qiskit.result.Counts]:
+        self._add_measurement_if_required(pe_circuit)
+
+        return pe_circuit
+
+    def _add_measurement_if_required(self, pe_circuit):
+        if not self._quantum_instance.is_statevector:
+            # Measure only the evaluation qubits.
+            regname = 'meas'
+            creg = ClassicalRegister(self._num_evaluation_qubits, regname)
+            pe_circuit.add_register(creg)
+            pe_circuit.barrier()
+            pe_circuit.measure(range(self._num_evaluation_qubits),
+                               range(self._num_evaluation_qubits))
+
+        return circuit
+
+    def _compute_phases(self,
+                        num_unitary_qubits: int,
+                        circuit_result: Result) -> Union[numpy.ndarray, qiskit.result.Counts]:
         """Compute frequencies/counts of phases from the result of running the QPE circuit.
 
         How the frequencies are computed depends on whether the backend computes amplitude or
@@ -122,6 +138,7 @@ def _compute_phases(self, circuit_result: Result) -> Union[numpy.ndarray,
         they can be easily understood when displaying or plotting a histogram.
 
         Args:
+            num_unitary_qubits: The number of qubits in the unitary.
             circuit_result: the result object returned by the backend that ran the QPE circuit.
 
         Returns:
@@ -133,7 +150,7 @@ def _compute_phases(self, circuit_result: Result) -> Union[numpy.ndarray,
             evaluation_density_matrix = qiskit.quantum_info.partial_trace(
                 state_vec,
                 range(self._num_evaluation_qubits,
-                      self._num_evaluation_qubits + self._num_unitary_qubits)
+                      self._num_evaluation_qubits + num_unitary_qubits)
             )
             phases = evaluation_density_matrix.probabilities()
         else:
@@ -152,9 +169,9 @@ def estimate(self,
                  state_preparation: Optional[QuantumCircuit] = None,
                  pe_circuit: Optional[QuantumCircuit] = None,
                  num_unitary_qubits: Optional[int] = None) -> PhaseEstimationResult:
-        """Run the circuit and return and return `PhaseEstimationResult`.
+        """Run the the phase estimation algorithm.
 
-           Args:
+        Args:
             unitary: The circuit representing the unitary operator whose eigenvalues (via phase)
                 will be measured. Exactly one of `pe_circuit` and `unitary` must be passed.
             state_preparation: The circuit that prepares the state whose eigenphase will be
@@ -166,35 +183,28 @@ def estimate(self,
                 if `pe_circuit` is passed. This parameter will be set from `unitary`
                 if `unitary` is passed.
 
+        Raises:
+            ValueError: If both `pe_circuit` and `unitary` are passed.
+            ValueError: If neither `pe_circuit` nor `unitary` are passed.
+
         Returns:
                An instance of qiskit.algorithms.phase_estimator_result.PhaseEstimationResult.
         """
+        num_evaluation_qubits = self._num_evaluation_qubits
 
         if unitary is not None:
             if pe_circuit is not None:
                 raise ValueError('Only one of `pe_circuit` and `unitary` may be passed.')
-            self._num_unitary_qubits = unitary.num_qubits
-            self._pe_circuit = circuit.library.PhaseEstimation(self._num_evaluation_qubits, unitary)
-            self._measurements_added = False
+            pe_circuit = self.construct_circuit(unitary, state_preparation)
+            num_unitary_qubits = unitary.num_qubits
 
-        if pe_circuit is not None:
-            if unitary is not None:
-                raise ValueError('Only one of `pe_circuit` and `unitary` may be passed.')
-            self._pe_circuit = pe_circuit
-            self._measurements_added = False
+        elif pe_circuit is not None:
+            self._add_measurement_if_required(pe_circuit)
 
-        if num_unitary_qubits is not None:
-            self._num_unitary_qubits = num_unitary_qubits
-
-        if state_preparation is not None:
-            self._pe_circuit.compose(
-                state_preparation,
-                qubits=range(self._num_evaluation_qubits,
-                             self._num_evaluation_qubits + self._num_unitary_qubits),
-                inplace=True,
-                front=True)
+        else:
+            raise ValueError('One of `pe_circuit` and `unitary` must be passed.')
 
-        circuit_result = self._quantum_instance.execute(self.construct_circuit())
-        phases = self._compute_phases(circuit_result)
-        return PhaseEstimationResult(self._num_evaluation_qubits, circuit_result=circuit_result,
+        circuit_result = self._quantum_instance.execute(pe_circuit)
+        phases = self._compute_phases(num_unitary_qubits, circuit_result)
+        return PhaseEstimationResult(num_evaluation_qubits, circuit_result=circuit_result,
                                      phases=phases)

qiskit/algorithms/phase_estimators/phase_estimation_result.py

@@ -36,12 +36,21 @@ def __init__(self, num_evaluation_qubits: int,
             circuit_result: result object returned by method running circuit.
             phases: ndarray or dict of phases and frequencies determined by QPE.
         """
-        # int: number of qubits in phase-readout register
-        self._num_evaluation_qubits = num_evaluation_qubits
+        super().__init__()
 
         self._phases = phases
+        # int: number of qubits in phase-readout register
+        self._num_evaluation_qubits = num_evaluation_qubits
         self._circuit_result = circuit_result
 
+    @property
+    def phases(self) -> Union[numpy.ndarray, dict]:
+        """Return all phases and their frequencies computed by QPE.
+
+        This is an array or dict whose values correspond to weights on bit strings.
+        """
+        return self._phases
+
     @property
     def circuit_result(self) -> Result:
         """Return the result object returned by running the QPE circuit (on hardware or simulator).

qiskit/algorithms/phase_estimators/phase_estimation_scale.py

@@ -114,7 +114,7 @@ def scale_phases(self, phases: Union[List, Dict], id_coefficient: float = 0.0
         return phases
 
     @classmethod
-    def from_pauli_sum(cls, pauli_sum: SummedOp):
+    def from_pauli_sum(cls, pauli_sum: SummedOp) -> 'PhaseEstimationScale':
         """Create a PhaseEstimationScale from a `SummedOp` representing a sum of Pauli Operators.
 
         It is assumed that the `SummedOp` `pauli_sum` is the sum of `PauliOp`s. The bound on

qiskit/algorithms/phase_estimators/phase_estimator.py

@@ -12,15 +12,11 @@
 
 """The Phase Estimator interface"""
 
-from typing import Union, Optional
-from abc import ABC, abstractmethod
-import numpy
-import qiskit.circuit as circuit
+from typing import Optional
+from abc import ABC, abstractmethod, abstractproperty
 from qiskit.circuit import QuantumCircuit
 from qiskit.algorithms import AlgorithmResult
 
-# pylint: disable=attribute-defined-outside-init
-
 
 class PhaseEstimator(ABC):
     """The Phase Estimator interface
@@ -36,31 +32,19 @@ def estimate(self,
                  state_preparation: Optional[QuantumCircuit] = None,
                  pe_circuit: Optional[QuantumCircuit] = None,
                  num_unitary_qubits: Optional[int] = None) -> 'PhaseEstimatorResult':
-        """
-        Estimate the phase.
-        """
-
-        return PhaseEstimatorResult()
+        """Estimate the phase."""
+        raise NotImplementedError
 
 
 class PhaseEstimatorResult(AlgorithmResult):
     """Phase Estimator Result."""
 
-    @property
-    def phases(self) -> Union[numpy.ndarray, dict]:
-        """Return all phases and their frequencies computed by QPE.
-
-        This is an array or dict whose values correspond to weights on bit strings.
-        """
-        # pylint: disable=no-member
-        return self._phases
-
-    @property
+    @abstractproperty
     def most_likely_phase(self) -> float:
         r"""Return the estimated phase as a number in :math:`[0.0, 1.0)`.
 
         1.0 corresponds to a phase of :math:`2\pi`. It is assumed that the input vector is an
         eigenvector of the unitary so that the peak of the probability density occurs at the bit
         string that most closely approximates the true phase.
         """
-        raise NotImplementedError()
+        raise NotImplementedError