Repository: PennyLaneAI/pennylane
Pull request #2180 (2022-02-08 11:07:54)

Add `parameter_frequencies` to operations.
This is the first PR in a series to introduce the general parameter-shift functionalities to `qml.gradients.param_shift`.

Here we only add the operation property `parameter_frequencies` as well as the specific frequencies to all operations where applicable.

github-actions[bot]: Hello. You may have forgotten to update the changelog!
Please edit [`doc/releases/changelog-dev.md`](/PennyLaneAI/pennylane/blob/master/doc/releases/changelog-dev.md) with:
* A one-to-two sentence description of the change. You may include a small working example for new features.
* A link back to this PR.
* Your name (or GitHub username) in the contributors section.
codecov[bot]: # [Codecov](https://codecov.io/gh/PennyLaneAI/pennylane/pull/2180?src=pr&el=h1&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI) Report
> Merging [#2180](https://codecov.io/gh/PennyLaneAI/pennylane/pull/2180?src=pr&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI) (5e7c9f0) into [master](https://codecov.io/gh/PennyLaneAI/pennylane/commit/c6d5cc4dfc1dbe15050c53301d8c3b245a8674ef?el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI) (c6d5cc4) will **decrease** coverage by `2.85%`.
> The diff coverage is `81.03%`.

[![Impacted file tree graph](https://codecov.io/gh/PennyLaneAI/pennylane/pull/2180/graphs/tree.svg?width=650&height=150&src=pr&token=9S3YzlUgGk&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI)](https://codecov.io/gh/PennyLaneAI/pennylane/pull/2180?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI)

```diff
@@            Coverage Diff             @@
##           master    #2180      +/-   ##
==========================================
- Coverage   99.20%   96.35%   -2.86%     
==========================================
  Files         231      231              
  Lines       17784    17842      +58     
==========================================
- Hits        17643    17191     -452     
- Misses        141      651     +510     
```


| [Impacted Files](https://codecov.io/gh/PennyLaneAI/pennylane/pull/2180?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI) | Coverage Δ | |
|---|---|---|
| [pennylane/operation.py](https://codecov.io/gh/PennyLaneAI/pennylane/pull/2180/diff?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI#diff-cGVubnlsYW5lL29wZXJhdGlvbi5weQ==) | `94.91% <25.00%> (-1.04%)` | :arrow_down: |
| [pennylane/ops/qubit/parametric\_ops.py](https://codecov.io/gh/PennyLaneAI/pennylane/pull/2180/diff?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI#diff-cGVubnlsYW5lL29wcy9xdWJpdC9wYXJhbWV0cmljX29wcy5weQ==) | `99.16% <88.88%> (-0.66%)` | :arrow_down: |
| [pennylane/ops/qubit/qchem\_ops.py](https://codecov.io/gh/PennyLaneAI/pennylane/pull/2180/diff?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI#diff-cGVubnlsYW5lL29wcy9xdWJpdC9xY2hlbV9vcHMucHk=) | `99.54% <92.85%> (-0.46%)` | :arrow_down: |
| [pennylane/interfaces/batch/jax.py](https://codecov.io/gh/PennyLaneAI/pennylane/pull/2180/diff?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI#diff-cGVubnlsYW5lL2ludGVyZmFjZXMvYmF0Y2gvamF4LnB5) | `0.00% <0.00%> (-100.00%)` | :arrow_down: |
| [pennylane/interfaces/batch/jax\_jit.py](https://codecov.io/gh/PennyLaneAI/pennylane/pull/2180/diff?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI#diff-cGVubnlsYW5lL2ludGVyZmFjZXMvYmF0Y2gvamF4X2ppdC5weQ==) | `0.00% <0.00%> (-100.00%)` | :arrow_down: |
| [pennylane/interfaces/jax.py](https://codecov.io/gh/PennyLaneAI/pennylane/pull/2180/diff?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI#diff-cGVubnlsYW5lL2ludGVyZmFjZXMvamF4LnB5) | `0.00% <0.00%> (-97.44%)` | :arrow_down: |
| [pennylane/devices/default\_qubit\_jax.py](https://codecov.io/gh/PennyLaneAI/pennylane/pull/2180/diff?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI#diff-cGVubnlsYW5lL2RldmljZXMvZGVmYXVsdF9xdWJpdF9qYXgucHk=) | `0.00% <0.00%> (-96.63%)` | :arrow_down: |
| [pennylane/qnn/cost.py](https://codecov.io/gh/PennyLaneAI/pennylane/pull/2180/diff?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI#diff-cGVubnlsYW5lL3Fubi9jb3N0LnB5) | `42.85% <0.00%> (-57.15%)` | :arrow_down: |
| [pennylane/transforms/classical\_jacobian.py](https://codecov.io/gh/PennyLaneAI/pennylane/pull/2180/diff?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI#diff-cGVubnlsYW5lL3RyYW5zZm9ybXMvY2xhc3NpY2FsX2phY29iaWFuLnB5) | `77.08% <0.00%> (-22.92%)` | :arrow_down: |
| [pennylane/math/utils.py](https://codecov.io/gh/PennyLaneAI/pennylane/pull/2180/diff?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI#diff-cGVubnlsYW5lL21hdGgvdXRpbHMucHk=) | `84.41% <0.00%> (-15.59%)` | :arrow_down: |
| ... and [16 more](https://codecov.io/gh/PennyLaneAI/pennylane/pull/2180/diff?src=pr&el=tree-more&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI) | |

------

[Continue to review full report at Codecov](https://codecov.io/gh/PennyLaneAI/pennylane/pull/2180?src=pr&el=continue&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI).
> **Legend** - [Click here to learn more](https://docs.codecov.io/docs/codecov-delta?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI)
> `Δ = absolute <relative> (impact)`, `ø = not affected`, `? = missing data`
> Powered by [Codecov](https://codecov.io/gh/PennyLaneAI/pennylane/pull/2180?src=pr&el=footer&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI). Last update [c6d5cc4...5e7c9f0](https://codecov.io/gh/PennyLaneAI/pennylane/pull/2180?src=pr&el=lastupdated&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI). Read the [comment docs](https://docs.codecov.io/docs/pull-request-comments?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI).

josh146: [sc-13053]

Patch:
diff --git a/doc/releases/changelog-dev.md b/doc/releases/changelog-dev.md
index afed6a8ec..c7ee27f1d 100644
--- a/doc/releases/changelog-dev.md
+++ b/doc/releases/changelog-dev.md
@@ -4,6 +4,17 @@
 
 <h3>New features since last release</h3>
 
+* Parametric operations now have the `parameter_frequencies`
+  method that returns the frequencies with which a parameter
+  enters a circuit when using the operation.
+  [(#2180)](https://github.com/PennyLaneAI/pennylane/pull/2180)
+
+  The frequencies can be used for circuit analysis, optimization
+  via the `RotosolveOptimizer` and differentiation with the
+  parameter-shift rule. They assume that the circuit returns
+  expectation values or probabilities, for a variance
+  measurement the frequencies will differ.
+
 * Continued development of the circuit-cutting compiler:
 
   A method for converting a quantum tape to a directed multigraph that is amenable
@@ -34,6 +45,11 @@
 
 <h3>Bug fixes</h3>
 
+* The operation `OrbitalRotation` previously was wrongfully registered to satisfy
+  the four-term parameter shift rule, it now will be decomposed instead when
+  using the parameter-shift rule.
+  [(#2180)](https://github.com/PennyLaneAI/pennylane/pull/2180)
+
 <h3>Documentation</h3>
 
 * The ``pennylane.numpy`` subpackage is now included in the PennyLane
diff --git a/pennylane/operation.py b/pennylane/operation.py
index c643ff7b3..8f262d4c4 100644
--- a/pennylane/operation.py
+++ b/pennylane/operation.py
@@ -1,5 +1,4 @@
 # Copyright 2018-2021 Xanadu Quantum Technologies Inc.
-
 # Licensed under the Apache License, Version 2.0 (the "License");
 # you may not use this file except in compliance with the License.
 # You may obtain a copy of the License at
@@ -229,6 +228,16 @@ def classproperty(func):
     return ClassPropertyDescriptor(func)
 
 
+# =============================================================================
+# Error classes
+# =============================================================================
+
+
+class OperatorPropertyUndefined(Exception):
+    """Generic exception to be used for undefined
+    Operator properties or methods."""
+
+
 # =============================================================================
 # Base Operator class
 # =============================================================================
@@ -677,7 +686,7 @@ class Operation(Operator):
         multiplier = 0.5 / np.sin(shift)
         a = 1
 
-        # We set the default recipe following:
+        # We set the following default recipe:
         # ∂f(x) = c*f(a*x+s) - c*f(a*x-s)
         # where we express a positive and a negative shift by default
         default_param_shift = [[multiplier, a, shift], [-multiplier, a, -shift]]
@@ -691,7 +700,7 @@ class Operation(Operator):
         A length-2 list ``[generator, scaling_factor]``, where
 
         * ``generator`` is an existing PennyLane
-          operation class or :math:`2\times 2` Hermitian array
+          operation class or a Hermitian array
           that acts as the generator of the current operation
 
         * ``scaling_factor`` represents a scaling factor applied
@@ -709,6 +718,57 @@ class Operation(Operator):
         """
         return [None, 1]
 
+    @property
+    def parameter_frequencies(self):
+        r"""Returns the frequencies for each operator parameter with respect
+        to an expectation value of the form
+        :math:`\langle \psi | U(\mathbf{p})^\dagger \hat{O} U(\mathbf{p})|\psi\rangle`.
+
+        These frequencies encode the behaviour of the operator :math:`U(\mathbf{p})`
+        on the value of the expectation value as the parameters are modified.
+        For more details, please see the :mod:`.pennylane.fourier` module.
+
+        Returns:
+            list[tuple[int or float]]: Tuple of frequencies for each parameter.
+            Note that only non-negative frequency values are returned.
+
+        **Example**
+
+        >>> op = qml.CRot(0.4, 0.1, 0.3, wires=[0, 1])
+        >>> op.parameter_frequencies
+        [(0.5, 1), (0.5, 1), (0.5, 1)]
+
+        For operators that define a generator, the parameter frequencies are directly
+        related to the eigenvalues of the generator:
+
+        >>> op = qml.ControlledPhaseShift(0.1, wires=[0, 1])
+        >>> op.parameter_frequencies
+        [(1,)]
+        >>> gen_eigvals = tuple(np.linalg.eigvals(op.generator[0] * op.generator[1]))
+        >>> qml.gradients.eigvals_to_frequencies(gen_eigvals)
+        (tensor(1., requires_grad=True),)
+
+        For more details on this relationship, see :func:`.eigvals_to_frequencies`.
+        """
+        if self.num_params == 1:
+            # if the operator has a single parameter, we can query the
+            # generator, and if defined, use its eigenvalues.
+            op, coeff = self.generator
+            if op is None:
+                raise OperatorPropertyUndefined(
+                    f"Operation {self.name} does not have parameter frequencies."
+                )
+
+            if not isinstance(op, np.ndarray):
+                op = op.matrix
+            gen_eigvals = tuple(np.linalg.eigvals(op))
+            coeff = np.abs(coeff)
+            return tuple(coeff * val for val in qml.gradients.eigvals_to_frequencies(gen_eigvals))
+
+        raise OperatorPropertyUndefined(
+            f"Operation {self.name} does not have parameter frequencies."
+        )
+
     @property
     def inverse(self):
         """Boolean determining if the inverse of the operation was requested."""
diff --git a/pennylane/ops/qubit/parametric_ops.py b/pennylane/ops/qubit/parametric_ops.py
index f5399e81a..8bb1563af 100644
--- a/pennylane/ops/qubit/parametric_ops.py
+++ b/pennylane/ops/qubit/parametric_ops.py
@@ -55,6 +55,7 @@ class RX(Operation):
     basis = "X"
     grad_method = "A"
     generator = [PauliX, -1 / 2]
+    parameter_frequencies = [(1,)]
 
     @classmethod
     def _matrix(cls, *params):
@@ -109,6 +110,7 @@ class RY(Operation):
     basis = "Y"
     grad_method = "A"
     generator = [PauliY, -1 / 2]
+    parameter_frequencies = [(1,)]
 
     @classmethod
     def _matrix(cls, *params):
@@ -157,6 +159,7 @@ class RZ(Operation):
     basis = "Z"
     grad_method = "A"
     generator = [PauliZ, -1 / 2]
+    parameter_frequencies = [(1,)]
 
     @classmethod
     def _matrix(cls, *params):
@@ -216,6 +219,7 @@ class PhaseShift(Operation):
     basis = "Z"
     grad_method = "A"
     generator = [np.array([[0, 0], [0, 1]]), 1]
+    parameter_frequencies = [(1,)]
 
     def label(self, decimals=None, base_label=None):
         return super().label(decimals=decimals, base_label=base_label or "Rϕ")
@@ -287,6 +291,7 @@ class ControlledPhaseShift(Operation):
     basis = "Z"
     grad_method = "A"
     generator = [np.array([[0, 0, 0, 0], [0, 0, 0, 0], [0, 0, 0, 0], [0, 0, 0, 1]]), 1]
+    parameter_frequencies = [(1,)]
 
     def label(self, decimals=None, base_label=None):
         return super().label(decimals=decimals, base_label=base_label or "Rϕ")
@@ -368,6 +373,7 @@ class Rot(Operation):
     num_wires = 1
     num_params = 3
     grad_method = "A"
+    parameter_frequencies = [(1,), (1,), (1,)]
 
     @classmethod
     def _matrix(cls, *params):
@@ -449,6 +455,7 @@ class MultiRZ(Operation):
     num_wires = AnyWires
     num_params = 1
     grad_method = "A"
+    parameter_frequencies = [(1,)]
 
     @classmethod
     def _matrix(cls, theta, n):
@@ -629,6 +636,8 @@ class PauliRot(Operation):
         """
         return all(pauli in PauliRot._ALLOWED_CHARACTERS for pauli in pauli_word)
 
+    parameter_frequencies = [(1,)]
+
     @classmethod
     def _matrix(cls, *params):
         pauli_word = params[1]
@@ -822,6 +831,7 @@ class CRX(Operation):
         np.array([[0, 0, 0, 0], [0, 0, 0, 0], [0, 0, 0, 1], [0, 0, 1, 0]]),
         -1 / 2,
     ]
+    parameter_frequencies = [(0.5, 1.0)]
 
     def label(self, decimals=None, base_label=None):
         return super().label(decimals=decimals, base_label=base_label or "RX")
@@ -920,6 +930,7 @@ class CRY(Operation):
         np.array([[0, 0, 0, 0], [0, 0, 0, 0], [0, 0, 0, -1j], [0, 0, 1j, 0]]),
         -1 / 2,
     ]
+    parameter_frequencies = [(0.5, 1.0)]
 
     def label(self, decimals=None, base_label=None):
         return super().label(decimals=decimals, base_label=base_label or "RY")
@@ -1012,6 +1023,7 @@ class CRZ(Operation):
         np.array([[0, 0, 0, 0], [0, 0, 0, 0], [0, 0, 1, 0], [0, 0, 0, -1]]),
         -1 / 2,
     ]
+    parameter_frequencies = [(0.5, 1.0)]
 
     def label(self, decimals=None, base_label=None):
         return super().label(decimals=decimals, base_label=base_label or "RZ")
@@ -1097,6 +1109,7 @@ class CRot(Operation):
     num_params = 3
     grad_method = "A"
     grad_recipe = four_term_grad_recipe * 3
+    parameter_frequencies = [(0.5, 1.0), (0.5, 1.0), (0.5, 1.0)]
 
     def label(self, decimals=None, base_label=None):
         return super().label(decimals=decimals, base_label=base_label or "Rot")
@@ -1189,6 +1202,7 @@ class U1(Operation):
     num_params = 1
     grad_method = "A"
     generator = [np.array([[0, 0], [0, 1]]), 1]
+    parameter_frequencies = [(1,)]
 
     @classmethod
     def _matrix(cls, *params):
@@ -1246,6 +1260,7 @@ class U2(Operation):
     num_wires = 1
     num_params = 2
     grad_method = "A"
+    parameter_frequencies = [(1,), (1,)]
 
     @classmethod
     def _matrix(cls, *params):
@@ -1319,6 +1334,7 @@ class U3(Operation):
     num_wires = 1
     num_params = 3
     grad_method = "A"
+    parameter_frequencies = [(1,), (1,), (1,)]
 
     @classmethod
     def _matrix(cls, *params):
@@ -1392,6 +1408,7 @@ class IsingXX(Operation):
         np.array([[0, 0, 0, 1], [0, 0, 1, 0], [0, 1, 0, 0], [1, 0, 0, 0]]),
         -1 / 2,
     ]
+    parameter_frequencies = [(1,)]
 
     @classmethod
     def _matrix(cls, *params):
@@ -1452,6 +1469,7 @@ class IsingYY(Operation):
         np.array([[0, 0, 0, -1], [0, 0, 1, 0], [0, 1, 0, 0], [-1, 0, 0, 0]]),
         -1 / 2,
     ]
+    parameter_frequencies = [(1,)]
 
     @staticmethod
     def decomposition(phi, wires):
@@ -1510,6 +1528,7 @@ class IsingZZ(Operation):
         np.array([[1, 0, 0, 0], [0, -1, 0, 0], [0, 0, -1, 0], [0, 0, 0, 1]]),
         -1 / 2,
     ]
+    parameter_frequencies = [(1,)]
 
     @staticmethod
     def decomposition(phi, wires):
diff --git a/pennylane/ops/qubit/qchem_ops.py b/pennylane/ops/qubit/qchem_ops.py
index 46a47d6e6..d3dddd0bd 100644
--- a/pennylane/ops/qubit/qchem_ops.py
+++ b/pennylane/ops/qubit/qchem_ops.py
@@ -85,6 +85,7 @@ class SingleExcitation(Operation):
         np.array([[0, 0, 0, 0], [0, 0, -1j, 0], [0, 1j, 0, 0], [0, 0, 0, 0]]),
         -1 / 2,
     ]
+    parameter_frequencies = [(0.5, 1.0)]
 
     @classmethod
     def _matrix(cls, *params):
@@ -114,23 +115,23 @@ class SingleExcitation(Operation):
         return super().label(decimals=decimals, base_label=base_label or "G")
 
 
-class SingleExcitationMinus(Operation):
-    r"""SingleExcitationMinus(phi, wires)
-    Single excitation rotation with negative phase-shift outside the rotation subspace.
+class SingleExcitationPlus(Operation):
+    r"""SingleExcitationPlus(phi, wires)
+    Single excitation rotation with positive phase-shift outside the rotation subspace.
 
-    .. math:: U_-(\phi) = \begin{bmatrix}
-                e^{-i\phi/2} & 0 & 0 & 0 \\
+    .. math:: U_+(\phi) = \begin{bmatrix}
+                e^{i\phi/2} & 0 & 0 & 0 \\
                 0 & \cos(\phi/2) & -\sin(\phi/2) & 0 \\
                 0 & \sin(\phi/2) & \cos(\phi/2) & 0 \\
-                0 & 0 & 0 & e^{-i\phi/2}
+                0 & 0 & 0 & e^{i\phi/2}
             \end{bmatrix}.
 
     **Details:**
 
     * Number of wires: 2
     * Number of parameters: 1
-    * Gradient recipe: :math:`\frac{d}{d\phi}f(U_-(\phi)) = \frac{1}{2}\left[f(U_-(\phi+\pi/2)) - f(U_-(\phi-\pi/2))\right]`
-      where :math:`f` is an expectation value depending on :math:`U_-(\phi)`.
+    * Gradient recipe: :math:`\frac{d}{d\phi}f(U_+(\phi)) = \frac{1}{2}\left[f(U_+(\phi+\pi/2)) - f(U_+(\phi-\pi/2))\right]`
+      where :math:`f` is an expectation value depending on :math:`U_+(\phi)`.
 
     Args:
         phi (float): rotation angle :math:`\phi`
@@ -141,9 +142,10 @@ class SingleExcitationMinus(Operation):
     num_params = 1
     grad_method = "A"
     generator = [
-        np.array([[1, 0, 0, 0], [0, 0, -1j, 0], [0, 1j, 0, 0], [0, 0, 0, 1]]),
+        np.array([[-1, 0, 0, 0], [0, 0, -1j, 0], [0, 1j, 0, 0], [0, 0, 0, -1]]),
         -1 / 2,
     ]
+    parameter_frequencies = [(1,)]
 
     @classmethod
     def _matrix(cls, *params):
@@ -159,7 +161,7 @@ class SingleExcitationMinus(Operation):
             c = qml.math.cast_like(c, 1j)
             s = qml.math.cast_like(s, 1j)
 
-        e = qml.math.exp(-1j * theta / 2)
+        e = qml.math.exp(1j * theta / 2)
         mat = qml.math.diag([e, 0, 0, e]) + qml.math.diag([0, c, c, 0])
         off_diag = qml.math.convert_like(np.diag([0, 1, -1, 0])[::-1].copy(), theta)
         return mat + s * qml.math.cast_like(off_diag, s)
@@ -169,10 +171,10 @@ class SingleExcitationMinus(Operation):
         decomp_ops = [
             qml.PauliX(wires=wires[0]),
             qml.PauliX(wires=wires[1]),
-            qml.ControlledPhaseShift(-theta / 2, wires=[wires[1], wires[0]]),
+            qml.ControlledPhaseShift(theta / 2, wires=[wires[1], wires[0]]),
             qml.PauliX(wires=wires[0]),
             qml.PauliX(wires=wires[1]),
-            qml.ControlledPhaseShift(-theta / 2, wires=[wires[0], wires[1]]),
+            qml.ControlledPhaseShift(theta / 2, wires=[wires[0], wires[1]]),
             qml.CNOT(wires=[wires[0], wires[1]]),
             qml.CRY(theta, wires=[wires[1], wires[0]]),
             qml.CNOT(wires=[wires[0], wires[1]]),
@@ -181,29 +183,29 @@ class SingleExcitationMinus(Operation):
 
     def adjoint(self):
         (phi,) = self.parameters
-        return SingleExcitationMinus(-phi, wires=self.wires)
+        return SingleExcitationPlus(-phi, wires=self.wires)
 
     def label(self, decimals=None, base_label=None):
-        return super().label(decimals=decimals, base_label=base_label or "G₋")
+        return super().label(decimals=decimals, base_label=base_label or "G₊")
 
 
-class SingleExcitationPlus(Operation):
-    r"""SingleExcitationPlus(phi, wires)
-    Single excitation rotation with positive phase-shift outside the rotation subspace.
+class SingleExcitationMinus(Operation):
+    r"""SingleExcitationMinus(phi, wires)
+    Single excitation rotation with negative phase-shift outside the rotation subspace.
 
-    .. math:: U_+(\phi) = \begin{bmatrix}
-                e^{i\phi/2} & 0 & 0 & 0 \\
+    .. math:: U_-(\phi) = \begin{bmatrix}
+                e^{-i\phi/2} & 0 & 0 & 0 \\
                 0 & \cos(\phi/2) & -\sin(\phi/2) & 0 \\
                 0 & \sin(\phi/2) & \cos(\phi/2) & 0 \\
-                0 & 0 & 0 & e^{i\phi/2}
+                0 & 0 & 0 & e^{-i\phi/2}
             \end{bmatrix}.
 
     **Details:**
 
     * Number of wires: 2
     * Number of parameters: 1
-    * Gradient recipe: :math:`\frac{d}{d\phi}f(U_+(\phi)) = \frac{1}{2}\left[f(U_+(\phi+\pi/2)) - f(U_+(\phi-\pi/2))\right]`
-      where :math:`f` is an expectation value depending on :math:`U_+(\phi)`.
+    * Gradient recipe: :math:`\frac{d}{d\phi}f(U_-(\phi)) = \frac{1}{2}\left[f(U_-(\phi+\pi/2)) - f(U_-(\phi-\pi/2))\right]`
+      where :math:`f` is an expectation value depending on :math:`U_-(\phi)`.
 
     Args:
         phi (float): rotation angle :math:`\phi`
@@ -214,9 +216,10 @@ class SingleExcitationPlus(Operation):
     num_params = 1
     grad_method = "A"
     generator = [
-        np.array([[-1, 0, 0, 0], [0, 0, -1j, 0], [0, 1j, 0, 0], [0, 0, 0, -1]]),
+        np.array([[1, 0, 0, 0], [0, 0, -1j, 0], [0, 1j, 0, 0], [0, 0, 0, 1]]),
         -1 / 2,
     ]
+    parameter_frequencies = [(1,)]
 
     @classmethod
     def _matrix(cls, *params):
@@ -232,7 +235,7 @@ class SingleExcitationPlus(Operation):
             c = qml.math.cast_like(c, 1j)
             s = qml.math.cast_like(s, 1j)
 
-        e = qml.math.exp(1j * theta / 2)
+        e = qml.math.exp(-1j * theta / 2)
         mat = qml.math.diag([e, 0, 0, e]) + qml.math.diag([0, c, c, 0])
         off_diag = qml.math.convert_like(np.diag([0, 1, -1, 0])[::-1].copy(), theta)
         return mat + s * qml.math.cast_like(off_diag, s)
@@ -242,10 +245,10 @@ class SingleExcitationPlus(Operation):
         decomp_ops = [
             qml.PauliX(wires=wires[0]),
             qml.PauliX(wires=wires[1]),
-            qml.ControlledPhaseShift(theta / 2, wires=[wires[1], wires[0]]),
+            qml.ControlledPhaseShift(-theta / 2, wires=[wires[1], wires[0]]),
             qml.PauliX(wires=wires[0]),
             qml.PauliX(wires=wires[1]),
-            qml.ControlledPhaseShift(theta / 2, wires=[wires[0], wires[1]]),
+            qml.ControlledPhaseShift(-theta / 2, wires=[wires[0], wires[1]]),
             qml.CNOT(wires=[wires[0], wires[1]]),
             qml.CRY(theta, wires=[wires[1], wires[0]]),
             qml.CNOT(wires=[wires[0], wires[1]]),
@@ -254,10 +257,10 @@ class SingleExcitationPlus(Operation):
 
     def adjoint(self):
         (phi,) = self.parameters
-        return SingleExcitationPlus(-phi, wires=self.wires)
+        return SingleExcitationMinus(-phi, wires=self.wires)
 
     def label(self, decimals=None, base_label=None):
-        return super().label(decimals=decimals, base_label=base_label or "G₊")
+        return super().label(decimals=decimals, base_label=base_label or "G₋")
 
 
 class DoubleExcitation(Operation):
@@ -316,6 +319,7 @@ class DoubleExcitation(Operation):
     G[3, 12] = -1j  # 3 (dec) = 0011 (bin)
     G[12, 3] = 1j  # 12 (dec) = 1100 (bin)
     generator = [G, -1 / 2]
+    parameter_frequencies = [(0.5, 1.0)]
 
     @classmethod
     def _matrix(cls, *params):
@@ -410,6 +414,7 @@ class DoubleExcitationPlus(Operation):
     G[3, 12] = -1j  # 3 (dec) = 0011 (bin)
     G[12, 3] = 1j  # 12 (dec) = 1100 (bin)
     generator = [G, -1 / 2]
+    parameter_frequencies = [(1,)]
 
     @classmethod
     def _matrix(cls, *params):
@@ -479,6 +484,7 @@ class DoubleExcitationMinus(Operation):
     G[3, 12] = -1j  # 3 (dec) = 0011 (bin)
     G[12, 3] = 1j  # 12 (dec) = 1100 (bin)
     generator = [G, -1 / 2]
+    parameter_frequencies = [(1,)]
 
     @classmethod
     def _matrix(cls, *params):
@@ -535,8 +541,9 @@ class OrbitalRotation(Operation):
 
     * Number of wires: 4
     * Number of parameters: 1
-    * Gradient recipe: The ``OrbitalRotation`` operator satisfies the four-term parameter-shift rule
-      (see Appendix F, https://arxiv.org/abs/2104.05695)
+    * Gradient recipe: The ``OrbitalRotation`` operator has 4 equidistant frequencies
+      :math:`\{0.5, 1, 1.5, 2\}`, and thus permits an 8-term parameter-shift rule.
+      (see https://arxiv.org/abs/2107.12390).
 
     Args:
         phi (float): rotation angle :math:`\phi`
@@ -563,7 +570,7 @@ class OrbitalRotation(Operation):
     num_wires = 4
     num_params = 1
     grad_method = "A"
-    grad_recipe = four_term_grad_recipe
+
     generator = [
         qml.math.array(
             [
@@ -587,12 +594,19 @@ class OrbitalRotation(Operation):
         ),
         -1 / 2,
     ]
+    parameter_frequencies = [(0.5, 1.0, 1.5, 2.0)]
+
+    @property
+    def grad_recipe(self):
+        coeffs, shifts = qml.gradients.generate_shift_rule(self.parameter_frequencies[0])
+        return [np.stack([coeffs, np.ones_like(coeffs), shifts]).T]
 
     @classmethod
     def _matrix(cls, *params):
-        # This matrix is the "sign flipped" version of that on p18 of https://arxiv.org/abs/2104.05695,
-        # where the sign flip is to adjust for the opposite convention used by authors for naming wires.
-        # Additionally, there was a typo in the sign of a matrix element "s" at [2, 8], which is fixed here.
+        r"""This matrix is the "sign flipped" version of that on p18 of https://arxiv.org/abs/2104.05695,
+        where the sign flip is to adjust for the opposite convention used by authors for naming wires.
+        Additionally, there was a typo in the sign of a matrix element "s" at [2, 8], which is fixed here.
+        """
 
         phi = params[0]
         c = qml.math.cos(phi / 2)


Test patch:
diff --git a/tests/ops/qubit/test_parametric_ops.py b/tests/ops/qubit/test_parametric_ops.py
index e35614a0e..ecc5851d5 100644
--- a/tests/ops/qubit/test_parametric_ops.py
+++ b/tests/ops/qubit/test_parametric_ops.py
@@ -28,23 +28,11 @@ PARAMETRIZED_OPERATIONS = [
     qml.RX(0.123, wires=0),
     qml.RY(1.434, wires=0),
     qml.RZ(2.774, wires=0),
-    qml.S(wires=0),
-    qml.SX(wires=0),
-    qml.T(wires=0),
-    qml.CNOT(wires=[0, 1]),
-    qml.CZ(wires=[0, 1]),
-    qml.CY(wires=[0, 1]),
-    qml.SWAP(wires=[0, 1]),
-    qml.ISWAP(wires=[0, 1]),
-    qml.SISWAP(wires=[0, 1]),
-    qml.SQISW(wires=[0, 1]),
-    qml.CSWAP(wires=[0, 1, 2]),
     qml.PauliRot(0.123, "Y", wires=0),
     qml.IsingXX(0.123, wires=[0, 1]),
     qml.IsingYY(0.123, wires=[0, 1]),
     qml.IsingZZ(0.123, wires=[0, 1]),
     qml.Rot(0.123, 0.456, 0.789, wires=0),
-    qml.Toffoli(wires=[0, 1, 2]),
     qml.PhaseShift(2.133, wires=0),
     qml.ControlledPhaseShift(1.777, wires=[0, 2]),
     qml.CPhase(1.777, wires=[0, 2]),
@@ -55,10 +43,6 @@ PARAMETRIZED_OPERATIONS = [
     qml.U1(0.123, wires=0),
     qml.U2(3.556, 2.134, wires=0),
     qml.U3(2.009, 1.894, 0.7789, wires=0),
-    qml.Hadamard(wires=0),
-    qml.PauliX(wires=0),
-    qml.PauliZ(wires=0),
-    qml.PauliY(wires=0),
     qml.CRot(0.123, 0.456, 0.789, wires=[0, 1]),
     qml.QubitUnitary(np.eye(2) * 1j, wires=0),
     qml.DiagonalQubitUnitary(np.array([1.0, 1.0j]), wires=1),
@@ -71,12 +55,34 @@ PARAMETRIZED_OPERATIONS = [
     qml.DoubleExcitation(0.123, wires=[0, 1, 2, 3]),
     qml.DoubleExcitationPlus(0.123, wires=[0, 1, 2, 3]),
     qml.DoubleExcitationMinus(0.123, wires=[0, 1, 2, 3]),
+]
+
+NON_PARAMETRIZED_OPERATIONS = [
+    qml.S(wires=0),
+    qml.SX(wires=0),
+    qml.T(wires=0),
+    qml.CNOT(wires=[0, 1]),
+    qml.CZ(wires=[0, 1]),
+    qml.CY(wires=[0, 1]),
+    qml.SWAP(wires=[0, 1]),
+    qml.ISWAP(wires=[0, 1]),
+    qml.SISWAP(wires=[0, 1]),
+    qml.SQISW(wires=[0, 1]),
+    qml.CSWAP(wires=[0, 1, 2]),
+    qml.Toffoli(wires=[0, 1, 2]),
+    qml.Hadamard(wires=0),
+    qml.PauliX(wires=0),
+    qml.PauliZ(wires=0),
+    qml.PauliY(wires=0),
+    qml.MultiControlledX(control_wires=[0, 1], wires=2, control_values="01"),
     qml.QubitSum(wires=[0, 1, 2]),
 ]
 
+ALL_OPERATIONS = NON_PARAMETRIZED_OPERATIONS + PARAMETRIZED_OPERATIONS
+
 
 class TestOperations:
-    @pytest.mark.parametrize("op", PARAMETRIZED_OPERATIONS)
+    @pytest.mark.parametrize("op", ALL_OPERATIONS)
     def test_parametrized_op_copy(self, op, tol):
         """Tests that copied parametrized ops function as expected"""
         copied_op = copy.copy(op)
@@ -87,7 +93,7 @@ class TestOperations:
         np.testing.assert_allclose(op.matrix, copied_op2.matrix, atol=tol)
         op.inv()
 
-    @pytest.mark.parametrize("op", PARAMETRIZED_OPERATIONS)
+    @pytest.mark.parametrize("op", ALL_OPERATIONS)
     def test_adjoint_unitaries(self, op, tol):
         op_d = op.adjoint()
         res1 = np.dot(op.matrix, op_d.matrix)
@@ -97,6 +103,32 @@ class TestOperations:
         assert op.wires == op_d.wires
 
 
+class TestParameterFrequencies:
+    @pytest.mark.parametrize("op", PARAMETRIZED_OPERATIONS)
+    def test_parameter_frequencies_match_generator(self, op, tol):
+        if op.generator[0] is None:
+            # For operations without generator, we only can try and check that the function
+            # executes properly, not its return value.
+            try:
+                op.parameter_frequencies
+            except qml.operation.OperatorPropertyUndefined:
+                pytest.skip(f"Operation {op.name} does not have parameter frequencies defined.")
+            pytest.skip(f"Operation {op.name} does not have a generator defined to test against.")
+        gen, coeff = op.generator
+        if isinstance(gen, np.ndarray):
+            matrix = gen
+        elif hasattr(gen, "matrix"):
+            matrix = gen.matrix
+        else:
+            raise ValueError
+
+        gen_eigvals = tuple(np.linalg.eigvalsh(matrix))
+        freqs_from_gen = np.abs(coeff) * np.array(qml.gradients.eigvals_to_frequencies(gen_eigvals))
+
+        freqs = op.parameter_frequencies
+        assert np.allclose(freqs, freqs_from_gen, atol=tol)
+
+
 class TestDecompositions:
     def test_phase_decomposition(self, tol):
         """Tests that the decomposition of the Phase gate is correct"""
diff --git a/tests/ops/qubit/test_qchem_ops.py b/tests/ops/qubit/test_qchem_ops.py
index 68a39a3fa..f335c338d 100644
--- a/tests/ops/qubit/test_qchem_ops.py
+++ b/tests/ops/qubit/test_qchem_ops.py
@@ -36,6 +36,38 @@ from gate_data import (
 )
 
 
+PARAMETRIZED_QCHEM_OPERATIONS = [
+    qml.SingleExcitation(0.14, wires=[0, 1]),
+    qml.SingleExcitationMinus(0.14, wires=[0, 1]),
+    qml.SingleExcitationPlus(0.14, wires=[0, 1]),
+    qml.DoubleExcitation(0.14, wires=[0, 1, 2, 3]),
+    qml.DoubleExcitationMinus(0.14, wires=[0, 1, 2, 3]),
+    qml.DoubleExcitationPlus(0.14, wires=[0, 1, 2, 3]),
+    qml.OrbitalRotation(0.14, wires=[0, 1, 2, 3]),
+]
+
+
+class TestParameterFrequencies:
+    @pytest.mark.parametrize("op", PARAMETRIZED_QCHEM_OPERATIONS)
+    def test_parameter_frequencies_match_generator(self, op, tol):
+        gen, coeff = op.generator
+        if isinstance(gen, np.ndarray):
+            matrix = gen
+        elif hasattr(gen, "matrix"):
+            matrix = gen.matrix
+        else:
+            raise ValueError
+
+        gen_eigvals = tuple(np.linalg.eigvalsh(coeff * matrix))
+        freqs_from_gen = np.array(
+            sorted(set(np.round(qml.gradients.eigvals_to_frequencies(gen_eigvals), 9)))
+        )
+        freqs_from_gen = freqs_from_gen[freqs_from_gen > 0]
+
+        freqs = op.parameter_frequencies
+        assert np.allclose(freqs, freqs_from_gen, atol=tol)
+
+
 class TestDecomposition:
     @pytest.mark.parametrize("phi", [-0.1, 0.2, 0.5])
     def test_single_excitation_plus_decomp(self, phi):
@@ -567,6 +599,21 @@ class TestDoubleExcitation:
 class TestOrbitalRotation:
     """Test OrbitalRotation gate operation"""
 
+    def grad_circuit_0(self, phi):
+        qml.PauliX(1)
+        qml.Hadamard(2)
+        qml.OrbitalRotation(phi, wires=[0, 1, 2, 3])
+        return qml.expval(qml.PauliX(0) @ qml.PauliX(1) @ qml.PauliX(3))
+
+    def grad_circuit_1(self, phi):
+        qml.PauliX(0)
+        qml.PauliX(1)
+        qml.OrbitalRotation(phi, wires=[0, 1, 2, 3])
+        return qml.expval(qml.PauliZ(0) @ qml.PauliZ(1) @ qml.PauliZ(3))
+
+    def expected_grad_fn(self, phi):
+        return -0.55 * np.sin(3 * phi / 2) * 3 / 2 - 0.7 * np.sin(phi) + 0.55 / 2 * np.sin(phi / 2)
+
     @pytest.mark.parametrize("phi", [-0.1, 0.2, np.pi / 4])
     def test_orbital_rotation_matrix(self, phi):
         """Tests that the OrbitalRotation operation calculates the correct matrix"""
@@ -812,7 +859,7 @@ class TestOrbitalRotation:
         "phi",
         [
             pnp.array(-0.1, requires_grad=True),
-            pnp.array(0.1, requires_grad=True),
+            pnp.array(0.1421, requires_grad=True),
         ],
     )
     def test_autograd_grad(self, phi, diff_method):
@@ -823,20 +870,20 @@ class TestOrbitalRotation:
 
         dev = qml.device("default.qubit.autograd", wires=4)
 
-        @qml.qnode(dev, interface="autograd", diff_method=diff_method)
-        def circuit(phi):
-            qml.PauliX(wires=0)
-            qml.PauliX(wires=1)
-            qml.OrbitalRotation(phi, wires=[0, 1, 2, 3])
+        circuit_0 = qml.QNode(
+            self.grad_circuit_0, dev, interface="autograd", diff_method=diff_method
+        )
+        circuit_1 = qml.QNode(
+            self.grad_circuit_1, dev, interface="autograd", diff_method=diff_method
+        )
+        total = lambda phi: 1.1 * circuit_0(phi) + 0.7 * circuit_1(phi)
 
-            return qml.expval(qml.PauliZ(0))
-
-        assert np.allclose(qml.grad(circuit)(phi), np.sin(phi))
+        assert np.allclose(qml.grad(total)(phi), self.expected_grad_fn(phi))
 
     @pytest.mark.parametrize("diff_method", ["parameter-shift", "backprop"])
     @pytest.mark.parametrize(
         ("phi"),
-        [(-0.1), (0.1)],
+        [-0.1, 0.1421],
     )
     def test_tf_grad(self, phi, diff_method):
         """Tests that gradients are computed correctly using the
@@ -845,28 +892,22 @@ class TestOrbitalRotation:
         tf = pytest.importorskip("tensorflow")
         dev = qml.device("default.qubit.tf", wires=4)
 
-        @qml.qnode(dev, interface="tf", diff_method=diff_method)
-        def circuit(phi):
-            qml.PauliX(wires=0)
-            qml.PauliX(wires=1)
-            qml.OrbitalRotation(phi, wires=[0, 1, 2, 3])
-            return qml.expval(qml.PauliZ(0))
+        circuit_0 = qml.QNode(self.grad_circuit_0, dev, interface="tf", diff_method=diff_method)
+        circuit_1 = qml.QNode(self.grad_circuit_1, dev, interface="tf", diff_method=diff_method)
+        total = lambda phi: 1.1 * circuit_0(phi) + 0.7 * circuit_1(phi)
 
         phi_t = tf.Variable(phi, dtype=tf.float64)
         with tf.GradientTape() as tape:
-            res = circuit(phi_t)
+            res = total(phi_t)
 
         grad = tape.gradient(res, phi_t)
 
-        assert np.allclose(grad, np.sin(phi))
+        assert np.allclose(grad, self.expected_grad_fn(phi))
 
     @pytest.mark.parametrize("diff_method", ["parameter-shift", "backprop"])
     @pytest.mark.parametrize(
         ("phi"),
-        [
-            (-0.1),
-            (0.1),
-        ],
+        [-0.1, 0.1421],
     )
     def test_jax_grad(self, phi, diff_method):
         """Tests that gradients and operations are computed correctly using the
@@ -879,22 +920,18 @@ class TestOrbitalRotation:
 
         dev = qml.device("default.qubit.jax", wires=4)
 
-        @qml.qnode(dev, interface="jax", diff_method=diff_method)
-        def circuit(phi):
-            qml.PauliX(wires=0)
-            qml.PauliX(wires=1)
-            qml.OrbitalRotation(phi, wires=[0, 1, 2, 3])
-            return qml.expval(qml.PauliZ(0))
+        circuit_0 = qml.QNode(self.grad_circuit_0, dev, interface="jax", diff_method=diff_method)
+        circuit_1 = qml.QNode(self.grad_circuit_1, dev, interface="jax", diff_method=diff_method)
+        total = lambda phi: 1.1 * circuit_0(phi) + 0.7 * circuit_1(phi)
 
-        assert np.allclose(jax.grad(circuit)(phi), np.sin(phi))
+        phi_j = jax.numpy.array(phi)
+
+        assert np.allclose(jax.grad(total)(phi_j), self.expected_grad_fn(phi))
 
     @pytest.mark.parametrize("diff_method", ["parameter-shift", "backprop"])
     @pytest.mark.parametrize(
         ("phi"),
-        [
-            (-0.1),
-            (0.1),
-        ],
+        [-0.1, 0.1421],
     )
     def test_torch_grad(self, phi, diff_method):
         """Tests that gradients and operations are computed correctly using the
@@ -904,19 +941,15 @@ class TestOrbitalRotation:
 
         dev = qml.device("default.qubit.torch", wires=4)
 
-        @qml.qnode(dev, interface="torch", diff_method=diff_method)
-        def circuit(phi):
-            qml.PauliX(wires=0)
-            qml.PauliX(wires=1)
-            qml.OrbitalRotation(phi, wires=[0, 1, 2, 3])
-            return qml.expval(qml.PauliZ(0))
+        circuit_0 = qml.QNode(self.grad_circuit_0, dev, interface="torch", diff_method=diff_method)
+        circuit_1 = qml.QNode(self.grad_circuit_1, dev, interface="torch", diff_method=diff_method)
+        total = lambda phi: 1.1 * circuit_0(phi) + 0.7 * circuit_1(phi)
 
         phi_t = torch.tensor(phi, dtype=torch.complex128, requires_grad=True)
-
-        result = circuit(phi_t)
+        result = total(phi_t)
         result.backward()
 
-        assert np.allclose(phi_t.grad, np.sin(phi))
+        assert np.allclose(phi_t.grad, self.expected_grad_fn(phi))
 
 
 label_data = [
diff --git a/tests/test_operation.py b/tests/test_operation.py
index e94e11888..a72ee228a 100644
--- a/tests/test_operation.py
+++ b/tests/test_operation.py
@@ -153,6 +153,59 @@ class TestOperationConstruction:
         op = DummyOp(x, wires=0)
         assert op.grad_recipe == ([[1.0, 1.0, x], [1.0, 0.0, -x]],)
 
+    def test_frequencies_default_single_param(self):
+        """Test that an operation with default parameter frequencies
+        and a single parameter works correctly."""
+
+        class DummyOp(qml.operation.Operation):
+            r"""Dummy custom operation"""
+            num_wires = 1
+            grad_method = "A"
+            generator = [qml.PauliX, -0.2]
+
+        x = 0.654
+        op = DummyOp(x, wires=0)
+        assert op.parameter_frequencies == (0.4,)
+
+    def test_frequencies_default_multi_param(self):
+        """Test that an operation with default parameter frequencies and multiple
+        parameters raises an error when calling parameter_frequencies."""
+
+        class DummyOp(qml.operation.Operation):
+            r"""Dummy custom operation"""
+            num_params = 3
+            num_wires = 1
+            grad_method = "A"
+
+        x = [0.654, 2.31, 0.1]
+        op = DummyOp(*x, wires=0)
+        with pytest.raises(
+            qml.operation.OperatorPropertyUndefined, match="does not have parameter"
+        ):
+            op.parameter_frequencies
+
+    @pytest.mark.parametrize("num_param", [1, 2])
+    def test_frequencies_parameter_dependent(self, num_param):
+        """Test that an operation with parameter frequencies that depend on
+        its instantiated parameter values works correctly"""
+
+        class DummyOp(qml.operation.Operation):
+            r"""Dummy custom operation"""
+            num_params = num_param
+            num_wires = 1
+            grad_method = "A"
+
+            @property
+            def parameter_frequencies(self):
+                x = self.data
+                return [(0.2, _x) for _x in x]
+
+        x = [0.654, 2.31][:num_param]
+        op = DummyOp(*x, wires=0)
+        f = op.parameter_frequencies
+        for i in range(num_param):
+            assert f[i] == (0.2, x[i])
+
     def test_no_wires_passed(self):
         """Test exception raised if no wires are passed"""