Qutrit Rotations #2 (TRY)

doc/releases/changelog-dev.md

@@ -12,6 +12,10 @@
   given subspace.
   [(#2845)](https://github.com/PennyLaneAI/pennylane/pull/2845)
 
+* Added the `TRY` qutrit rotation operator, which allows applying a Y rotation on a
+  given subspace.
+  [(#2846)](https://github.com/PennyLaneAI/pennylane/pull/2846)
+
 <h3>Improvements 🛠</h3>
 
 * `DiagonalQubitUnitary` now decomposes into `RZ`, `IsingZZ` and `MultiRZ` gates

pennylane/devices/default_qutrit.py

@@ -62,6 +62,7 @@ class DefaultQutrit(QutritDevice):
         "TSWAP",
         "THadamard",
         "TRX",
+        "TRY",
     }
 
     # Identity is supported as an observable for qml.state() to work correctly. However, any

pennylane/ops/qutrit/__init__.py

@@ -37,6 +37,7 @@
     "TSWAP",
     "THadamard",
     "TRX",
+    "TRY",
 }
 
 obs = {

pennylane/ops/qutrit/parametric_ops.py

@@ -99,13 +99,13 @@ def __init__(self, phi, wires, subspace=(0, 1), do_queue=True, id=None):
 
     @property
     def subspace(self):
-        """The single-qutrit basis states the operator acts on.
+        """The single-qutrit basis states which the operator acts on
 
-        The subspace defines which two basis states the opreation acts on. The basis state
-        not included in the subspace remains unaffected.
+        This subspace determines which two single-qutrit basis states the operator acts on.
+        The remaining basis state is not affected by the operator.
 
         Returns:
-            tuple[int]: subspace on which the operator acts
+            tuple[int]: subspace on which operator acts
         """
         return self._subspace
 
@@ -169,3 +169,147 @@ def adjoint(self):
 
     def pow(self, z):
         return [TRX(self.data[0] * z, wires=self.wires, subspace=self.subspace)]
+
+
+class TRY(Operation):
+    r"""
+    The single qutrit Y rotation
+
+    Performs the RY operation on the specified 2D subspace. The subspace is
+    given as a keyword argument and determines which two of three single-qutrit
+    basis states the operation applies to.
+
+    The construction of this operator is based on section 3 of
+    `Di et al. (2012) <https://arxiv.org/abs/1105.5485>`_.
+
+    .. math:: TR_y^{jk}(\phi) = \exp(-i\phi\sigma_y^{jk}/2),
+                \sigma_y^{jk} = -i |j\rangle\langle k| + i |k\rangle\langle j|,
+                j, k \in \{0, 1, 2\}, j \neq k
+
+    **Details:**
+
+    * Number of wires: 1
+    * Number of parameters: 1
+    * Number of dimensions per parameter: (0,)
+
+    Args:
+        phi (float): rotation angle :math:`\phi`
+        wires (Sequence[int] or int): the wire the operation acts on
+        subspace (Sequence[int]): the 2D subspace on which to apply operation
+        do_queue (bool): Indicates whether the operator should be
+            immediately pushed into the Operator queue (optional)
+        id (str or None): String representing the operation (optional)
+
+    **Example**
+
+    The specified subspace will determine which basis states the operation actually
+    applies to:
+
+    >>> qml.TRY(0.5, wires=0, subspace=[0, 1]).matrix()
+    array([[ 0.96891242+0.j, -0.24740396-0.j,  0.        +0.j],
+           [ 0.24740396+0.j,  0.96891242+0.j,  0.        +0.j],
+           [ 0.        +0.j,  0.        +0.j,  1.        +0.j]])
+
+    >>> qml.TRY(0.5, wires=0, subspace=[0, 2]).matrix()
+    array([[ 0.96891242+0.j,  0.        +0.j, -0.24740396-0.j],
+           [ 0.        +0.j,  1.        +0.j,  0.        +0.j],
+           [ 0.24740396+0.j,  0.        +0.j,  0.96891242+0.j]])
+
+    >>> qml.TRY(0.5, wires=0, subspace=[1, 2]).matrix()
+    array([[ 1.        +0.j,  0.        +0.j,  0.        +0.j],
+           [ 0.        +0.j,  0.96891242+0.j, -0.24740396-0.j],
+           [ 0.        +0.j,  0.24740396+0.j,  0.96891242+0.j]])
+    """
+    num_wires = 1
+    num_params = 1
+    """int: Number of trainable parameters that the operator depends on."""
+
+    ndim_params = (0,)
+    """tuple[int]: Number of dimensions per trainable parameter that the operator depends on."""
+
+    grad_method = "A"
+    parameter_frequencies = [(0.5, 1)]
+
+    # Internal dictionary to map subpsaces to Gell-Mann observable for the generator
+    _index_dict = {(0, 1): 2, (0, 2): 5, (1, 2): 7}
+
+    def generator(self):
+        return qml.s_prod(-0.5, qml.GellMann(self.wires, index=self._index_dict[self.subspace]))
+
+    def __init__(self, phi, wires, subspace=(0, 1), do_queue=True, id=None):
+        self._subspace = self.validate_subspace(subspace)
+        self._hyperparameters = {
+            "subspace": self._subspace,
+        }
+        super().__init__(phi, wires=wires, do_queue=do_queue, id=id)
+
+    @property
+    def subspace(self):
+        """The single-qutrit basis states which the operator acts on
+
+        This subspace determines which two single-qutrit basis states the operator acts on.
+        The remaining basis state is not affected by the operator.
+
+        Returns:
+            tuple[int]: subspace on which operator acts
+        """
+        return self._subspace
+
+    @staticmethod
+    def compute_matrix(theta, subspace=(0, 1)):  # pylint: disable=arguments-differ
+        r"""Representation of the operator as a canonical matrix in the computational basis (static method).
+
+        The canonical matrix is the textbook matrix representation that does not consider wires.
+        Implicitly, this assumes that the wires of the operator correspond to the global wire order.
+
+        .. seealso:: :meth:`~.TRY.matrix`
+
+        Args:
+            theta (tensor_like or float): rotation angle
+            subspace (Sequence[int]): the 2D subspace on which to apply operation
+
+        Returns:
+            tensor_like: canonical matrix
+
+        **Example**
+
+        >>> qml.TRY.compute_matrix(torch.tensor(0.5), subspace=[0, 2])
+        tensor([[ 0.9689+0.j,  0.0000+0.j, -0.2474-0.j],
+                [ 0.0000+0.j,  1.0000+0.j,  0.0000+0.j],
+                [ 0.2474+0.j,  0.0000+0.j,  0.9689+0.j]])
+        """
+        c = qml.math.cos(theta / 2)
+        s = qml.math.sin(theta / 2)
+
+        if qml.math.get_interface(theta) == "tensorflow":
+            c = qml.math.cast_like(c, 1j)
+            s = qml.math.cast_like(s, 1j)
+
+        # The following avoids casting an imaginary quantity to reals when backpropagating
+        c = (1 + 0j) * c
+        s = (1 + 0j) * s
+
+        shape = qml.math.shape(theta)
+        is_broadcasted = len(shape) != 0 and shape[0] > 1
+        # Construct identity matrices and cast to complex type
+        mat = (
+            qml.math.tensordot([1] * qml.math.shape(theta)[0], qml.math.eye(3), axes=0)
+            if is_broadcasted
+            else qml.math.eye(3)
+        )
+        mat = qml.math.cast_like(mat, s)
+        slices = tuple(itertools.product(subspace, subspace))
+        if is_broadcasted:
+            slices = [(Ellipsis, *s) for s in slices]
+
+        # Put rotation terms in the appropriate indices using the slices
+        mat[slices[0]] = mat[slices[3]] = c
+        mat[slices[1]] = -s
+        mat[slices[2]] = s
+        return qml.math.convert_like(mat, theta)
+
+    def adjoint(self):
+        return TRY(-self.data[0], wires=self.wires, subspace=self.subspace)
+
+    def pow(self, z):
+        return [TRY(self.data[0] * z, wires=self.wires, subspace=self.subspace)]