Repository: PennyLaneAI/pennylane
Pull request #5199 (2024-02-13 21:24:42)

Make non-parametric controlled ops inherit from Controlled
**Context:**
This is the second part of a rework for controlled operations (more to come). This PR simply moves the rest of non-parametric controlled operations from `ops/qubit` to `ops/op_math`, and makes them inherit from `ControlledOp`

**Description of the Change:**
1. `CCZ`, `CSWAP`, `CH`, `CNOT`, `Toffoli` now inherit from `ControlledOp`.
2. Removed some duplicate code in `controlled_ops` involving `matrix` and `decomposition`

**Benefits:**
Cleaner code, more consistent behaviour

**Related GitHub Issues:**
https://github.com/PennyLaneAI/pennylane/issues/1447

**Related Shortcut Story:**
[sc-37951]

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://app.codecov.io/gh/PennyLaneAI/pennylane/pull/5199?src=pr&el=h1&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI) Report
Attention: `1 lines` in your changes are missing coverage. Please review.
> Comparison is base [(`256b1cc`)](https://app.codecov.io/gh/PennyLaneAI/pennylane/commit/256b1ccefe2420d9951ed276b16b230502d6f33a?el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI) 99.68% compared to head [(`eec3b0b`)](https://app.codecov.io/gh/PennyLaneAI/pennylane/pull/5199?src=pr&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI) 99.67%.
> Report is 1 commits behind head on master.

| [Files](https://app.codecov.io/gh/PennyLaneAI/pennylane/pull/5199?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI) | Patch % | Lines |
|---|---|---|
| [.../transforms/decompositions/clifford\_t\_transform.py](https://app.codecov.io/gh/PennyLaneAI/pennylane/pull/5199?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI#diff-cGVubnlsYW5lL3RyYW5zZm9ybXMvZGVjb21wb3NpdGlvbnMvY2xpZmZvcmRfdF90cmFuc2Zvcm0ucHk=) | 75.00% | [1 Missing :warning: ](https://app.codecov.io/gh/PennyLaneAI/pennylane/pull/5199?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI) |

<details><summary>Additional details and impacted files</summary>


```diff
@@            Coverage Diff             @@
##           master    #5199      +/-   ##
==========================================
- Coverage   99.68%   99.67%   -0.02%     
==========================================
  Files         396      396              
  Lines       36640    36350     -290     
==========================================
- Hits        36526    36233     -293     
- Misses        114      117       +3     
```



</details>

[:umbrella: View full report in Codecov by Sentry](https://app.codecov.io/gh/PennyLaneAI/pennylane/pull/5199?src=pr&el=continue&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI).   
:loudspeaker: Have feedback on the report? [Share it here](https://about.codecov.io/codecov-pr-comment-feedback/?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI).


Patch:
diff --git a/doc/introduction/compiling_circuits.rst b/doc/introduction/compiling_circuits.rst
index 26fee8008..4deed05b0 100644
--- a/doc/introduction/compiling_circuits.rst
+++ b/doc/introduction/compiling_circuits.rst
@@ -236,7 +236,7 @@ We define the custom decompositions like so, and pass them to a device:
 
 .. code-block:: python
 
-    def custom_cnot(wires):
+    def custom_cnot(wires, **_):
         return [
             qml.Hadamard(wires=wires[1]),
             qml.CZ(wires=[wires[0], wires[1]]),
@@ -248,6 +248,8 @@ We define the custom decompositions like so, and pass them to a device:
     decomp_dev = qml.device("default.qubit", wires=3, custom_decomps=custom_decomps)
     decomp_qnode = qml.QNode(circuit, decomp_dev)
 
+Note that custom decomposition functions should accept keyword arguments even when it is not used.
+
 Now when we draw or run a QNode on this device, the gates will be expanded
 according to our specifications:
 
diff --git a/doc/releases/changelog-dev.md b/doc/releases/changelog-dev.md
index 4374d5079..b5370480d 100644
--- a/doc/releases/changelog-dev.md
+++ b/doc/releases/changelog-dev.md
@@ -338,8 +338,9 @@
   at different stages.
   [(#5084)](https://github.com/PennyLaneAI/pennylane/pull/5084)
 
-* `CRX`, `CRY`, `CRZ`, `CROT`, and `ControlledPhaseShift` (i.e. `CPhaseShift`) now inherit from `ControlledOp`, giving them additional properties such as `control_wire` and `control_values`. Calling `qml.ctrl` on `RX`, `RY`, `RZ`, `Rot`, and `PhaseShift` with a single control wire will return gates of types `CRX`, `CRY`, etc. as opposed to a general `Controlled` operator.
+* All custom controlled operations such as `CRX`, `CZ`, `CNOT`, `ControlledPhaseShift` now inherit from `ControlledOp`, giving them additional properties such as `control_wire` and `control_values`. Calling `qml.ctrl` on `RX`, `RY`, `RZ`, `Rot`, and `PhaseShift` with a single control wire will return gates of types `CRX`, `CRY`, etc. as opposed to a general `Controlled` operator.
   [(#5069)](https://github.com/PennyLaneAI/pennylane/pull/5069)
+  [(#5199)](https://github.com/PennyLaneAI/pennylane/pull/5199)
 
 * CI will now fail if coverage data fails to upload to codecov. Previously, it would silently pass
   and the codecov check itself would never execute.
diff --git a/pennylane/__init__.py b/pennylane/__init__.py
index 8a346ac0a..79a28c0d3 100644
--- a/pennylane/__init__.py
+++ b/pennylane/__init__.py
@@ -295,7 +295,7 @@ def device(name, *args, **kwargs):
 
     .. code-block:: python
 
-        def ion_trap_cnot(wires):
+        def ion_trap_cnot(wires, **_):
             return [
                 qml.RY(np.pi/2, wires=wires[0]),
                 qml.IsingXX(np.pi/2, wires=wires),
@@ -307,6 +307,8 @@ def device(name, *args, **kwargs):
     Next, we create a device, and a QNode for testing. When constructing the
     QNode, we can set the expansion strategy to ``"device"`` to ensure the
     decomposition is applied and will be viewable when we draw the circuit.
+    Note that custom decompositions should accept keyword arguments even when
+    it is not used.
 
     .. code-block:: python
 
diff --git a/pennylane/drawer/utils.py b/pennylane/drawer/utils.py
index 90c05ec08..3a18ddea6 100644
--- a/pennylane/drawer/utils.py
+++ b/pennylane/drawer/utils.py
@@ -14,7 +14,6 @@
 """
 This module contains some useful utility functions for circuit drawing.
 """
-import pennylane as qml
 from pennylane.ops import Controlled, Conditional
 from pennylane.measurements import MeasurementProcess, MidMeasureMP, MeasurementValue
 
@@ -96,12 +95,6 @@ def convert_wire_order(tape, wire_order=None, show_all_wires=False):
     return {wire: ind for ind, wire in enumerate(wire_order)}
 
 
-def _is_controlled(op):
-    """Checks whether an operation is a controlled operation."""
-    # TODO: remove the other ones once they all inherit from Controlled [sc-37951]
-    return isinstance(op, (Controlled, qml.CCZ, qml.CSWAP, qml.CH, qml.CNOT, qml.Toffoli))
-
-
 def unwrap_controls(op):
     """Unwraps nested controlled operations for drawing.
 
@@ -129,7 +122,7 @@ def unwrap_controls(op):
 
         while hasattr(next_ctrl, "base"):
 
-            if _is_controlled(next_ctrl.base):
+            if isinstance(next_ctrl.base, Controlled):
                 base_control_wires = getattr(next_ctrl.base, "control_wires", [])
                 control_wires += base_control_wires
 
diff --git a/pennylane/ops/functions/bind_new_parameters.py b/pennylane/ops/functions/bind_new_parameters.py
index dd8bff6f5..482548951 100644
--- a/pennylane/ops/functions/bind_new_parameters.py
+++ b/pennylane/ops/functions/bind_new_parameters.py
@@ -113,10 +113,13 @@ def bind_new_parameters_composite_op(op: CompositeOp, params: Sequence[TensorLik
 
 @bind_new_parameters.register(qml.CY)
 @bind_new_parameters.register(qml.CZ)
+@bind_new_parameters.register(qml.CH)
+@bind_new_parameters.register(qml.CCZ)
+@bind_new_parameters.register(qml.CSWAP)
+@bind_new_parameters.register(qml.CNOT)
+@bind_new_parameters.register(qml.Toffoli)
 @bind_new_parameters.register(qml.MultiControlledX)
-def bind_new_parameters_copy(
-    op: Union[qml.CY, qml.CZ, qml.MultiControlledX], params: Sequence[TensorLike]
-):  # pylint:disable=unused-argument
+def bind_new_parameters_copy(op, params: Sequence[TensorLike]):  # pylint:disable=unused-argument
     return copy.copy(op)
 
 
diff --git a/pennylane/ops/op_math/__init__.py b/pennylane/ops/op_math/__init__.py
index 20089d30d..cc0b050ed 100644
--- a/pennylane/ops/op_math/__init__.py
+++ b/pennylane/ops/op_math/__init__.py
@@ -64,6 +64,12 @@ Controlled Operator Classes
     ~ControlledQubitUnitary
     ~CY
     ~CZ
+    ~CH
+    ~CCZ
+    ~CSWAP
+    ~CNOT
+    ~Toffoli
+    ~MultiControlledX
     ~CRX
     ~CRY
     ~CRZ
@@ -109,6 +115,11 @@ from .controlled_ops import (
     CRZ,
     CY,
     CZ,
+    CH,
+    CCZ,
+    CSWAP,
+    CNOT,
+    Toffoli,
     MultiControlledX,
 )
 from .decompositions import one_qubit_decomposition, two_qubit_decomposition, sk_decomposition
@@ -125,6 +136,11 @@ controlled_qubit_ops = {
     "ControlledQubitUnitary",
     "CY",
     "CZ",
+    "CH",
+    "CCZ",
+    "CSWAP",
+    "CNOT",
+    "Toffoli",
     "MultiControlledX",
     "CRX",
     "CRY",
diff --git a/pennylane/ops/op_math/controlled.py b/pennylane/ops/op_math/controlled.py
index f9ee20dee..ccb9c3b67 100644
--- a/pennylane/ops/op_math/controlled.py
+++ b/pennylane/ops/op_math/controlled.py
@@ -267,16 +267,6 @@ def _handle_pauli_x_based_controlled_ops(op, control, control_values, work_wires
             wires=control + op.wires, control_values=control_values, work_wires=work_wires
         )
 
-    # TODO: remove special handling of CNOT and Toffoli when they inherit from Controlled
-    if isinstance(op, qml.CNOT):
-        return qml.MultiControlledX(
-            wires=control + op.wires, control_values=control_values + [1], work_wires=work_wires
-        )
-    if isinstance(op, qml.Toffoli):
-        return qml.MultiControlledX(
-            wires=control + op.wires, control_values=control_values + [1, 1], work_wires=work_wires
-        )
-
     work_wires = work_wires or []
     return qml.MultiControlledX(
         wires=control + op.wires,
@@ -542,7 +532,7 @@ class Controlled(SymbolicOp):
     def label(self, decimals=None, base_label=None, cache=None):
         return self.base.label(decimals=decimals, base_label=base_label, cache=cache)
 
-    def matrix(self, wire_order=None):
+    def _compute_matrix_from_base(self):
         base_matrix = self.base.matrix()
         interface = qmlmath.get_interface(base_matrix)
 
@@ -562,17 +552,19 @@ class Controlled(SymbolicOp):
 
         shape = qml.math.shape(base_matrix)
         if len(shape) == 3:  # stack if batching
-            canonical_matrix = qml.math.stack(
+            return qml.math.stack(
                 [qml.math.block_diag([left_pad, _U, right_pad]) for _U in base_matrix]
             )
-        else:
-            canonical_matrix = qmlmath.block_diag([left_pad, base_matrix, right_pad])
 
-        if wire_order is None or self.wires == Wires(wire_order):
-            return qml.math.expand_matrix(
-                canonical_matrix, wires=self.active_wires, wire_order=self.wires
-            )
+        return qmlmath.block_diag([left_pad, base_matrix, right_pad])
 
+    def matrix(self, wire_order=None):
+        if self.compute_matrix is not Operator.compute_matrix:
+            canonical_matrix = self.compute_matrix(*self.data)
+        else:
+            canonical_matrix = self._compute_matrix_from_base()
+
+        wire_order = wire_order or self.wires
         return qml.math.expand_matrix(
             canonical_matrix, wires=self.active_wires, wire_order=wire_order
         )
@@ -622,12 +614,12 @@ class Controlled(SymbolicOp):
 
     @property
     def has_decomposition(self):
+        if self.compute_decomposition is not Operator.compute_decomposition:
+            return True
         if not all(self.control_values):
             return True
         if len(self.control_wires) == 1 and hasattr(self.base, "_controlled"):
             return True
-        if isinstance(self.base, _get_pauli_x_based_ops()):
-            return True
         if _is_single_qubit_special_unitary(self.base):
             return True
         if self.base.has_decomposition:
@@ -636,6 +628,10 @@ class Controlled(SymbolicOp):
         return False
 
     def decomposition(self):
+
+        if self.compute_decomposition is not Operator.compute_decomposition:
+            return self.compute_decomposition(*self.data, self.wires)
+
         if all(self.control_values):
             decomp = _decompose_no_control_values(self)
             if decomp is None:
diff --git a/pennylane/ops/op_math/controlled_ops.py b/pennylane/ops/op_math/controlled_ops.py
index 5601351bb..df862bc42 100644
--- a/pennylane/ops/op_math/controlled_ops.py
+++ b/pennylane/ops/op_math/controlled_ops.py
@@ -23,7 +23,10 @@ import numpy as np
 from scipy.linalg import block_diag
 
 import pennylane as qml
-from pennylane.operation import AnyWires, Wires
+from pennylane.operation import (
+    AnyWires,
+    Wires,
+)
 from pennylane.ops.qubit.matrix_ops import QubitUnitary
 from pennylane.ops.qubit.parametric_ops_single_qubit import stack_last
 from .controlled import ControlledOp
@@ -169,6 +172,111 @@ class ControlledQubitUnitary(ControlledOp):
         return True
 
 
+class CH(ControlledOp):
+    r"""CH(wires)
+    The controlled-Hadamard operator
+
+    .. math:: CH = \begin{bmatrix}
+            1 & 0 & 0 & 0 \\
+            0 & 1 & 0 & 0 \\
+            0 & 0 & \frac{1}{\sqrt{2}} & \frac{1}{\sqrt{2}} \\
+            0 & 0 & \frac{1}{\sqrt{2}} & -\frac{1}{\sqrt{2}}
+        \end{bmatrix}.
+
+    .. note:: The first wire provided corresponds to the **control qubit**.
+
+    **Details:**
+
+    * Number of wires: 2
+    * Number of parameters: 0
+
+    Args:
+        wires (Sequence[int]): the wires the operation acts on
+    """
+
+    num_wires = 2
+    """int: Number of wires that the operation acts on."""
+
+    num_params = 0
+    """int: Number of trainable parameters that the operator depends on."""
+
+    ndim_params = ()
+    """tuple[int]: Number of dimensions per trainable parameter that the operator depends on."""
+
+    name = "CH"
+
+    def _flatten(self):
+        return tuple(), (self.wires,)
+
+    @classmethod
+    def _unflatten(cls, data, metadata):
+        return cls(metadata[0])
+
+    def __init__(self, wires, id=None):
+        control_wires = wires[:1]
+        target_wires = wires[1:]
+        super().__init__(qml.Hadamard(wires=target_wires), control_wires, id=id)
+
+    @staticmethod
+    @lru_cache()
+    def compute_matrix():  # 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:`~.CH.matrix`
+
+
+        Returns:
+            ndarray: matrix
+
+        **Example**
+
+        >>> print(qml.CH.compute_matrix())
+        [[ 1.          0.          0.          0.        ]
+         [ 0.          1.          0.          0.        ]
+         [ 0.          0.          0.70710678  0.70710678]
+         [ 0.          0.          0.70710678 -0.70710678]]
+        """
+        return np.array(
+            [
+                [1, 0, 0, 0],
+                [0, 1, 0, 0],
+                [0, 0, INV_SQRT2, INV_SQRT2],
+                [0, 0, INV_SQRT2, -INV_SQRT2],
+            ]
+        )
+
+    @staticmethod
+    def compute_decomposition(wires):  # pylint: disable=arguments-differ
+        r"""Representation of the operator as a product of other operators (static method).
+
+
+        .. math:: O = O_1 O_2 \dots O_n.
+
+
+        .. seealso:: :meth:`~.CH.decomposition`.
+
+        Args:
+            wires (Iterable, Wires): wires that the operator acts on
+
+        Returns:
+            list[Operator]: decomposition into lower level operations
+
+        **Example:**
+
+        >>> print(qml.CH.compute_decomposition([0, 1]))
+        [qml.RY(-0.78539816339, wires=[1]), CZ(wires=[0, 1]), qml.RY(0.78539816339, wires=[1])]
+
+        """
+        return [
+            qml.RY(-np.pi / 4, wires=wires[1]),
+            qml.CZ(wires=wires),
+            qml.RY(+np.pi / 4, wires=wires[1]),
+        ]
+
+
 class CY(ControlledOp):
     r"""CY(wires)
     The controlled-Y operator
@@ -270,9 +378,6 @@ class CY(ControlledOp):
         """
         return [qml.CRY(np.pi, wires=wires), qml.S(wires=wires[0])]
 
-    def decomposition(self):
-        return self.compute_decomposition(self.wires)
-
 
 class CZ(ControlledOp):
     r"""CZ(wires)
@@ -314,9 +419,9 @@ class CZ(ControlledOp):
     def _unflatten(cls, data, metadata):
         return cls(metadata[0])
 
-    def __init__(self, wires):
+    def __init__(self, wires, id=None):
         control_wire, wire = wires
-        super().__init__(qml.PauliZ(wires=wire), control_wire)
+        super().__init__(qml.PauliZ(wires=wire), control_wire, id=id)
 
     @staticmethod
     @lru_cache()
@@ -348,8 +453,489 @@ class CZ(ControlledOp):
     def compute_decomposition(wires):  # pylint: disable=arguments-differ
         return [qml.ControlledPhaseShift(np.pi, wires=wires)]
 
-    def decomposition(self):
-        return self.compute_decomposition(self.wires)
+
+class CSWAP(ControlledOp):
+    r"""CSWAP(wires)
+    The controlled-swap operator
+
+    .. math:: CSWAP = \begin{bmatrix}
+            1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\
+            0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 \\
+            0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 \\
+            0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 \\
+            0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 \\
+            0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 \\
+            0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 \\
+            0 & 0 & 0 & 0 & 0 & 0 & 0 & 1
+        \end{bmatrix}.
+
+    .. note:: The first wire provided corresponds to the **control qubit**.
+
+    **Details:**
+
+    * Number of wires: 3
+    * Number of parameters: 0
+
+    Args:
+        wires (Sequence[int]): the wires the operation acts on
+    """
+
+    num_wires = 3
+    """int : Number of wires that the operation acts on."""
+
+    num_params = 0
+    """int: Number of trainable parameters that the operator depends on."""
+
+    ndim_params = ()
+    """tuple[int]: Number of dimensions per trainable parameter that the operator depends on."""
+
+    name = "CSWAP"
+
+    def _flatten(self):
+        return tuple(), (self.wires,)
+
+    @classmethod
+    def _unflatten(cls, data, metadata):
+        return cls(metadata[0])
+
+    def __init__(self, wires, id=None):
+        control_wires = wires[:1]
+        target_wires = wires[1:]
+        super().__init__(qml.SWAP(wires=target_wires), control_wires, id=id)
+
+    @staticmethod
+    @lru_cache()
+    def compute_matrix():  # 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:`~.CSWAP.matrix`
+
+        Returns:
+            ndarray: matrix
+
+        **Example**
+
+        >>> print(qml.CSWAP.compute_matrix())
+        [[1 0 0 0 0 0 0 0]
+         [0 1 0 0 0 0 0 0]
+         [0 0 1 0 0 0 0 0]
+         [0 0 0 1 0 0 0 0]
+         [0 0 0 0 1 0 0 0]
+         [0 0 0 0 0 0 1 0]
+         [0 0 0 0 0 1 0 0]
+         [0 0 0 0 0 0 0 1]]
+        """
+        return np.array(
+            [
+                [1, 0, 0, 0, 0, 0, 0, 0],
+                [0, 1, 0, 0, 0, 0, 0, 0],
+                [0, 0, 1, 0, 0, 0, 0, 0],
+                [0, 0, 0, 1, 0, 0, 0, 0],
+                [0, 0, 0, 0, 1, 0, 0, 0],
+                [0, 0, 0, 0, 0, 0, 1, 0],
+                [0, 0, 0, 0, 0, 1, 0, 0],
+                [0, 0, 0, 0, 0, 0, 0, 1],
+            ]
+        )
+
+    @staticmethod
+    def compute_decomposition(wires):  # pylint: disable=arguments-differ
+        r"""Representation of the operator as a product of other operators (static method).
+
+        .. math:: O = O_1 O_2 \dots O_n.
+
+
+        .. seealso:: :meth:`~.CSWAP.decomposition`.
+
+        Args:
+            wires (Iterable, Wires): wires that the operator acts on
+
+        Returns:
+            list[Operator]: decomposition into lower level operations
+
+        **Example:**
+
+        >>> print(qml.CSWAP.compute_decomposition((0,1,2)))
+        [Toffoli(wires=[0, 2, 1]), Toffoli(wires=[0, 1, 2]), Toffoli(wires=[0, 2, 1])]
+
+        """
+        decomp_ops = [
+            qml.Toffoli(wires=[wires[0], wires[2], wires[1]]),
+            qml.Toffoli(wires=[wires[0], wires[1], wires[2]]),
+            qml.Toffoli(wires=[wires[0], wires[2], wires[1]]),
+        ]
+        return decomp_ops
+
+
+class CCZ(ControlledOp):
+    r"""CCZ(wires)
+    CCZ (controlled-controlled-Z) gate.
+
+    .. math::
+
+        CCZ =
+        \begin{pmatrix}
+        1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\
+        0 & 1 & 0 & 0 & 0 & 0 & 0 & 0\\
+        0 & 0 & 1 & 0 & 0 & 0 & 0 & 0\\
+        0 & 0 & 0 & 1 & 0 & 0 & 0 & 0\\
+        0 & 0 & 0 & 0 & 1 & 0 & 0 & 0\\
+        0 & 0 & 0 & 0 & 0 & 1 & 0 & 0\\
+        0 & 0 & 0 & 0 & 0 & 0 & 1 & 0\\
+        0 & 0 & 0 & 0 & 0 & 0 & 0 & -1
+        \end{pmatrix}
+
+    **Details:**
+
+    * Number of wires: 3
+    * Number of parameters: 0
+
+    Args:
+        wires (Sequence[int]): the subsystem the gate acts on
+    """
+
+    def _flatten(self):
+        return tuple(), (self.wires,)
+
+    @classmethod
+    def _unflatten(cls, data, metadata):
+        return cls(metadata[0])
+
+    num_wires = 3
+    """int: Number of wires that the operator acts on."""
+
+    num_params = 0
+    """int: Number of trainable parameters that the operator depends on."""
+
+    ndim_params = ()
+    """tuple[int]: Number of dimensions per trainable parameter that the operator depends on."""
+
+    name = "CCZ"
+
+    def __init__(self, wires, id=None):
+        control_wires = wires[:2]
+        target_wires = wires[2:]
+        super().__init__(qml.PauliZ(wires=target_wires), control_wires, id=id)
+
+    @staticmethod
+    @lru_cache()
+    def compute_matrix():  # 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:`~.CCZ.matrix`
+
+
+        Returns:
+            ndarray: matrix
+
+        **Example**
+
+        >>> print(qml.CCZ.compute_matrix())
+        [[1 0 0 0 0 0 0 0]
+         [0 1 0 0 0 0 0 0]
+         [0 0 1 0 0 0 0 0]
+         [0 0 0 1 0 0 0 0]
+         [0 0 0 0 1 0 0 0]
+         [0 0 0 0 0 1 0 0]
+         [0 0 0 0 0 0 1 0]
+         [0 0 0 0 0 0 0 -1]]
+        """
+        return np.array(
+            [
+                [1, 0, 0, 0, 0, 0, 0, 0],
+                [0, 1, 0, 0, 0, 0, 0, 0],
+                [0, 0, 1, 0, 0, 0, 0, 0],
+                [0, 0, 0, 1, 0, 0, 0, 0],
+                [0, 0, 0, 0, 1, 0, 0, 0],
+                [0, 0, 0, 0, 0, 1, 0, 0],
+                [0, 0, 0, 0, 0, 0, 1, 0],
+                [0, 0, 0, 0, 0, 0, 0, -1],
+            ]
+        )
+
+    @staticmethod
+    def compute_decomposition(wires):  # pylint: disable=arguments-differ
+        r"""Representation of the operator as a product of other operators (static method).
+
+        .. math:: O = O_1 O_2 \dots O_n.
+
+
+        .. seealso:: :meth:`~.Toffoli.decomposition`.
+
+        Args:
+            wires (Iterable, Wires): wires that the operator acts on
+
+        Returns:
+            list[Operator]: decomposition into lower level operations
+
+        **Example:**
+
+        >>> qml.CCZ.compute_decomposition((0,1,2))
+        [CNOT(wires=[1, 2]),
+        Adjoint(T)(wires=[2]),
+        CNOT(wires=[0, 2]),
+        T(wires=[2]),
+        CNOT(wires=[1, 2]),
+        Adjoint(T)(wires=[2]),
+        CNOT(wires=[0, 2]),
+        T(wires=[2]),
+        T(wires=[1]),
+        CNOT(wires=[0, 1]),
+        Hadamard(wires=[2]),
+        T(wires=[0]),
+        Adjoint(T)(wires=[1]),
+        CNOT(wires=[0, 1]),
+        Hadamard(wires=[2]),]
+
+        """
+        return [
+            qml.CNOT(wires=[wires[1], wires[2]]),
+            qml.adjoint(qml.T(wires=wires[2])),
+            qml.CNOT(wires=[wires[0], wires[2]]),
+            qml.T(wires=wires[2]),
+            qml.CNOT(wires=[wires[1], wires[2]]),
+            qml.adjoint(qml.T(wires=wires[2])),
+            qml.CNOT(wires=[wires[0], wires[2]]),
+            qml.T(wires=wires[2]),
+            qml.T(wires=wires[1]),
+            qml.CNOT(wires=[wires[0], wires[1]]),
+            qml.Hadamard(wires=wires[2]),
+            qml.T(wires=wires[0]),
+            qml.adjoint(qml.T(wires=wires[1])),
+            qml.CNOT(wires=[wires[0], wires[1]]),
+            qml.Hadamard(wires=wires[2]),
+        ]
+
+
+class CNOT(ControlledOp):
+    r"""CNOT(wires)
+    The controlled-NOT operator
+
+    .. math:: CNOT = \begin{bmatrix}
+            1 & 0 & 0 & 0 \\
+            0 & 1 & 0 & 0\\
+            0 & 0 & 0 & 1\\
+            0 & 0 & 1 & 0
+        \end{bmatrix}.
+
+    .. note:: The first wire provided corresponds to the **control qubit**.
+
+    **Details:**
+
+    * Number of wires: 2
+    * Number of parameters: 0
+
+    Args:
+        wires (Sequence[int]): the wires the operation acts on
+    """
+
+    num_wires = 2
+    """int: Number of wires that the operator acts on."""
+
+    num_params = 0
+    """int: Number of trainable parameters that the operator depends on."""
+
+    ndim_params = ()
+    """tuple[int]: Number of dimensions per trainable parameter that the operator depends on."""
+
+    arithmetic_depth = 0
+
+    name = "CNOT"
+
+    def _flatten(self):
+        return tuple(), (self.wires,)
+
+    @classmethod
+    def _unflatten(cls, data, metadata):
+        return cls(metadata[0])
+
+    def __init__(self, wires, id=None):
+        control_wire, wire = wires
+        super().__init__(qml.PauliX(wires=wire), control_wire, id=id)
+
+    def __repr__(self):
+        return f"CNOT(wires={self.wires.tolist()})"
+
+    @staticmethod
+    @lru_cache()
+    def compute_matrix():  # 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:`~.CNOT.matrix`
+
+
+        Returns:
+            ndarray: matrix
+
+        **Example**
+
+        >>> print(qml.CNOT.compute_matrix())
+        [[1 0 0 0]
+         [0 1 0 0]
+         [0 0 0 1]
+         [0 0 1 0]]
+        """
+        return np.array([[1, 0, 0, 0], [0, 1, 0, 0], [0, 0, 0, 1], [0, 0, 1, 0]])
+
+    def _controlled(self, wire):
+        return qml.Toffoli(wires=wire + self.wires)
+
+
+class Toffoli(ControlledOp):
+    r"""Toffoli(wires)
+    Toffoli (controlled-controlled-X) gate.
+
+    .. math::
+
+        Toffoli =
+        \begin{pmatrix}
+        1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\
+        0 & 1 & 0 & 0 & 0 & 0 & 0 & 0\\
+        0 & 0 & 1 & 0 & 0 & 0 & 0 & 0\\
+        0 & 0 & 0 & 1 & 0 & 0 & 0 & 0\\
+        0 & 0 & 0 & 0 & 1 & 0 & 0 & 0\\
+        0 & 0 & 0 & 0 & 0 & 1 & 0 & 0\\
+        0 & 0 & 0 & 0 & 0 & 0 & 0 & 1\\
+        0 & 0 & 0 & 0 & 0 & 0 & 1 & 0
+        \end{pmatrix}
+
+    **Details:**
+
+    * Number of wires: 3
+    * Number of parameters: 0
+
+    Args:
+        wires (Sequence[int]): the subsystem the gate acts on
+    """
+
+    num_wires = 3
+    """int: Number of wires that the operator acts on."""
+
+    num_params = 0
+    """int: Number of trainable parameters that the operator depends on."""
+
+    ndim_params = ()
+    """tuple[int]: Number of dimensions per trainable parameter that the operator depends on."""
+
+    arithmetic_depth = 0
+
+    name = "Toffoli"
+
+    def _flatten(self):
+        return tuple(), (self.wires,)
+
+    @classmethod
+    def _unflatten(cls, _, metadata):
+        return cls(metadata[0])
+
+    def __init__(self, wires, id=None):
+        control_wires = wires[:2]
+        target_wires = wires[2:]
+        super().__init__(qml.PauliX(wires=target_wires), control_wires, id=id)
+
+    def __repr__(self):
+        return f"Toffoli(wires={self.wires.tolist()})"
+
+    @staticmethod
+    @lru_cache()
+    def compute_matrix():  # 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:`~.Toffoli.matrix`
+
+
+        Returns:
+            ndarray: matrix
+
+        **Example**
+
+        >>> print(qml.Toffoli.compute_matrix())
+        [[1 0 0 0 0 0 0 0]
+         [0 1 0 0 0 0 0 0]
+         [0 0 1 0 0 0 0 0]
+         [0 0 0 1 0 0 0 0]
+         [0 0 0 0 1 0 0 0]
+         [0 0 0 0 0 1 0 0]
+         [0 0 0 0 0 0 0 1]
+         [0 0 0 0 0 0 1 0]]
+        """
+        return np.array(
+            [
+                [1, 0, 0, 0, 0, 0, 0, 0],
+                [0, 1, 0, 0, 0, 0, 0, 0],
+                [0, 0, 1, 0, 0, 0, 0, 0],
+                [0, 0, 0, 1, 0, 0, 0, 0],
+                [0, 0, 0, 0, 1, 0, 0, 0],
+                [0, 0, 0, 0, 0, 1, 0, 0],
+                [0, 0, 0, 0, 0, 0, 0, 1],
+                [0, 0, 0, 0, 0, 0, 1, 0],
+            ]
+        )
+
+    @staticmethod
+    def compute_decomposition(wires):  # pylint: disable=arguments-differ
+        r"""Representation of the operator as a product of other operators (static method).
+
+        .. math:: O = O_1 O_2 \dots O_n.
+
+
+        .. seealso:: :meth:`~.Toffoli.decomposition`.
+
+        Args:
+            wires (Iterable, Wires): wires that the operator acts on
+
+        Returns:
+            list[Operator]: decomposition into lower level operations
+
+        **Example:**
+
+        >>> qml.Toffoli.compute_decomposition((0,1,2))
+        [Hadamard(wires=[2]),
+        CNOT(wires=[1, 2]),
+        Adjoint(T)(wires=[2]),
+        CNOT(wires=[0, 2]),
+        T(wires=[2]),
+        CNOT(wires=[1, 2]),
+        Adjoint(T)(wires=[2]),
+        CNOT(wires=[0, 2]),
+        T(wires=[2]),
+        T(wires=[1]),
+        CNOT(wires=[0, 1]),
+        Hadamard(wires=[2]),
+        T(wires=[0]),
+        Adjoint(T)(wires=[1]),
+        CNOT(wires=[0, 1])]
+
+        """
+        return [
+            qml.Hadamard(wires=wires[2]),
+            CNOT(wires=[wires[1], wires[2]]),
+            qml.adjoint(qml.T(wires=wires[2])),
+            CNOT(wires=[wires[0], wires[2]]),
+            qml.T(wires=wires[2]),
+            CNOT(wires=[wires[1], wires[2]]),
+            qml.adjoint(qml.T(wires=wires[2])),
+            CNOT(wires=[wires[0], wires[2]]),
+            qml.T(wires=wires[2]),
+            qml.T(wires=wires[1]),
+            CNOT(wires=[wires[0], wires[1]]),
+            qml.Hadamard(wires=wires[2]),
+            qml.T(wires=wires[0]),
+            qml.adjoint(qml.T(wires=wires[1])),
+            CNOT(wires=[wires[0], wires[1]]),
+        ]
 
 
 def _check_and_convert_control_values(control_values, control_wires):
@@ -725,20 +1311,6 @@ class CRX(ControlledOp):
             qml.RZ(-pi_half, wires=wires[1]),
         ]
 
-    def decomposition(self):
-        r"""Representation of the operator as a product of other operators.
-
-        .. math:: O = O_1 O_2 \dots O_n
-
-        A ``DecompositionUndefinedError`` is raised if no representation by decomposition is defined.
-
-        .. seealso:: :meth:`~.Operator.compute_decomposition`.
-
-        Returns:
-            list[Operator]: decomposition of the operator
-        """
-        return self.compute_decomposition(*self.parameters, wires=self.wires)
-
 
 class CRY(ControlledOp):
     r"""The controlled-RY operator
@@ -879,20 +1451,6 @@ class CRY(ControlledOp):
             qml.CNOT(wires=wires),
         ]
 
-    def decomposition(self):
-        r"""Representation of the operator as a product of other operators.
-
-        .. math:: O = O_1 O_2 \dots O_n
-
-        A ``DecompositionUndefinedError`` is raised if no representation by decomposition is defined.
-
-        .. seealso:: :meth:`~.Operator.compute_decomposition`.
-
-        Returns:
-            list[Operator]: decomposition of the operator
-        """
-        return self.compute_decomposition(*self.parameters, wires=self.wires)
-
 
 class CRZ(ControlledOp):
     r"""The controlled-RZ operator
@@ -1074,20 +1632,6 @@ class CRZ(ControlledOp):
             qml.CNOT(wires=wires),
         ]
 
-    def decomposition(self):
-        r"""Representation of the operator as a product of other operators.
-
-        .. math:: O = O_1 O_2 \dots O_n
-
-        A ``DecompositionUndefinedError`` is raised if no representation by decomposition is defined.
-
-        .. seealso:: :meth:`~.Operator.compute_decomposition`.
-
-        Returns:
-            list[Operator]: decomposition of the operator
-        """
-        return self.compute_decomposition(*self.parameters, wires=self.wires)
-
 
 class CRot(ControlledOp):
     r"""The controlled-Rot operator
@@ -1257,20 +1801,6 @@ class CRot(ControlledOp):
             qml.RZ(omega, wires=wires[1]),
         ]
 
-    def decomposition(self):
-        r"""Representation of the operator as a product of other operators.
-
-        .. math:: O = O_1 O_2 \dots O_n
-
-        A ``DecompositionUndefinedError`` is raised if no representation by decomposition is defined.
-
-        .. seealso:: :meth:`~.Operator.compute_decomposition`.
-
-        Returns:
-            list[Operator]: decomposition of the operator
-        """
-        return self.compute_decomposition(*self.parameters, wires=self.wires)
-
 
 class ControlledPhaseShift(ControlledOp):
     r"""A qubit controlled phase shift.
@@ -1436,19 +1966,5 @@ class ControlledPhaseShift(ControlledOp):
             qml.PhaseShift(phi / 2, wires=wires[1]),
         ]
 
-    def decomposition(self):
-        r"""Representation of the operator as a product of other operators.
-
-        .. math:: O = O_1 O_2 \dots O_n
-
-        A ``DecompositionUndefinedError`` is raised if no representation by decomposition is defined.
-
-        .. seealso:: :meth:`~.Operator.compute_decomposition`.
-
-        Returns:
-            list[Operator]: decomposition of the operator
-        """
-        return self.compute_decomposition(*self.parameters, wires=self.wires)
-
 
 CPhase = ControlledPhaseShift
diff --git a/pennylane/ops/qubit/__init__.py b/pennylane/ops/qubit/__init__.py
index 4f8902255..6a95e789f 100644
--- a/pennylane/ops/qubit/__init__.py
+++ b/pennylane/ops/qubit/__init__.py
@@ -57,17 +57,12 @@ __ops__ = {
     "S",
     "T",
     "SX",
-    "CNOT",
-    "CCZ",
-    "CH",
     "SWAP",
     "ISWAP",
     "SISWAP",
     "SQISW",
-    "CSWAP",
     "PSWAP",
     "ECR",
-    "Toffoli",
     "RX",
     "RY",
     "RZ",
diff --git a/pennylane/ops/qubit/non_parametric_ops.py b/pennylane/ops/qubit/non_parametric_ops.py
index cb0ec6d24..cfe8228a3 100644
--- a/pennylane/ops/qubit/non_parametric_ops.py
+++ b/pennylane/ops/qubit/non_parametric_ops.py
@@ -163,7 +163,7 @@ class Hadamard(Observable, Operation):
         ]
 
     def _controlled(self, wire):
-        return CH(wires=Wires(wire) + self.wires)
+        return qml.CH(wires=Wires(wire) + self.wires)
 
     def adjoint(self):
         return Hadamard(wires=self.wires)
@@ -333,7 +333,7 @@ class PauliX(Observable, Operation):
         return super().pow(z_mod2)
 
     def _controlled(self, wire):
-        return CNOT(wires=Wires(wire) + self.wires)
+        return qml.CNOT(wires=Wires(wire) + self.wires)
 
     def single_qubit_rot_angles(self):
         # X = RZ(-\pi/2) RY(\pi) RZ(\pi/2)
@@ -1012,188 +1012,6 @@ class SX(Operation):
         return [np.pi / 2, np.pi / 2, -np.pi / 2]
 
 
-class CNOT(Operation):
-    r"""CNOT(wires)
-    The controlled-NOT operator
-
-    .. math:: CNOT = \begin{bmatrix}
-            1 & 0 & 0 & 0 \\
-            0 & 1 & 0 & 0\\
-            0 & 0 & 0 & 1\\
-            0 & 0 & 1 & 0
-        \end{bmatrix}.
-
-    .. note:: The first wire provided corresponds to the **control qubit**.
-
-    **Details:**
-
-    * Number of wires: 2
-    * Number of parameters: 0
-
-    Args:
-        wires (Sequence[int]): the wires the operation acts on
-    """
-
-    num_wires = 2
-    num_params = 0
-    """int: Number of trainable parameters that the operator depends on."""
-
-    basis = "X"
-
-    batch_size = None
-
-    def label(self, decimals=None, base_label=None, cache=None):
-        return base_label or "X"
-
-    @staticmethod
-    @lru_cache()
-    def compute_matrix():  # 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:`~.CNOT.matrix`
-
-
-        Returns:
-            ndarray: matrix
-
-        **Example**
-
-        >>> print(qml.CNOT.compute_matrix())
-        [[1 0 0 0]
-         [0 1 0 0]
-         [0 0 0 1]
-         [0 0 1 0]]
-        """
-        return np.array([[1, 0, 0, 0], [0, 1, 0, 0], [0, 0, 0, 1], [0, 0, 1, 0]])
-
-    def adjoint(self):
-        return CNOT(wires=self.wires)
-
-    def pow(self, z):
-        return super().pow(z % 2)
-
-    def _controlled(self, wire):
-        return Toffoli(wires=wire + self.wires)
-
-    @property
-    def control_wires(self):
-        return Wires(self.wires[0])
-
-    @property
-    def is_hermitian(self):
-        return True
-
-
-class CH(Operation):
-    r"""CH(wires)
-    The controlled-Hadamard operator
-
-    .. math:: CH = \begin{bmatrix}
-            1 & 0 & 0 & 0 \\
-            0 & 1 & 0 & 0 \\
-            0 & 0 & \frac{1}{\sqrt{2}} & \frac{1}{\sqrt{2}} \\
-            0 & 0 & \frac{1}{\sqrt{2}} & -\frac{1}{\sqrt{2}}
-        \end{bmatrix}.
-
-    .. note:: The first wire provided corresponds to the **control qubit**.
-
-    **Details:**
-
-    * Number of wires: 2
-    * Number of parameters: 0
-
-    Args:
-        wires (Sequence[int]): the wires the operation acts on
-    """
-
-    num_wires = 2
-    num_params = 0
-    """int: Number of trainable parameters that the operator depends on."""
-
-    basis = "Hadamard"
-
-    batch_size = None
-
-    def label(self, decimals=None, base_label=None, cache=None):
-        return base_label or "H"
-
-    @staticmethod
-    @lru_cache()
-    def compute_matrix():  # 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:`~.CH.matrix`
-
-
-        Returns:
-            ndarray: matrix
-
-        **Example**
-
-        >>> print(qml.CH.compute_matrix())
-        [[ 1.          0.          0.          0.        ]
-         [ 0.          1.          0.          0.        ]
-         [ 0.          0.          0.70710678  0.70710678]
-         [ 0.          0.          0.70710678 -0.70710678]]
-        """
-        return np.array(
-            [
-                [1, 0, 0, 0],
-                [0, 1, 0, 0],
-                [0, 0, INV_SQRT2, INV_SQRT2],
-                [0, 0, INV_SQRT2, -INV_SQRT2],
-            ]
-        )
-
-    @staticmethod
-    def compute_decomposition(wires):
-        r"""Representation of the operator as a product of other operators (static method).
-
-
-        .. math:: O = O_1 O_2 \dots O_n.
-
-
-        .. seealso:: :meth:`~.CH.decomposition`.
-
-        Args:
-            wires (Iterable, Wires): wires that the operator acts on
-
-        Returns:
-            list[Operator]: decomposition into lower level operations
-
-        **Example:**
-
-        >>> print(qml.CH.compute_decomposition([0, 1]))
-        [qml.RY(-0.78539816339, wires=[1]), CZ(wires=[0, 1]), qml.RY(0.78539816339, wires=[1])]
-
-        """
-        return [
-            qml.RY(-np.pi / 4, wires=wires[1]),
-            qml.CZ(wires=wires),
-            qml.RY(+np.pi / 4, wires=wires[1]),
-        ]
-
-    def adjoint(self):
-        return CH(wires=self.wires)
-
-    def pow(self, z):
-        return super().pow(z % 2)
-
-    @property
-    def control_wires(self):
-        return Wires(self.wires[0])
-
-    @property
-    def is_hermitian(self):
-        return True
-
-
 class SWAP(Operation):
     r"""SWAP(wires)
     The swap operator
@@ -1277,7 +1095,7 @@ class SWAP(Operation):
         return SWAP(wires=self.wires)
 
     def _controlled(self, wire):
-        return CSWAP(wires=wire + self.wires)
+        return qml.CSWAP(wires=wire + self.wires)
 
     @property
     def is_hermitian(self):
@@ -1399,7 +1217,7 @@ class ECR(Operation):
         pi = np.pi
         return [
             PauliZ(wires=[wires[0]]),
-            CNOT(wires=[wires[0], wires[1]]),
+            qml.CNOT(wires=[wires[0], wires[1]]),
             SX(wires=[wires[1]]),
             qml.RX(pi / 2, wires=[wires[0]]),
             qml.RY(pi / 2, wires=[wires[0]]),
@@ -1518,8 +1336,8 @@ class ISWAP(Operation):
             S(wires=wires[0]),
             S(wires=wires[1]),
             Hadamard(wires=wires[0]),
-            CNOT(wires=[wires[0], wires[1]]),
-            CNOT(wires=[wires[1], wires[0]]),
+            qml.CNOT(wires=[wires[0], wires[1]]),
+            qml.CNOT(wires=[wires[1], wires[0]]),
             Hadamard(wires=wires[1]),
         ]
 
@@ -1648,14 +1466,14 @@ class SISWAP(Operation):
         return [
             SX(wires=wires[0]),
             qml.RZ(np.pi / 2, wires=wires[0]),
-            CNOT(wires=[wires[0], wires[1]]),
+            qml.CNOT(wires=[wires[0], wires[1]]),
             SX(wires=wires[0]),
             qml.RZ(7 * np.pi / 4, wires=wires[0]),
             SX(wires=wires[0]),
             qml.RZ(np.pi / 2, wires=wires[0]),
             SX(wires=wires[1]),
             qml.RZ(7 * np.pi / 4, wires=wires[1]),
-            CNOT(wires=[wires[0], wires[1]]),
+            qml.CNOT(wires=[wires[0], wires[1]]),
             SX(wires=wires[0]),
             SX(wires=wires[1]),
         ]
@@ -1666,410 +1484,3 @@ class SISWAP(Operation):
 
 
 SQISW = SISWAP
-
-
-class CSWAP(Operation):
-    r"""CSWAP(wires)
-    The controlled-swap operator
-
-    .. math:: CSWAP = \begin{bmatrix}
-            1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\
-            0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 \\
-            0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 \\
-            0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 \\
-            0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 \\
-            0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 \\
-            0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 \\
-            0 & 0 & 0 & 0 & 0 & 0 & 0 & 1
-        \end{bmatrix}.
-
-    .. note:: The first wire provided corresponds to the **control qubit**.
-
-    **Details:**
-
-    * Number of wires: 3
-    * Number of parameters: 0
-
-    Args:
-        wires (Sequence[int]): the wires the operation acts on
-    """
-
-    is_self_inverse = True
-    num_wires = 3
-    num_params = 0
-    """int: Number of trainable parameters that the operator depends on."""
-
-    batch_size = None
-
-    def label(self, decimals=None, base_label=None, cache=None):
-        return base_label or "SWAP"
-
-    @staticmethod
-    @lru_cache()
-    def compute_matrix():  # 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:`~.CSWAP.matrix`
-
-        Returns:
-            ndarray: matrix
-
-        **Example**
-
-        >>> print(qml.CSWAP.compute_matrix())
-        [[1 0 0 0 0 0 0 0]
-         [0 1 0 0 0 0 0 0]
-         [0 0 1 0 0 0 0 0]
-         [0 0 0 1 0 0 0 0]
-         [0 0 0 0 1 0 0 0]
-         [0 0 0 0 0 0 1 0]
-         [0 0 0 0 0 1 0 0]
-         [0 0 0 0 0 0 0 1]]
-        """
-        return np.array(
-            [
-                [1, 0, 0, 0, 0, 0, 0, 0],
-                [0, 1, 0, 0, 0, 0, 0, 0],
-                [0, 0, 1, 0, 0, 0, 0, 0],
-                [0, 0, 0, 1, 0, 0, 0, 0],
-                [0, 0, 0, 0, 1, 0, 0, 0],
-                [0, 0, 0, 0, 0, 0, 1, 0],
-                [0, 0, 0, 0, 0, 1, 0, 0],
-                [0, 0, 0, 0, 0, 0, 0, 1],
-            ]
-        )
-
-    @staticmethod
-    def compute_decomposition(wires):
-        r"""Representation of the operator as a product of other operators (static method).
-
-        .. math:: O = O_1 O_2 \dots O_n.
-
-
-        .. seealso:: :meth:`~.CSWAP.decomposition`.
-
-        Args:
-            wires (Iterable, Wires): wires that the operator acts on
-
-        Returns:
-            list[Operator]: decomposition into lower level operations
-
-        **Example:**
-
-        >>> print(qml.CSWAP.compute_decomposition((0,1,2)))
-        [Toffoli(wires=[0, 2, 1]), Toffoli(wires=[0, 1, 2]), Toffoli(wires=[0, 2, 1])]
-
-        """
-        decomp_ops = [
-            qml.Toffoli(wires=[wires[0], wires[2], wires[1]]),
-            qml.Toffoli(wires=[wires[0], wires[1], wires[2]]),
-            qml.Toffoli(wires=[wires[0], wires[2], wires[1]]),
-        ]
-        return decomp_ops
-
-    def pow(self, z):
-        return super().pow(z % 2)
-
-    def adjoint(self):
-        return CSWAP(wires=self.wires)
-
-    @property
-    def control_wires(self):
-        return Wires(self.wires[0])
-
-    @property
-    def is_hermitian(self):
-        return True
-
-
-class CCZ(Operation):
-    r"""CCZ(wires)
-    CCZ (controlled-controlled-Z) gate.
-
-    .. math::
-
-        CCZ =
-        \begin{pmatrix}
-        1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\
-        0 & 1 & 0 & 0 & 0 & 0 & 0 & 0\\
-        0 & 0 & 1 & 0 & 0 & 0 & 0 & 0\\
-        0 & 0 & 0 & 1 & 0 & 0 & 0 & 0\\
-        0 & 0 & 0 & 0 & 1 & 0 & 0 & 0\\
-        0 & 0 & 0 & 0 & 0 & 1 & 0 & 0\\
-        0 & 0 & 0 & 0 & 0 & 0 & 1 & 0\\
-        0 & 0 & 0 & 0 & 0 & 0 & 0 & -1
-        \end{pmatrix}
-
-    **Details:**
-
-    * Number of wires: 3
-    * Number of parameters: 0
-
-    Args:
-        wires (Sequence[int]): the subsystem the gate acts on
-    """
-
-    num_wires = 3
-    num_params = 0
-    """int: Number of trainable parameters that the operator depends on."""
-
-    basis = "Z"
-
-    batch_size = None
-
-    def label(self, decimals=None, base_label=None, cache=None):
-        return base_label or "Z"
-
-    @staticmethod
-    @lru_cache()
-    def compute_matrix():  # 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:`~.CCZ.matrix`
-
-
-        Returns:
-            ndarray: matrix
-
-        **Example**
-
-        >>> print(qml.CCZ.compute_matrix())
-        [[1 0 0 0 0 0 0 0]
-         [0 1 0 0 0 0 0 0]
-         [0 0 1 0 0 0 0 0]
-         [0 0 0 1 0 0 0 0]
-         [0 0 0 0 1 0 0 0]
-         [0 0 0 0 0 1 0 0]
-         [0 0 0 0 0 0 1 0]
-         [0 0 0 0 0 0 0 -1]]
-        """
-        return np.array(
-            [
-                [1, 0, 0, 0, 0, 0, 0, 0],
-                [0, 1, 0, 0, 0, 0, 0, 0],
-                [0, 0, 1, 0, 0, 0, 0, 0],
-                [0, 0, 0, 1, 0, 0, 0, 0],
-                [0, 0, 0, 0, 1, 0, 0, 0],
-                [0, 0, 0, 0, 0, 1, 0, 0],
-                [0, 0, 0, 0, 0, 0, 1, 0],
-                [0, 0, 0, 0, 0, 0, 0, -1],
-            ]
-        )
-
-    @staticmethod
-    def compute_decomposition(wires):
-        r"""Representation of the operator as a product of other operators (static method).
-
-        .. math:: O = O_1 O_2 \dots O_n.
-
-
-        .. seealso:: :meth:`~.Toffoli.decomposition`.
-
-        Args:
-            wires (Iterable, Wires): wires that the operator acts on
-
-        Returns:
-            list[Operator]: decomposition into lower level operations
-
-        **Example:**
-
-        >>> qml.CCZ.compute_decomposition((0,1,2))
-        [CNOT(wires=[1, 2]),
-        Adjoint(T)(wires=[2]),
-        CNOT(wires=[0, 2]),
-        T(wires=[2]),
-        CNOT(wires=[1, 2]),
-        Adjoint(T)(wires=[2]),
-        CNOT(wires=[0, 2]),
-        T(wires=[2]),
-        T(wires=[1]),
-        CNOT(wires=[0, 1]),
-        Hadamard(wires=[2]),
-        T(wires=[0]),
-        Adjoint(T)(wires=[1]),
-        CNOT(wires=[0, 1]),
-        Hadamard(wires=[2]),]
-
-        """
-        return [
-            CNOT(wires=[wires[1], wires[2]]),
-            qml.adjoint(T(wires=wires[2])),
-            CNOT(wires=[wires[0], wires[2]]),
-            T(wires=wires[2]),
-            CNOT(wires=[wires[1], wires[2]]),
-            qml.adjoint(T(wires=wires[2])),
-            CNOT(wires=[wires[0], wires[2]]),
-            T(wires=wires[2]),
-            T(wires=wires[1]),
-            CNOT(wires=[wires[0], wires[1]]),
-            Hadamard(wires=wires[2]),
-            T(wires=wires[0]),
-            qml.adjoint(T(wires=wires[1])),
-            CNOT(wires=[wires[0], wires[1]]),
-            Hadamard(wires=wires[2]),
-        ]
-
-    def adjoint(self):
-        return CCZ(wires=self.wires)
-
-    def pow(self, z):
-        return super().pow(z % 2)
-
-    @property
-    def control_wires(self):
-        return Wires(self.wires[:2])
-
-    @property
-    def is_hermitian(self):
-        return True
-
-
-class Toffoli(Operation):
-    r"""Toffoli(wires)
-    Toffoli (controlled-controlled-X) gate.
-
-    .. math::
-
-        Toffoli =
-        \begin{pmatrix}
-        1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\
-        0 & 1 & 0 & 0 & 0 & 0 & 0 & 0\\
-        0 & 0 & 1 & 0 & 0 & 0 & 0 & 0\\
-        0 & 0 & 0 & 1 & 0 & 0 & 0 & 0\\
-        0 & 0 & 0 & 0 & 1 & 0 & 0 & 0\\
-        0 & 0 & 0 & 0 & 0 & 1 & 0 & 0\\
-        0 & 0 & 0 & 0 & 0 & 0 & 0 & 1\\
-        0 & 0 & 0 & 0 & 0 & 0 & 1 & 0
-        \end{pmatrix}
-
-    **Details:**
-
-    * Number of wires: 3
-    * Number of parameters: 0
-
-    Args:
-        wires (Sequence[int]): the subsystem the gate acts on
-    """
-
-    num_wires = 3
-    num_params = 0
-    """int: Number of trainable parameters that the operator depends on."""
-
-    basis = "X"
-
-    batch_size = None
-
-    def label(self, decimals=None, base_label=None, cache=None):
-        return base_label or "X"
-
-    @staticmethod
-    @lru_cache()
-    def compute_matrix():  # 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:`~.Toffoli.matrix`
-
-
-        Returns:
-            ndarray: matrix
-
-        **Example**
-
-        >>> print(qml.Toffoli.compute_matrix())
-        [[1 0 0 0 0 0 0 0]
-         [0 1 0 0 0 0 0 0]
-         [0 0 1 0 0 0 0 0]
-         [0 0 0 1 0 0 0 0]
-         [0 0 0 0 1 0 0 0]
-         [0 0 0 0 0 1 0 0]
-         [0 0 0 0 0 0 0 1]
-         [0 0 0 0 0 0 1 0]]
-        """
-        return np.array(
-            [
-                [1, 0, 0, 0, 0, 0, 0, 0],
-                [0, 1, 0, 0, 0, 0, 0, 0],
-                [0, 0, 1, 0, 0, 0, 0, 0],
-                [0, 0, 0, 1, 0, 0, 0, 0],
-                [0, 0, 0, 0, 1, 0, 0, 0],
-                [0, 0, 0, 0, 0, 1, 0, 0],
-                [0, 0, 0, 0, 0, 0, 0, 1],
-                [0, 0, 0, 0, 0, 0, 1, 0],
-            ]
-        )
-
-    @staticmethod
-    def compute_decomposition(wires):
-        r"""Representation of the operator as a product of other operators (static method).
-
-        .. math:: O = O_1 O_2 \dots O_n.
-
-
-        .. seealso:: :meth:`~.Toffoli.decomposition`.
-
-        Args:
-            wires (Iterable, Wires): wires that the operator acts on
-
-        Returns:
-            list[Operator]: decomposition into lower level operations
-
-        **Example:**
-
-        >>> qml.Toffoli.compute_decomposition((0,1,2))
-        [Hadamard(wires=[2]),
-        CNOT(wires=[1, 2]),
-        Adjoint(T)(wires=[2]),
-        CNOT(wires=[0, 2]),
-        T(wires=[2]),
-        CNOT(wires=[1, 2]),
-        Adjoint(T)(wires=[2]),
-        CNOT(wires=[0, 2]),
-        T(wires=[2]),
-        T(wires=[1]),
-        CNOT(wires=[0, 1]),
-        Hadamard(wires=[2]),
-        T(wires=[0]),
-        Adjoint(T)(wires=[1]),
-        CNOT(wires=[0, 1])]
-
-        """
-        return [
-            Hadamard(wires=wires[2]),
-            CNOT(wires=[wires[1], wires[2]]),
-            qml.adjoint(T(wires=wires[2])),
-            CNOT(wires=[wires[0], wires[2]]),
-            T(wires=wires[2]),
-            CNOT(wires=[wires[1], wires[2]]),
-            qml.adjoint(T(wires=wires[2])),
-            CNOT(wires=[wires[0], wires[2]]),
-            T(wires=wires[2]),
-            T(wires=wires[1]),
-            CNOT(wires=[wires[0], wires[1]]),
-            Hadamard(wires=wires[2]),
-            T(wires=wires[0]),
-            qml.adjoint(T(wires=wires[1])),
-            CNOT(wires=[wires[0], wires[1]]),
-        ]
-
-    def adjoint(self):
-        return Toffoli(wires=self.wires)
-
-    def pow(self, z):
-        return super().pow(z % 2)
-
-    @property
-    def control_wires(self):
-        return Wires(self.wires[:2])
-
-    @property
-    def is_hermitian(self):
-        return True
diff --git a/pennylane/transforms/decompositions/clifford_t_transform.py b/pennylane/transforms/decompositions/clifford_t_transform.py
index 852ea83a1..6620e95ba 100644
--- a/pennylane/transforms/decompositions/clifford_t_transform.py
+++ b/pennylane/transforms/decompositions/clifford_t_transform.py
@@ -19,6 +19,7 @@ from itertools import product
 from typing import Sequence, Callable
 
 import pennylane as qml
+from pennylane.ops import Adjoint
 from pennylane.queuing import QueuingManager
 from pennylane.transforms.core import transform
 from pennylane.tape import QuantumTape
@@ -130,7 +131,10 @@ def check_clifford_t(op, use_decomposition=False):
         Bool that represents whether the provided operator is Clifford+T or not.
     """
     # Get the base operation for an adjointed operation
-    base = getattr(op, "base", None)
+    if isinstance(op, Adjoint):
+        base = op.base
+    else:
+        base = None
 
     # Save time and check from the pre-computed list
     if isinstance(op, _CLIFFORD_T_GATES) or isinstance(base, _CLIFFORD_T_GATES):


Test patch:
diff --git a/tests/drawer/test_drawer_utils.py b/tests/drawer/test_drawer_utils.py
index 902b18e92..0458d0a39 100644
--- a/tests/drawer/test_drawer_utils.py
+++ b/tests/drawer/test_drawer_utils.py
@@ -137,8 +137,8 @@ class TestUnwrapControls:
         "op,expected_control_wires,expected_control_values",
         [
             (qml.PauliX(wires="a"), Wires([]), None),
-            (qml.CNOT(wires=["a", "b"]), Wires("a"), None),
-            (qml.ctrl(qml.PauliX(wires="b"), control="a"), Wires("a"), None),
+            (qml.CNOT(wires=["a", "b"]), Wires("a"), [True]),
+            (qml.ctrl(qml.PauliX(wires="b"), control="a"), Wires("a"), [True]),
             (
                 qml.ctrl(qml.PauliX(wires="b"), control=["a", "c", "d"]),
                 Wires(["a", "c", "d"]),
diff --git a/tests/ops/op_math/test_pow_op.py b/tests/ops/op_math/test_pow_op.py
index ba0016ba6..934380030 100644
--- a/tests/ops/op_math/test_pow_op.py
+++ b/tests/ops/op_math/test_pow_op.py
@@ -468,8 +468,7 @@ class TestSimplify:
         final_op = qml.CH([1, 0], id=3)
         simplified_op = pow_op.simplify()
 
-        # TODO: uncomment this check when all controlled operations inherit from ControlledOp
-        # assert isinstance(simplified_op, ControlledOp)
+        assert isinstance(simplified_op, ControlledOp)
         assert final_op.data == simplified_op.data
         assert final_op.wires == simplified_op.wires
         assert final_op.arithmetic_depth == simplified_op.arithmetic_depth
diff --git a/tests/ops/qubit/test_non_parametric_ops.py b/tests/ops/qubit/test_non_parametric_ops.py
index eca827f4b..95bd2e579 100644
--- a/tests/ops/qubit/test_non_parametric_ops.py
+++ b/tests/ops/qubit/test_non_parametric_ops.py
@@ -935,7 +935,7 @@ class TestPowMethod:
     @pytest.mark.parametrize("op", period_two_ops)
     def test_period_two_noninteger_power(self, op):
         """Test that ops with a period of 2 raised to a non-integer power raise an error."""
-        if op.__class__ in [qml.PauliZ, qml.CZ]:
+        if op.__class__ in [qml.PauliZ, qml.CZ, qml.CCZ]:
             pytest.skip("PauliZ can be raised to any power.")
         with pytest.raises(qml.operation.PowUndefinedError):
             op.pow(1.234)
diff --git a/tests/templates/test_subroutines/test_single_excitation.py b/tests/templates/test_subroutines/test_single_excitation.py
index 76eda25a2..1a898036b 100644
--- a/tests/templates/test_subroutines/test_single_excitation.py
+++ b/tests/templates/test_subroutines/test_single_excitation.py
@@ -102,7 +102,7 @@ class TestDecomposition:
             assert isinstance(res_gate, exp_gate)
 
             exp_wires = gate[2]
-            res_wires = queue[idx]._wires
+            res_wires = queue[idx].wires
             assert res_wires.tolist() == exp_wires
 
             exp_weight = gate[3]
diff --git a/tests/templates/test_subroutines/test_uccsd.py b/tests/templates/test_subroutines/test_uccsd.py
index 22ed28965..827abb4e1 100644
--- a/tests/templates/test_subroutines/test_uccsd.py
+++ b/tests/templates/test_subroutines/test_uccsd.py
@@ -139,7 +139,7 @@ class TestDecomposition:
             assert isinstance(res_gate, exp_gate)
 
             exp_wires = gate[2]
-            res_wires = queue[idx]._wires
+            res_wires = queue[idx].wires
             assert res_wires.tolist() == exp_wires
 
             exp_weight = gate[3]
diff --git a/tests/transforms/test_cliffordt_transform.py b/tests/transforms/test_cliffordt_transform.py
index 9020a3425..fc014ae1f 100644
--- a/tests/transforms/test_cliffordt_transform.py
+++ b/tests/transforms/test_cliffordt_transform.py
@@ -90,6 +90,7 @@ class TestCliffordCompile:
             (qml.DoubleExcitation(2.0, wires=[0, 1, 2, 3]), False),
             (qml.PauliX(wires=[1]), True),
             (qml.RX(3 * PI, wires=[1]), True),
+            (qml.adjoint(qml.RX(3 * PI, wires=[1])), True),
             (qml.PhaseShift(2 * PI, wires=["a"]), True),
             (qml.ECR(wires=["e", "f"]), True),
             (qml.CH(wires=["a", "b"]), False),
diff --git a/tests/transforms/test_commutation_dag.py b/tests/transforms/test_commutation_dag.py
index 58b32a79b..a33bf3114 100644
--- a/tests/transforms/test_commutation_dag.py
+++ b/tests/transforms/test_commutation_dag.py
@@ -142,19 +142,24 @@ class TestCommutationDAG:
         "Test a the DAG and its attributes for a more complicated circuit."
         # pylint: disable=protected-access
 
+        op_wires = [
+            (qml.CNOT, [3, 0]),
+            (qml.PauliX, [4]),
+            (qml.PauliZ, [0]),
+            (qml.CNOT, [4, 2]),
+            (qml.CNOT, [0, 1]),
+            (qml.CNOT, [3, 4]),
+            (qml.CNOT, [1, 2]),
+            (qml.PauliX, [1]),
+            (qml.CNOT, [1, 0]),
+            (qml.PauliX, [1]),
+            (qml.CNOT, [1, 2]),
+            (qml.CNOT, [0, 3]),
+        ]
+
         def circuit():
-            qml.CNOT(wires=[3, 0])
-            qml.PauliX(wires=4)
-            qml.PauliZ(wires=0)
-            qml.CNOT(wires=[4, 2])
-            qml.CNOT(wires=[0, 1])
-            qml.CNOT(wires=[3, 4])
-            qml.CNOT(wires=[1, 2])
-            qml.PauliX(wires=1)
-            qml.CNOT(wires=[1, 0])
-            qml.PauliX(wires=1)
-            qml.CNOT(wires=[1, 2])
-            qml.CNOT(wires=[0, 3])
+            for op, ws in op_wires:
+                op(wires=ws)
 
         dag = qml.transforms.commutation_dag(circuit)()
 
@@ -162,23 +167,7 @@ class TestCommutationDAG:
         consecutive_wires = Wires(range(len(wires)))
         wires_map = OrderedDict(zip(wires, consecutive_wires))
 
-        nodes = [
-            qml.CNOT(wires=[3, 0]),
-            qml.PauliX(wires=4),
-            qml.PauliZ(wires=0),
-            qml.CNOT(wires=[4, 2]),
-            qml.CNOT(wires=[0, 1]),
-            qml.CNOT(wires=[3, 4]),
-            qml.CNOT(wires=[1, 2]),
-            qml.PauliX(wires=1),
-            qml.CNOT(wires=[1, 0]),
-            qml.PauliX(wires=1),
-            qml.CNOT(wires=[1, 2]),
-            qml.CNOT(wires=[0, 3]),
-        ]
-
-        for node in nodes:
-            node._wires = Wires([wires_map[wire] for wire in node.wires.tolist()])
+        nodes = list(op(wires=[wires_map[w] for w in ws]) for op, ws in op_wires)
 
         edges = [
             (0, 2, {"commute": False}),
diff --git a/tests/transforms/test_tape_expand.py b/tests/transforms/test_tape_expand.py
index 7dfd56ca3..861ae5337 100644
--- a/tests/transforms/test_tape_expand.py
+++ b/tests/transforms/test_tape_expand.py
@@ -388,7 +388,7 @@ class TestExpandInvalidTrainable:
 
 
 # Custom decomposition functions for testing.
-def custom_cnot(wires):
+def custom_cnot(wires, **_):
     return [
         qml.Hadamard(wires=wires[1]),
         qml.CZ(wires=[wires[0], wires[1]]),