Porting to QubitDevice (#38)

* Fixes the plugin to support pyquil 2.16

* bump pyQuil version

* added rewiring

* Fixed qubit measuring order

* store compiled program

* remove partial rewiring

* return partial

* fix wiring key error

* First draft of porting to QubitDevice

* Modifying BasisState such that it uses the wires passed (instead of enumerate); rename CCNOT to Toffoli; fix test for test_qvm

* Fixing qvm and qpu and their tests

* Fixes for the wavelength simulator and tests

* Fixing WavefunctionDevice and NumpyWaveFunction Device; added functionality to raise an error for QVM analytic=True case (test included)

* Increasing shots for QPU tests

* Modifying comments, removing active_wires from test

* Update pennylane_forest/qpu.py

Co-Authored-By: Josh Izaac <josh146@gmail.com>

* Implementing feedback

* Docstring

* Using keyword arguments for apply, fixing matrices that were used previously from default.qubit, increasing number of shots such that the stochastic tests will more likely pass

* Modifying the requirements file for the checks

* Decreasing number of shots and adding flaky tests to stochastic test cases (each is tried 10 times and they need to pass at least once)

* Rename apply_wiring to remap_wires

* Remove no longer used instance attribute _eigs from QVMDevice

* Refactor QVMDevice such that generate_samples returns samples instead of returning and assigning it to dev._samples; adjusting the tests accordingly

* Remove unused imports

* Modfying qpu test to be a parametrized test

* Adjust flaky runs

* Remove unnecessary else

* Restricting the qpu lattices that can be used

* Decresase the number of required pass runs for test

Co-authored-by: Josh Izaac <josh146@gmail.com>

pennylane_forest/device.py

@@ -44,7 +44,7 @@
 # following gates are not supported by PennyLane
 from pyquil.gates import S, T, CPHASE00, CPHASE01, CPHASE10, CPHASE, CCNOT, CSWAP, ISWAP, PSWAP
 
-from pennylane import Device
+from pennylane import QubitDevice
 
 from ._version import __version__
 
@@ -64,7 +64,7 @@ def basis_state(par, *wires):
         list: list of PauliX matrix operators acting on each wire
     """
     # pylint: disable=unused-argument
-    return [X(w) for w, p in enumerate(par) if p == 1]
+    return [X(w) for w, p in zip(wires, par) if p == 1]
 
 
 def qubit_unitary(par, *wires):
@@ -151,15 +151,15 @@ def controlled_phase(phi, q, *wires):
     # the following gates are provided by the PL-Forest plugin
     "S": S,
     "T": T,
-    "CCNOT": CCNOT,
+    "Toffoli": CCNOT,
     "CPHASE": controlled_phase,
     "CSWAP": CSWAP,
     "ISWAP": ISWAP,
     "PSWAP": PSWAP,
 }
 
 
-class ForestDevice(Device):
+class ForestDevice(QubitDevice):
     r"""Abstract Forest device for PennyLane.
 
     Args:
@@ -184,10 +184,10 @@ class ForestDevice(Device):
     author = "Josh Izaac"
 
     _operation_map = pyquil_operation_map
-    _capabilities = {"model": "qubit"}
+    _capabilities = {"model": "qubit", "tensor_observables": True}
 
     def __init__(self, wires, shots=1000, analytic=False,  **kwargs):
-        super().__init__(wires, shots)
+        super().__init__(wires, shots, analytic=analytic)
         self.analytic = analytic
         self.forest_url = kwargs.get("forest_url", pyquil_config.forest_url)
         self.qvm_url = kwargs.get("qvm_url", pyquil_config.qvm_url)
@@ -224,31 +224,46 @@ def program(self):
         """View the last evaluated Quil program"""
         return self.prog
 
-    def apply(self, operation, wires, par):
-        # pylint: disable=attribute-defined-outside-init
+    def remap_wires(self, wires):
+        """Use the wiring specified for the device if applicable.
+
+        Returns:
+            list: wires as integers corresponding to the wiring if applicable
+        """
         if hasattr(self, "wiring"):
-            qubits = [int(self.wiring[i]) for i in wires]
-        else:
-            qubits = [int(w) for w in wires]
+            return [int(self.wiring[i]) for i in wires]
+
+        return [int(w) for w in wires]
+
+    def apply(self, operations, **kwargs):
+        # pylint: disable=attribute-defined-outside-init
+        rotations = kwargs.get("rotations", [])
+
+        # Storing the active wires
+        self._active_wires = ForestDevice.active_wires(operations + rotations)
+
+        # Apply the circuit operations
+        for i, operation in enumerate(operations):
+            # number of wires on device
+            wires = self.remap_wires(operation.wires)
+            par = operation.parameters
 
-        self.prog += self._operation_map[operation](*par, *qubits)
+            if i > 0 and operation.name in ("QubitStateVector", "BasisState"):
+                raise DeviceError("Operation {} cannot be used after other Operations have already been applied "
+                                  "on a {} device.".format(operation.name, self.short_name))
 
-        # keep track of the active wires. This is required, as the
-        # pyQuil wavefunction simulator creates qubits dynamically.
-        if wires:
-            self.active_wires = self.active_wires.union(set(wires))
-        else:
-            self.active_wires = set(range(self.num_wires))
+            self.prog += self._operation_map[operation.name](*par, *wires)
 
-    @abc.abstractmethod
-    def pre_measure(self):  # pragma no cover
-        """Run the QVM or QPU"""
-        raise NotImplementedError
+        # Apply the circuit rotations
+        for operation in rotations:
+            wires = self.remap_wires(operation.wires)
+            par = operation.parameters
+            self.prog += self._operation_map[operation.name](*par, *wires)
 
     def reset(self):
         self.prog = Program()
-        self.active_wires = set()
-        self.state = None
+        self._active_wires = set()
+        self._state = None
 
     @property
     def operations(self):
@@ -327,3 +342,29 @@ def mat_vec_product(self, mat, vec, wires):
         state_multi_index = np.transpose(tdot, inv_perm)
 
         return np.reshape(state_multi_index, 2 ** self.num_wires)
+
+    def probability(self, wires=None):
+        """Return the (marginal) probability of each computational basis
+        state from the last run of the device.
+
+        If no wires are specified, then all the basis states representable by
+        the device are considered and no marginalization takes place.
+
+        .. warning:: This method will have to be redefined for hardware devices, since it uses
+            the ``device._state`` attribute. This attribute might not be available for such devices.
+
+        Args:
+            wires (Sequence[int]): Sequence of wires to return
+                marginal probabilities for. Wires not provided
+                are traced out of the system.
+
+        Returns:
+            List[float]: list of the probabilities
+        """
+        if self._state is None:
+            return None
+
+        wires = wires or range(self.num_wires)
+        wires = self.remap_wires(wires)
+        prob = self.marginal_prob(np.abs(self._state) ** 2, wires)
+        return prob

pennylane_forest/numpy_wavefunction.py

@@ -26,11 +26,11 @@
 from pyquil.pyqvm import PyQVM
 from pyquil.simulation import NumpyWavefunctionSimulator
 
-from .wavefunction import WavefunctionDevice
+from .device import ForestDevice
 from ._version import __version__
 
 
-class NumpyWavefunctionDevice(WavefunctionDevice):
+class NumpyWavefunctionDevice(ForestDevice):
     r"""NumpyWavefunction simulator device for PennyLane.
 
     Args:
@@ -44,18 +44,18 @@ class NumpyWavefunctionDevice(WavefunctionDevice):
     observables = {"PauliX", "PauliY", "PauliZ", "Hadamard", "Hermitian", "Identity"}
 
     def __init__(self, wires, *, shots=1000, analytic=True, **kwargs):
-        super(WavefunctionDevice, self).__init__(wires, shots, analytic, **kwargs)
+        super().__init__(wires, shots, analytic, **kwargs)
         self.qc = PyQVM(n_qubits=wires, quantum_simulator_type=NumpyWavefunctionSimulator)
-        self.state = None
+        self._state = None
 
-    def pre_apply(self):
+    def apply(self, operations, **kwargs):
         self.reset()
         self.qc.wf_simulator.reset()
+        super().apply(operations, **kwargs)
 
-    def pre_measure(self):
         # TODO: currently, the PyQVM considers qubit 0 as the leftmost bit and therefore
         # returns amplitudes in the opposite of the Rigetti Lisp QVM (which considers qubit
         # 0 as the rightmost bit). This may change in the future, so in the future this
         # might need to get udpated to be similar to the pre_measure function of
         # pennylane_forest/wavefunction.py
-        self.state = self.qc.execute(self.prog).wf_simulator.wf.flatten()
+        self._state = self.qc.execute(self.prog).wf_simulator.wf.flatten()

pennylane_forest/qpu.py

@@ -104,14 +104,8 @@ def __init__(self, device, *, shots=1024, active_reset=True, load_qc=True, reado
         self.calibrate_readout = calibrate_readout
         self.wiring = {i: q for i, q in enumerate(self.qc.qubits())}
 
-    def pre_rotations(self, observable, wires):
-        """
-        This is used in `QVM.pre_measure`. Since the pre-rotations are handled
-        by `measure_observables` (see `expval`), we simply don't do anything here
-        """
-        pass
-
-    def expval(self, observable, wires, par):
+    def expval(self, observable):
+        wires = observable.wires
         # Single-qubit observable
         if len(wires) == 1:
             # identify Experiment Settings for each of the possible single-qubit observables
@@ -126,7 +120,7 @@ def expval(self, observable, wires, par):
                              ExperimentSetting(TensorProductState(), float(np.sqrt(1/2)) * sZ(qubit))]
             }
             # expectation values for single-qubit observables
-            if observable in ["PauliX", "PauliY", "PauliZ", "Identity", "Hadamard"]:
+            if observable.name in ["PauliX", "PauliY", "PauliZ", "Identity", "Hadamard"]:
                 prep_prog = Program()
                 for instr in self.program.instructions:
                     if isinstance(instr, Gate):
@@ -142,32 +136,20 @@ def expval(self, observable, wires, par):
                 if self.readout_error is not None:
                     prep_prog.define_noisy_readout(qubit, p00=self.readout_error[0],
                                                           p11=self.readout_error[1])
-                tomo_expt = Experiment(settings=d_expt_settings[observable], program=prep_prog)
+
+                # All observables are rotated and can be measured in the PauliZ basis
+                tomo_expt = Experiment(settings=d_expt_settings["PauliZ"], program=prep_prog)
                 grouped_tomo_expt = group_experiments(tomo_expt)
                 meas_obs = list(measure_observables(self.qc, grouped_tomo_expt,
                                                     active_reset=self.active_reset,
                                                     symmetrize_readout=self.symmetrize_readout,
                                                     calibrate_readout=self.calibrate_readout))
                 return np.sum([expt_result.expectation for expt_result in meas_obs])
 
-            elif observable == 'Hermitian':
+            elif observable.name == 'Hermitian':
                 # <H> = \sum_i w_i p_i
                 Hkey = tuple(par[0].flatten().tolist())
                 w = self._eigs[Hkey]['eigval']
                 return w[0]*p0 + w[1]*p1
 
-        # Multi-qubit observable
-        # ----------------------
-        # Currently, we only support qml.expval.Hermitian(A, wires),
-        # where A is a 2^N x 2^N matrix acting on N wires.
-        #
-        # Eventually, we will also support tensor products of Pauli
-        # matrices in the PennyLane UI.
-
-        probs = self.probabilities(wires)
-
-        if observable == 'Hermitian':
-            Hkey = tuple(par[0].flatten().tolist())
-            w = self._eigs[Hkey]['eigval']
-            # <A> = \sum_i w_i p_i
-            return w @ probs
+        return super().expval(observable)