Repository: PennyLaneAI/pennylane
Pull request #3449 (2022-11-30 21:37:12)

[Return-types Fix] Update metric tensor transforms
**Context:**
Part of the ongoing effort to update the `batch_transforms` for the new return types system.

**Description of the Change:**
Update the `metric_tensor`, `adjoint_metric_tensor`, and `classical_fisher` transforms.


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/3449?src=pr&el=h1&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI) Report
> Merging [#3449](https://codecov.io/gh/PennyLaneAI/pennylane/pull/3449?src=pr&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI) (cf65aaa) into [master](https://codecov.io/gh/PennyLaneAI/pennylane/commit/287156194ac0511015dd7cbbe7f6d0a07e1d3564?el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI) (2871561) will **increase** coverage by `0.00%`.
> The diff coverage is `100.00%`.

```diff
@@           Coverage Diff           @@
##           master    #3449   +/-   ##
=======================================
  Coverage   99.68%   99.68%           
=======================================
  Files         299      299           
  Lines       26755    26759    +4     
=======================================
+ Hits        26670    26674    +4     
  Misses         85       85           
```


| [Impacted Files](https://codecov.io/gh/PennyLaneAI/pennylane/pull/3449?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI) | Coverage Δ | |
|---|---|---|
| [pennylane/qinfo/transforms.py](https://codecov.io/gh/PennyLaneAI/pennylane/pull/3449/diff?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI#diff-cGVubnlsYW5lL3FpbmZvL3RyYW5zZm9ybXMucHk=) | `100.00% <100.00%> (ø)` | |
| [pennylane/transforms/metric\_tensor.py](https://codecov.io/gh/PennyLaneAI/pennylane/pull/3449/diff?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI#diff-cGVubnlsYW5lL3RyYW5zZm9ybXMvbWV0cmljX3RlbnNvci5weQ==) | `100.00% <100.00%> (ø)` | |

Help us with your feedback. Take ten seconds to tell us [how you rate us](https://about.codecov.io/nps?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI). Have a feature suggestion? [Share it here.](https://app.codecov.io/gh/feedback/?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=PennyLaneAI)

eddddddy: @rmoyard Yes sorry, I forgot to mention that none of the test files were modified.
eddddddy: Not a lot of changes, so I think 1 reviewer is enough.

Patch:
diff --git a/doc/releases/changelog-dev.md b/doc/releases/changelog-dev.md
index 8582a09cb..b8b391df3 100644
--- a/doc/releases/changelog-dev.md
+++ b/doc/releases/changelog-dev.md
@@ -188,15 +188,15 @@
 * Made `gradients.finite_diff` more convenient to use with custom data type observables/devices.
   [(#3426)](https://github.com/PennyLaneAI/pennylane/pull/3426)
 
-* The `qml.ISWAP` gate is now natively supported on `default.mixed`, improving on its efficiency. 
+* The `qml.ISWAP` gate is now natively supported on `default.mixed`, improving on its efficiency.
   [(#3284)](https://github.com/PennyLaneAI/pennylane/pull/3284)
-  
+
 * Added more input validation to `hamiltonian_expand` such that Hamiltonian objects with no terms raise an error.
   [(#3339)](https://github.com/PennyLaneAI/pennylane/pull/3339)
 
 * Continuous integration checks are now performed for Python 3.11 and Torch v1.13. Python 3.7 is dropped.
   [(#3276)](https://github.com/PennyLaneAI/pennylane/pull/3276)
-  
+
 * `qml.Tracker` now also logs results in `tracker.history` when tracking execution of a circuit.
    [(#3306)](https://github.com/PennyLaneAI/pennylane/pull/3306)
 
@@ -351,6 +351,10 @@
 * Updated `qml.transforms.mitigate_with_zne` to support the new return types.
   [(#3415)](https://github.com/PennyLaneAI/pennylane/pull/3415)
 
+* Updated `qml.transforms.metric_tensor`, `qml.transforms.adjoint_metric_tensor`,
+  `qml.qinfo.classical_fisher`, and `qml.qinfo.quantum_fisher` to support the new return types.
+  [(#3449)](https://github.com/PennyLaneAI/pennylane/pull/3449)
+
 
 <h3>Breaking changes</h3>
 
diff --git a/pennylane/qinfo/transforms.py b/pennylane/qinfo/transforms.py
index 7155d52a1..ff451ce95 100644
--- a/pennylane/qinfo/transforms.py
+++ b/pennylane/qinfo/transforms.py
@@ -313,8 +313,12 @@ def _make_probs(tape, wires=None, post_processing_fn=None):
 
     if post_processing_fn is None:
 
-        def post_processing_fn(x):
-            return qml.math.squeeze(qml.math.stack(x))
+        def post_processing_fn(res):
+            if qml.active_return():
+                # only a single probs measurement, so no stacking needed
+                return res[0]
+
+            return qml.math.squeeze(qml.math.stack(res))
 
     return [new_tape], post_processing_fn
 
@@ -462,7 +466,7 @@ def classical_fisher(qnode, argnums=0):
             [2.16840434e-18, 2.81967252e-01]]]))
 
     """
-    new_qnode = _make_probs(qnode, post_processing_fn=lambda x: qml.math.squeeze(qml.math.stack(x)))
+    new_qnode = _make_probs(qnode)
 
     interface = qnode.interface
 
diff --git a/pennylane/transforms/metric_tensor.py b/pennylane/transforms/metric_tensor.py
index 99950d3ee..cb0535a2d 100644
--- a/pennylane/transforms/metric_tensor.py
+++ b/pennylane/transforms/metric_tensor.py
@@ -319,7 +319,9 @@ def qnode_execution_wrapper(self, qnode, targs, tkwargs):
     tkwargs.setdefault("device_wires", qnode.device.wires)
     mt_fn = self.default_qnode_wrapper(qnode, targs, tkwargs)
 
-    _expand_fn = lambda tape: self.expand_fn(tape, *targs, **tkwargs)
+    def _expand_fn(tape):
+        return self.expand_fn(tape, *targs, **tkwargs)
+
     cjac_fn = qml.transforms.classical_jacobian(qnode, expand_fn=_expand_fn)
 
     def wrapper(*args, **kwargs):
@@ -351,7 +353,7 @@ def qnode_execution_wrapper(self, qnode, targs, tkwargs):
                     "much more efficient to request the block-diagonal approximation directly!"
                 )
             tkwargs["approx"] = "block-diag"
-            return self.__call__(qnode, *targs, **tkwargs)(*args, **kwargs)
+            return self(qnode, *targs, **tkwargs)(*args, **kwargs)
 
         if not hybrid:
             return mt
@@ -611,6 +613,10 @@ def _metric_tensor_hadamard(tape, allow_nonunitary, aux_wire, device_wires):
         # Get full block-diagonal tensor
         diag_mt = diag_proc_fn(diag_res)
 
+        if qml.active_return():
+            # the off diag tapes only have a single expval measurement
+            off_diag_res = [qml.math.expand_dims(res, 0) for res in off_diag_res]
+
         # Prepare the mask to match the used interface
         mask = qml.math.convert_like(mask, diag_mt)
 


Test patch:
diff --git a/tests/returntypes/transforms/test_adjoint_metric_tensor_new.py b/tests/returntypes/transforms/test_adjoint_metric_tensor_new.py
new file mode 100644
index 000000000..dd1a12cb8
--- /dev/null
+++ b/tests/returntypes/transforms/test_adjoint_metric_tensor_new.py
@@ -0,0 +1,863 @@
+# 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
+
+#     http://www.apache.org/licenses/LICENSE-2.0
+
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+"""
+Unit tests for the adjoint_metric_tensor function.
+"""
+import pytest
+from pennylane import numpy as np
+import pennylane as qml
+
+from pennylane.transforms.adjoint_metric_tensor import _apply_operations
+
+
+class TestApplyOperations:
+    """Tests the application of operations via the helper function
+    _apply_operations used in the adjoint metric tensor."""
+
+    device = qml.device("default.qubit", wires=2)
+    x = 0.5
+
+    def test_simple_operation(self):
+        """Test that an operation is applied correctly."""
+        op = qml.RX(self.x, wires=0)
+        out = _apply_operations(self.device._state, op, self.device)
+        out = qml.math.reshape(out, 4)
+        exp = np.array([np.cos(self.x / 2), 0.0, -1j * np.sin(self.x / 2), 0.0])
+        assert np.allclose(out, exp)
+        assert not op.inverse
+
+    def test_simple_operation_inv(self):
+        """Test that an operation is applied correctly when using invert=True
+        but does not alter the operation (in particular its inverse flag) that is used."""
+        op = qml.RX(self.x, wires=0)
+        out = _apply_operations(self.device._state, op, self.device, invert=True)
+        out = qml.math.reshape(out, 4)
+        exp = np.array([np.cos(self.x / 2), 0.0, 1j * np.sin(self.x / 2), 0.0])
+        assert np.allclose(out, exp)
+        assert not op.inverse
+
+    def test_inv_operation(self):
+        """Test that an operation with inverse=True is applied correctly
+        but does not alter the operation (in particular its inverse flag) that is used."""
+        op = qml.RX(self.x, wires=0).inv()
+        out = _apply_operations(self.device._state, op, self.device)
+        out = qml.math.reshape(out, 4)
+        exp = np.array([np.cos(self.x / 2), 0.0, 1j * np.sin(self.x / 2), 0.0])
+        assert np.allclose(out, exp)
+        assert op.inverse
+
+    def test_inv_operation_inv(self):
+        """Test that an operation with inverse=True is applied correctly when using invert=True
+        but does not alter the operation (in particular its inverse flag) that is used."""
+        op = qml.RX(self.x, wires=0).inv()
+        out = _apply_operations(self.device._state, op, self.device, invert=True)
+        out = qml.math.reshape(out, 4)
+        exp = np.array([np.cos(self.x / 2), 0.0, -1j * np.sin(self.x / 2), 0.0])
+        assert np.allclose(out, exp)
+        assert op.inverse
+
+    def test_operation_group(self):
+        """Test that a group of operations with is applied correctly
+        but does not alter the operations (in particular their order and
+        inverse flags) that are used."""
+        op = [qml.RX(self.x, wires=0).inv(), qml.Hadamard(wires=1), qml.CNOT(wires=[1, 0])]
+        out = _apply_operations(self.device._state, op, self.device)
+        out = qml.math.reshape(out, 4)
+        exp = np.array(
+            [
+                np.cos(self.x / 2) / np.sqrt(2),
+                1j * np.sin(self.x / 2) / np.sqrt(2),
+                1j * np.sin(self.x / 2) / np.sqrt(2),
+                np.cos(self.x / 2) / np.sqrt(2),
+            ]
+        )
+        assert np.allclose(out, exp)
+        assert isinstance(op[0], qml.RX) and op[0].inverse
+        assert isinstance(op[1], qml.Hadamard) and not op[1].inverse
+        assert isinstance(op[2], qml.CNOT) and not op[2].inverse
+
+    def test_operation_group_inv(self):
+        """Test that a group of operations with is applied correctly when using invert=True
+        but does not alter the operations (in particular their order and
+        inverse flags) that are used."""
+        op = [qml.RX(self.x, wires=0).inv(), qml.Hadamard(wires=1), qml.CNOT(wires=[1, 0])]
+        out = _apply_operations(self.device._state, op, self.device, invert=True)
+        out = qml.math.reshape(out, 4)
+        exp = np.array(
+            [
+                np.cos(self.x / 2) / np.sqrt(2),
+                np.cos(self.x / 2) / np.sqrt(2),
+                -1j * np.sin(self.x / 2) / np.sqrt(2),
+                -1j * np.sin(self.x / 2) / np.sqrt(2),
+            ]
+        )
+        assert np.allclose(out, exp)
+        assert isinstance(op[0], qml.RX) and op[0].inverse
+        assert isinstance(op[1], qml.Hadamard) and not op[1].inverse
+        assert isinstance(op[2], qml.CNOT) and not op[2].inverse
+
+    def test_qubit_statevector(self):
+        """Test that a statevector preparation is applied correctly."""
+        state = np.array([0.4, 1.2 - 0.2j, 9.5, -0.3 + 1.1j])
+        state /= np.linalg.norm(state, ord=2)
+        op = qml.QubitStateVector(state, wires=self.device.wires)
+        out = _apply_operations(None, op, self.device, invert=False)
+        out = qml.math.reshape(out, 4)
+        assert np.allclose(out, state)
+
+    def test_error_qubit_statevector(self):
+        """Test that an error is raised for a statevector preparation with invert=True."""
+        state = np.array([0.4, 1.2 - 0.2j, 9.5, -0.3 + 1.1j])
+        state = np.array([0.4, 1.2 - 0.2j, 9.5, -0.3 + 1.1j])
+        state /= np.linalg.norm(state, ord=2)
+        op = qml.QubitStateVector(state, wires=self.device.wires)
+        with pytest.raises(ValueError, match="Can't invert state preparation."):
+            _apply_operations(None, op, self.device, invert=True)
+
+    def test_basisstate(self):
+        """Test that a basis state preparation is applied correctly."""
+        op = qml.BasisState(np.array([1, 0]), wires=self.device.wires)
+        out = _apply_operations(None, op, self.device, invert=False)
+        out = qml.math.reshape(out, 4)
+        exp = np.array([0.0, 0.0, 1.0, 0.0])
+        assert np.allclose(out, exp)
+
+    def test_error_basisstate(self):
+        """Test that an error is raised for a basis state preparation with invert=True."""
+        op = qml.BasisState(np.array([1, 0]), wires=self.device.wires)
+        with pytest.raises(ValueError, match="Can't invert state preparation."):
+            _apply_operations(None, op, self.device, invert=True)
+
+
+@pytest.mark.parametrize("invert", [False, True])
+class TestApplyOperationsDifferentiability:
+    """Tests the differentiability of applying operations via the helper function
+    _apply_operations used in the adjoint metric tensor."""
+
+    x = 0.5
+
+    @pytest.mark.autograd
+    def test_simple_operation_autograd(self, invert):
+        """Test differentiability for a simple operation with Autograd."""
+        device = qml.device("default.qubit.autograd", wires=2)
+        x = np.array(self.x, requires_grad=True)
+        r_fn = lambda x: qml.math.real(
+            _apply_operations(device._state, qml.RX(x, wires=0), device, invert)
+        )
+        i_fn = lambda x: qml.math.imag(
+            _apply_operations(device._state, qml.RX(x, wires=0), device, invert)
+        )
+        out = qml.jacobian(r_fn)(x) + 1j * qml.jacobian(i_fn)(x)
+        exp = (
+            np.array([[-np.sin(self.x / 2), 0.0], [-1j * (-1) ** invert * np.cos(self.x / 2), 0.0]])
+            / 2
+        )
+        assert np.allclose(out, exp)
+
+    @pytest.mark.jax
+    def test_simple_operation_jax(self, invert):
+        """Test differentiability for a simple operation with JAX."""
+        import jax
+
+        device = qml.device("default.qubit.jax", wires=2)
+        x = jax.numpy.array(self.x)
+        r_fn = lambda x: qml.math.real(
+            _apply_operations(device._state, qml.RX(x, wires=0), device, invert)
+        )
+        i_fn = lambda x: qml.math.imag(
+            _apply_operations(device._state, qml.RX(x, wires=0), device, invert)
+        )
+        out = jax.jacobian(r_fn)(x) + 1j * jax.jacobian(i_fn)(x)
+        exp = (
+            np.array([[-np.sin(self.x / 2), 0.0], [-1j * (-1) ** invert * np.cos(self.x / 2), 0.0]])
+            / 2
+        )
+        assert np.allclose(out, exp)
+
+    @pytest.mark.tf
+    def test_simple_operation_tf(self, invert):
+        """Test differentiability for a simple operation with TensorFlow."""
+        import tensorflow as tf
+
+        device = qml.device("default.qubit.tf", wires=2)
+        x = tf.Variable(self.x, dtype=tf.float64)
+        r_fn = lambda x: qml.math.real(
+            _apply_operations(device._state, qml.RX(x, wires=0), device, invert)
+        )
+        i_fn = lambda x: qml.math.imag(
+            _apply_operations(device._state, qml.RX(x, wires=0), device, invert)
+        )
+        with tf.GradientTape(persistent=True) as tape:
+            r_state = r_fn(x)
+            i_state = i_fn(x)
+        out = qml.math.complex(tape.jacobian(r_state, x), tape.jacobian(i_state, x))
+        exp = (
+            np.array([[-np.sin(self.x / 2), 0.0], [-1j * (-1) ** invert * np.cos(self.x / 2), 0.0]])
+            / 2
+        )
+        assert np.allclose(out, exp)
+
+    @pytest.mark.torch
+    def test_simple_operation_torch(self, invert):
+        """Test differentiability for a simple operation with Torch."""
+        import torch
+
+        jac_fn = torch.autograd.functional.jacobian
+        device = qml.device("default.qubit.torch", wires=2)
+        x = torch.tensor(self.x, requires_grad=True)
+        r_fn = lambda x: qml.math.real(
+            _apply_operations(device._state, qml.RX(x, wires=0), device, invert)
+        )
+        i_fn = lambda x: qml.math.imag(
+            _apply_operations(device._state, qml.RX(x, wires=0), device, invert)
+        )
+        out = jac_fn(r_fn, x) + 1j * jac_fn(i_fn, x)
+        exp = (
+            np.array([[-np.sin(self.x / 2), 0.0], [-1j * (-1) ** invert * np.cos(self.x / 2), 0.0]])
+            / 2
+        )
+        assert np.allclose(out, exp)
+
+
+fixed_pars = [-0.2, 0.2, 0.5, 0.3, 0.7]
+
+
+def fubini_ansatz0(params, wires=None):
+    qml.RX(params[0], wires=0)
+    qml.RY(fixed_pars[0], wires=0)
+    qml.CNOT(wires=[wires[0], wires[1]])
+    qml.RZ(params[1], wires=0)
+    qml.CNOT(wires=[wires[0], wires[1]])
+
+
+def fubini_ansatz1(params, wires=None):
+    qml.RX(fixed_pars[1], wires=0)
+    for wire in wires:
+        qml.Rot(*params[0][wire], wires=wire)
+    qml.CNOT(wires=[0, 1])
+    qml.RY(fixed_pars[1], wires=0).inv()
+    qml.CNOT(wires=[1, 2])
+    for wire in wires:
+        qml.Rot(*params[1][wire], wires=wire)
+    qml.CNOT(wires=[1, 2])
+    qml.RX(fixed_pars[2], wires=1)
+
+
+def fubini_ansatz2(params, wires=None):
+    params0 = params[0]
+    params1 = params[1]
+    qml.RX(fixed_pars[1], wires=0)
+    qml.Rot(*fixed_pars[2:5], wires=1)
+    qml.CNOT(wires=[0, 1])
+    qml.RY(params0, wires=0)
+    qml.RY(params0, wires=1)
+    qml.CNOT(wires=[0, 1])
+    qml.RX(params1, wires=0).inv()
+    qml.RX(params1, wires=1)
+
+
+def fubini_ansatz3(params, wires=None):
+    params0 = params[0]
+    params1 = params[1]
+    params2 = params[2]
+    qml.RX(fixed_pars[1], wires=0)
+    qml.RX(fixed_pars[3], wires=1)
+    qml.CNOT(wires=[0, 1])
+    qml.CNOT(wires=[1, 2])
+    qml.RX(params0, wires=0)
+    qml.RX(params0, wires=1)
+    qml.CNOT(wires=[0, 1])
+    qml.CNOT(wires=[1, 2])
+    qml.CNOT(wires=[2, 0])
+    qml.RY(params1, wires=0)
+    qml.RY(params1, wires=1)
+    qml.RY(params1, wires=2)
+    qml.RZ(params2, wires=0)
+    qml.RZ(params2, wires=1)
+    qml.RZ(params2, wires=2)
+
+
+def fubini_ansatz4(params00, params_rest, wires=None):
+    params01 = params_rest[0]
+    params10 = params_rest[1]
+    params11 = params_rest[2]
+    qml.RY(fixed_pars[3], wires=0)
+    qml.RY(fixed_pars[2], wires=1)
+    qml.CNOT(wires=[0, 1])
+    qml.CNOT(wires=[1, 2])
+    qml.RY(fixed_pars[4], wires=0)
+    qml.RX(params00, wires=0)
+    qml.CNOT(wires=[0, 1])
+    qml.RX(params01, wires=1)
+    qml.RZ(params10, wires=1)
+    qml.CNOT(wires=[0, 1])
+    qml.RZ(params11, wires=1)
+
+
+def fubini_ansatz5(params, wires=None):
+    fubini_ansatz4(params[0], [params[0], params[1], params[1]], wires=wires)
+
+
+def fubini_ansatz6(params, wires=None):
+    fubini_ansatz4(params[0], [params[0], params[1], -params[1]], wires=wires)
+
+
+def fubini_ansatz7(params0, params1, wires=None):
+    fubini_ansatz4(params0[0], [params0[1], params1[0], params1[1]], wires=wires)
+
+
+def fubini_ansatz8(x, wires=None):
+    qml.RX(fixed_pars[0], wires=0)
+    qml.RX(x, wires=0)
+
+
+def fubini_ansatz9(params, wires=None):
+    params0 = params[0]
+    params1 = params[1]
+    qml.RX(fixed_pars[1], wires=[0])
+    qml.RY(fixed_pars[3], wires=[0])
+    qml.RZ(fixed_pars[2], wires=[0])
+    qml.RX(fixed_pars[2], wires=[1])
+    qml.RY(fixed_pars[2], wires=[1])
+    qml.RZ(fixed_pars[4], wires=[1])
+    qml.CNOT(wires=[0, 1])
+    qml.RX(fixed_pars[0], wires=[0])
+    qml.RY(fixed_pars[1], wires=[0])
+    qml.RZ(fixed_pars[3], wires=[0])
+    qml.RX(fixed_pars[1], wires=[1])
+    qml.RY(fixed_pars[2], wires=[1])
+    qml.RZ(fixed_pars[0], wires=[1])
+    qml.CNOT(wires=[0, 1])
+    qml.RX(params0, wires=[0])
+    qml.RX(params0, wires=[1])
+    qml.CNOT(wires=[0, 1])
+    qml.RY(fixed_pars[4], wires=[1])
+    qml.RY(params1, wires=[0])
+    qml.RY(params1, wires=[1])
+    qml.CNOT(wires=[0, 1])
+    qml.RX(fixed_pars[2], wires=[1])
+
+
+def fubini_ansatz10(weights, wires=None):
+    qml.templates.BasicEntanglerLayers(weights, wires=[0, 1])
+
+
+B = np.array(
+    [
+        [
+            [0.73, 0.49, 0.04],
+            [0.29, 0.45, 0.59],
+            [0.64, 0.06, 0.26],
+        ],
+        [
+            [0.93, 0.14, 0.46],
+            [0.31, 0.83, 0.79],
+            [0.25, 0.40, 0.16],
+        ],
+    ],
+    requires_grad=True,
+)
+fubini_ansatze_tape = [fubini_ansatz0, fubini_ansatz1, fubini_ansatz8]
+fubini_params_tape = [
+    (np.array([0.3434, -0.7245345], requires_grad=True),),
+    (B,),
+    (np.array(-0.1735, requires_grad=True),),
+]
+
+fubini_ansatze = [
+    fubini_ansatz0,
+    fubini_ansatz1,
+    fubini_ansatz2,
+    fubini_ansatz3,
+    fubini_ansatz4,
+    fubini_ansatz5,
+    fubini_ansatz6,
+    fubini_ansatz7,
+    fubini_ansatz8,
+    fubini_ansatz9,
+    fubini_ansatz10,
+]
+
+fubini_params = [
+    (np.array([0.3434, -0.7245345], requires_grad=True),),
+    (B,),
+    (np.array([-0.1111, -0.2222], requires_grad=True),),
+    (np.array([-0.1111, -0.2222, 0.4554], requires_grad=True),),
+    (
+        np.array(-0.1735, requires_grad=True),
+        np.array([-0.1735, -0.2846, -0.2846], requires_grad=True),
+    ),
+    (np.array([-0.1735, -0.2846], requires_grad=True),),
+    (np.array([-0.1735, -0.2846], requires_grad=True),),
+    (
+        np.array([-0.1735, -0.2846], requires_grad=True),
+        np.array([0.9812, -0.1492], requires_grad=True),
+    ),
+    (np.array(-0.1735, requires_grad=True),),
+    (np.array([-0.1111, 0.3333], requires_grad=True),),
+    (np.array([[0.21, 9.29], [-0.2, 0.12], [0.3, -2.1]], requires_grad=True),),
+]
+
+
+def autodiff_metric_tensor(ansatz, num_wires):
+    """Compute the metric tensor by full state vector
+    differentiation via autograd."""
+    dev = qml.device("default.qubit", wires=num_wires)
+
+    @qml.qnode(dev)
+    def qnode(*params):
+        ansatz(*params, wires=dev.wires)
+        return qml.state()
+
+    def mt(*params):
+        state = qnode(*params)
+        rqnode = lambda *params: np.real(qnode(*params))
+        iqnode = lambda *params: np.imag(qnode(*params))
+        rjac = qml.jacobian(rqnode)(*params)
+        ijac = qml.jacobian(iqnode)(*params)
+
+        if isinstance(rjac, tuple):
+            out = []
+            for rc, ic in zip(rjac, ijac):
+                c = rc + 1j * ic
+                psidpsi = np.tensordot(np.conj(state), c, axes=([0], [0]))
+                out.append(
+                    np.real(
+                        np.tensordot(np.conj(c), c, axes=([0], [0]))
+                        - np.tensordot(np.conj(psidpsi), psidpsi, axes=0)
+                    )
+                )
+            return tuple(out)
+
+        jac = rjac + 1j * ijac
+        psidpsi = np.tensordot(np.conj(state), jac, axes=([0], [0]))
+        return np.real(
+            np.tensordot(np.conj(jac), jac, axes=([0], [0]))
+            - np.tensordot(np.conj(psidpsi), psidpsi, axes=0)
+        )
+
+    return mt
+
+
+@pytest.mark.parametrize("ansatz, params", list(zip(fubini_ansatze_tape, fubini_params_tape)))
+class TestAdjointMetricTensorTape:
+    """Test the adjoint method for the metric tensor when calling it directly on
+    a tape.
+    """
+
+    num_wires = 3
+
+    @pytest.mark.autograd
+    def test_correct_output_tape_autograd(self, ansatz, params):
+        """Test that the output is correct when using Autograd and
+        calling the adjoint metric tensor directly on a tape."""
+        expected = autodiff_metric_tensor(ansatz, self.num_wires)(*params)
+        dev = qml.device("default.qubit.autograd", wires=self.num_wires)
+
+        @qml.qnode(dev, interface="autograd")
+        def circuit(*params):
+            """Circuit with dummy output to create a QNode."""
+            ansatz(*params, dev.wires)
+            return qml.expval(qml.PauliZ(0))
+
+        circuit(*params)
+        mt = qml.adjoint_metric_tensor(circuit.qtape, dev)
+        expected = qml.math.reshape(expected, qml.math.shape(mt))
+        assert qml.math.allclose(mt, expected)
+
+        mt = qml.adjoint_metric_tensor(circuit, hybrid=False)(*params)
+        assert qml.math.allclose(mt, expected)
+
+    @pytest.mark.jax
+    @pytest.mark.skip("JAX does not support forward pass executiong of the metric tensor.")
+    def test_correct_output_tape_jax(self, ansatz, params):
+        """Test that the output is correct when using JAX and
+        calling the adjoint metric tensor directly on a tape."""
+
+        import jax
+        from jax.config import config
+
+        config.update("jax_enable_x64", True)
+
+        expected = autodiff_metric_tensor(ansatz, self.num_wires)(*params)
+        j_params = tuple(jax.numpy.array(p) for p in params)
+        dev = qml.device("default.qubit.jax", wires=self.num_wires)
+
+        @qml.qnode(dev, interface="jax")
+        def circuit(*params):
+            """Circuit with dummy output to create a QNode."""
+            ansatz(*params, dev.wires)
+            return qml.expval(qml.PauliZ(0))
+
+        circuit(*j_params)
+        mt = qml.adjoint_metric_tensor(circuit.qtape, dev)
+        expected = qml.math.reshape(expected, qml.math.shape(mt))
+        assert qml.math.allclose(mt, expected)
+
+        mt = qml.adjoint_metric_tensor(circuit, hybrid=False)(*j_params)
+        assert qml.math.allclose(mt, expected)
+
+    @pytest.mark.torch
+    def test_correct_output_tape_torch(self, ansatz, params):
+        """Test that the output is correct when using Torch and
+        calling the adjoint metric tensor directly on a tape."""
+
+        import torch
+
+        expected = autodiff_metric_tensor(ansatz, self.num_wires)(*params)
+        t_params = tuple(torch.tensor(p, requires_grad=True) for p in params)
+        dev = qml.device("default.qubit.torch", wires=self.num_wires)
+
+        @qml.qnode(dev, interface="torch")
+        def circuit(*params):
+            """Circuit with dummy output to create a QNode."""
+            ansatz(*params, dev.wires)
+            return qml.expval(qml.PauliZ(0))
+
+        circuit(*t_params)
+        mt = qml.adjoint_metric_tensor(circuit.qtape, dev)
+        expected = qml.math.reshape(expected, qml.math.shape(mt))
+        assert qml.math.allclose(mt.detach().numpy(), expected)
+
+        mt = qml.adjoint_metric_tensor(circuit, hybrid=False)(*t_params)
+        assert qml.math.allclose(mt, expected)
+
+    @pytest.mark.tf
+    def test_correct_output_tape_tf(self, ansatz, params):
+        """Test that the output is correct when using TensorFlow and
+        calling the adjoint metric tensor directly on a tape."""
+
+        import tensorflow as tf
+
+        expected = autodiff_metric_tensor(ansatz, self.num_wires)(*params)
+        t_params = tuple(tf.Variable(p) for p in params)
+        dev = qml.device("default.qubit.tf", wires=self.num_wires)
+
+        @qml.qnode(dev, interface="tf")
+        def circuit(*params):
+            """Circuit with dummy output to create a QNode."""
+            ansatz(*params, dev.wires)
+            return qml.expval(qml.PauliZ(0))
+
+        with tf.GradientTape() as t:
+            circuit(*t_params)
+            mt = qml.adjoint_metric_tensor(circuit.qtape, dev)
+
+        expected = qml.math.reshape(expected, qml.math.shape(mt))
+        assert qml.math.allclose(mt, expected)
+
+        with tf.GradientTape() as t:
+            mt = qml.adjoint_metric_tensor(circuit, hybrid=False)(*t_params)
+        assert qml.math.allclose(mt, expected)
+
+
+class TestAdjointMetricTensorQNode:
+    """Test the adjoint method for the metric tensor when calling it on
+    a QNode.
+    """
+
+    num_wires = 3
+
+    @pytest.mark.autograd
+    @pytest.mark.parametrize("ansatz, params", list(zip(fubini_ansatze, fubini_params)))
+    def test_correct_output_qnode_autograd(self, ansatz, params):
+        """Test that the output is correct when using Autograd and
+        calling the adjoint metric tensor on a QNode."""
+        expected = autodiff_metric_tensor(ansatz, self.num_wires)(*params)
+        dev = qml.device("default.qubit", wires=self.num_wires)
+
+        @qml.qnode(dev, interface="autograd")
+        def circuit(*params):
+            """Circuit with dummy output to create a QNode."""
+            ansatz(*params, dev.wires)
+            return qml.expval(qml.PauliZ(0))
+
+        mt = qml.adjoint_metric_tensor(circuit)(*params)
+
+        if isinstance(mt, tuple):
+            assert all(qml.math.allclose(_mt, _exp) for _mt, _exp in zip(mt, expected))
+        else:
+            assert qml.math.allclose(mt, expected)
+
+    @pytest.mark.jax
+    @pytest.mark.skip("JAX does not support forward pass executiong of the metric tensor.")
+    @pytest.mark.parametrize("ansatz, params", list(zip(fubini_ansatze, fubini_params)))
+    def test_correct_output_qnode_jax(self, ansatz, params):
+        """Test that the output is correct when using JAX and
+        calling the adjoint metric tensor on a QNode."""
+
+        import jax
+        from jax.config import config
+
+        config.update("jax_enable_x64", True)
+
+        expected = autodiff_metric_tensor(ansatz, self.num_wires)(*params)
+        j_params = tuple(jax.numpy.array(p) for p in params)
+        dev = qml.device("default.qubit", wires=self.num_wires)
+
+        @qml.qnode(dev, interface="jax")
+        def circuit(*params):
+            """Circuit with dummy output to create a QNode."""
+            ansatz(*params, dev.wires)
+            return qml.expval(qml.PauliZ(0))
+
+        mt = qml.adjoint_metric_tensor(circuit)(*j_params)
+
+        if isinstance(mt, tuple):
+            assert all(qml.math.allclose(_mt, _exp) for _mt, _exp in zip(mt, expected))
+        else:
+            assert qml.math.allclose(mt, expected)
+
+    @pytest.mark.torch
+    @pytest.mark.parametrize("ansatz, params", list(zip(fubini_ansatze, fubini_params)))
+    def test_correct_output_qnode_torch(self, ansatz, params):
+        """Test that the output is correct when using Torch and
+        calling the adjoint metric tensor on a QNode."""
+
+        import torch
+
+        expected = autodiff_metric_tensor(ansatz, self.num_wires)(*params)
+        t_params = tuple(torch.tensor(p, requires_grad=True, dtype=torch.float64) for p in params)
+        dev = qml.device("default.qubit", wires=self.num_wires)
+
+        @qml.qnode(dev, interface="torch")
+        def circuit(*params):
+            """Circuit with dummy output to create a QNode."""
+            ansatz(*params, dev.wires)
+            return qml.expval(qml.PauliZ(0))
+
+        mt = qml.adjoint_metric_tensor(circuit)(*t_params)
+
+        if isinstance(mt, tuple):
+            assert all(qml.math.allclose(_mt, _exp) for _mt, _exp in zip(mt, expected))
+        else:
+            assert qml.math.allclose(mt, expected)
+
+    @pytest.mark.tf
+    @pytest.mark.parametrize("ansatz, params", list(zip(fubini_ansatze, fubini_params)))
+    def test_correct_output_qnode_tf(self, ansatz, params):
+        """Test that the output is correct when using TensorFlow and
+        calling the adjoint metric tensor on a QNode."""
+
+        import tensorflow as tf
+
+        expected = autodiff_metric_tensor(ansatz, self.num_wires)(*params)
+        t_params = tuple(tf.Variable(p, dtype=tf.float64) for p in params)
+        dev = qml.device("default.qubit", wires=self.num_wires)
+
+        @qml.qnode(dev, interface="tf")
+        def circuit(*params):
+            """Circuit with dummy output to create a QNode."""
+            ansatz(*params, dev.wires)
+            return qml.expval(qml.PauliZ(0))
+
+        with tf.GradientTape() as t:
+            mt = qml.adjoint_metric_tensor(circuit)(*t_params)
+
+        if isinstance(mt, tuple):
+            assert all(qml.math.allclose(_mt, _exp) for _mt, _exp in zip(mt, expected))
+        else:
+            assert qml.math.allclose(mt, expected)
+
+    @pytest.mark.autograd
+    def test_autograd_with_other_device(self):
+        """Test passing an extra device to the QNode wrapper."""
+        ansatz = fubini_ansatz2
+        params = fubini_params[2]
+
+        exp_fn = autodiff_metric_tensor(ansatz, self.num_wires)
+        expected = qml.jacobian(exp_fn)(*params)
+        dev = qml.device("default.qubit", wires=self.num_wires)
+        dev2 = qml.device("default.qubit.autograd", wires=self.num_wires)
+
+        @qml.qnode(dev, interface="autograd")
+        def circuit(*params):
+            """Circuit with dummy output to create a QNode."""
+            ansatz(*params, dev.wires)
+            return qml.expval(qml.PauliZ(0))
+
+        mt = qml.jacobian(qml.adjoint_metric_tensor(circuit, device=dev2))(*params)
+
+        if isinstance(mt, tuple):
+            assert all(qml.math.allclose(_mt, _exp) for _mt, _exp in zip(mt, expected))
+        else:
+            assert qml.math.allclose(mt, expected)
+
+
+diff_fubini_ansatze = [
+    fubini_ansatz0,
+    fubini_ansatz2,
+    fubini_ansatz10,
+]
+
+diff_fubini_params = [
+    fubini_params[0],
+    fubini_params[2],
+    fubini_params[10],
+]
+
+
+@pytest.mark.parametrize("ansatz, params", list(zip(diff_fubini_ansatze, diff_fubini_params)))
+class TestAdjointMetricTensorDifferentiability:
+    """Test the differentiability of the adjoint method for the metric
+    tensor when calling it on a QNode.
+    """
+
+    num_wires = 3
+
+    @pytest.mark.autograd
+    def test_autograd(self, ansatz, params):
+        """Test that the derivative is correct when using Autograd and
+        calling the adjoint metric tensor on a QNode."""
+        exp_fn = autodiff_metric_tensor(ansatz, self.num_wires)
+        expected = qml.jacobian(exp_fn)(*params)
+        dev = qml.device("default.qubit", wires=self.num_wires)
+
+        @qml.qnode(dev, interface="autograd")
+        def circuit(*params):
+            """Circuit with dummy output to create a QNode."""
+            ansatz(*params, dev.wires)
+            return qml.expval(qml.PauliZ(0))
+
+        mt_jac = qml.jacobian(qml.adjoint_metric_tensor(circuit))(*params)
+
+        if isinstance(mt_jac, tuple):
+            assert all(qml.math.allclose(_mt, _exp) for _mt, _exp in zip(mt_jac, expected))
+        else:
+            assert qml.math.allclose(mt_jac, expected)
+
+    @pytest.mark.jax
+    def test_correct_output_qnode_jax(self, ansatz, params):
+        """Test that the derivative is correct when using JAX and
+        calling the adjoint metric tensor on a QNode."""
+
+        import jax
+        from jax.config import config
+
+        config.update("jax_enable_x64", True)
+
+        expected = qml.jacobian(autodiff_metric_tensor(ansatz, self.num_wires))(*params)
+        j_params = tuple(jax.numpy.array(p) for p in params)
+        dev = qml.device("default.qubit", wires=self.num_wires)
+
+        @qml.qnode(dev, interface="jax")
+        def circuit(*params):
+            """Circuit with dummy output to create a QNode."""
+            ansatz(*params, dev.wires)
+            return qml.expval(qml.PauliZ(0))
+
+        mt_fn = qml.adjoint_metric_tensor(circuit, hybrid=True)
+        argnums = list(range(len(params)))
+        mt_jac = jax.jacobian(mt_fn, argnums=argnums)(*j_params)
+
+        if isinstance(mt_jac, tuple):
+            if not isinstance(expected, tuple) and len(mt_jac) == 1:
+                expected = (expected,)
+            assert all(qml.math.allclose(_mt, _exp) for _mt, _exp in zip(mt_jac, expected))
+        else:
+            assert qml.math.allclose(mt_jac, expected)
+
+    @pytest.mark.torch
+    def test_correct_output_qnode_torch(self, ansatz, params):
+        """Test that the derivative is correct when using Torch and
+        calling the adjoint metric tensor on a QNode."""
+
+        import torch
+
+        expected = qml.jacobian(autodiff_metric_tensor(ansatz, self.num_wires))(*params)
+        t_params = tuple(torch.tensor(p, requires_grad=True, dtype=torch.float64) for p in params)
+        dev = qml.device("default.qubit", wires=self.num_wires)
+
+        @qml.qnode(dev, interface="torch")
+        def circuit(*params):
+            """Circuit with dummy output to create a QNode."""
+            ansatz(*params, dev.wires)
+            return qml.expval(qml.PauliZ(0))
+
+        mt_fn = qml.adjoint_metric_tensor(circuit)
+        mt_jac = torch.autograd.functional.jacobian(mt_fn, *t_params)
+
+        if isinstance(mt_jac, tuple):
+            assert all(qml.math.allclose(_mt, _exp) for _mt, _exp in zip(mt_jac, expected))
+        else:
+            assert qml.math.allclose(mt_jac, expected)
+
+    @pytest.mark.tf
+    def test_correct_output_qnode_tf(self, ansatz, params):
+        """Test that the derivative is correct when using TensorFlow and
+        calling the adjoint metric tensor on a QNode."""
+
+        import tensorflow as tf
+
+        expected = qml.jacobian(autodiff_metric_tensor(ansatz, self.num_wires))(*params)
+        t_params = tuple(tf.Variable(p, dtype=tf.float64) for p in params)
+        dev = qml.device("default.qubit", wires=self.num_wires)
+
+        @qml.qnode(dev, interface="tf")
+        def circuit(*params):
+            """Circuit with dummy output to create a QNode."""
+            ansatz(*params, dev.wires)
+            return qml.expval(qml.PauliZ(0))
+
+        with tf.GradientTape() as t:
+            mt = qml.adjoint_metric_tensor(circuit)(*t_params)
+
+        mt_jac = t.jacobian(mt, t_params)
+        if isinstance(mt_jac, tuple):
+            if not isinstance(expected, tuple) and len(mt_jac) == 1:
+                expected = (expected,)
+            assert all(qml.math.allclose(_mt, _exp) for _mt, _exp in zip(mt_jac, expected))
+        else:
+            assert qml.math.allclose(mt_jac, expected)
+
+
+class TestErrors:
+    """Test that errors are raised correctly."""
+
+    def test_error_wrong_object_passed(self):
+        """Test that an error is raised if neither a tape nor a QNode is passed."""
+
+        def ansatz(x, y):
+            qml.RX(x, wires=0)
+            qml.RY(y, wires=1)
+
+        dev = qml.device("default.qubit", wires=2)
+
+        with pytest.raises(qml.QuantumFunctionError, match="The passed object is not a "):
+            qml.adjoint_metric_tensor(ansatz, device=dev)
+
+    def test_error_finite_shots(self):
+        """Test that an error is raised if the device has a finite number of shots set."""
+        with qml.tape.QuantumTape() as tape:
+            qml.RX(0.2, wires=0)
+            qml.RY(1.9, wires=1)
+        dev = qml.device("default.qubit", wires=2, shots=1)
+
+        with pytest.raises(ValueError, match="The adjoint method for the metric tensor"):
+            qml.adjoint_metric_tensor(tape, device=dev)
+
+    def test_warning_multiple_devices(self):
+        """Test that a warning is issued if an ExpvalCost with multiple
+        devices is passed."""
+        dev1 = qml.device("default.qubit", wires=2)
+        dev2 = qml.device("default.qubit", wires=1)
+        H = qml.Hamiltonian([0.2, 0.9], [qml.PauliZ(0), qml.PauliY(0)])
+
+        def ansatz(x, wires):
+            qml.RX(x, wires=0)
+
+        with pytest.warns(UserWarning, match="is deprecated,"):
+            cost = qml.ExpvalCost(ansatz, H, [dev1, dev2])
+        with pytest.warns(UserWarning, match="ExpvalCost was instantiated"):
+            mt = qml.adjoint_metric_tensor(cost)
diff --git a/tests/returntypes/transforms/test_fisher_new.py b/tests/returntypes/transforms/test_fisher_new.py
new file mode 100644
index 000000000..54445d129
--- /dev/null
+++ b/tests/returntypes/transforms/test_fisher_new.py
@@ -0,0 +1,570 @@
+# Copyright 2022 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
+
+#     http://www.apache.org/licenses/LICENSE-2.0
+
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+"""
+Tests for the classical fisher information matrix in the pennylane.qinfo
+"""
+# pylint: disable=no-self-use, import-outside-toplevel, no-member, import-error, too-few-public-methods, bad-continuation
+import pytest
+
+import pennylane as qml
+import pennylane.numpy as pnp
+import numpy as np
+
+
+from pennylane.qinfo import classical_fisher, quantum_fisher
+from pennylane.qinfo.transforms import _make_probs, _compute_cfim
+
+
+class TestMakeProbs:
+    """Class to test the private function _make_probs."""
+
+    def test_make_probs_makes_probs(self):
+        """Testing the correctness of _make_probs."""
+        dev = qml.device("default.qubit", wires=3)
+
+        @qml.qnode(dev)
+        def qnode(x):
+            qml.RX(x, 0)
+            qml.CNOT(wires=[0, 1])
+            return qml.expval(qml.PauliX(0))
+
+        x = pnp.array([0.5])
+        new_qnode = _make_probs(qnode)
+        tape, _ = new_qnode.construct(x, {})
+        assert tape[0].observables[0].return_type == qml.measurements.Probability
+
+    def test_make_probs(self):
+        """Testing the private _make_probs transform"""
+        with qml.tape.QuantumTape() as tape:
+            qml.PauliX(0)
+            qml.PauliZ(1)
+            qml.PauliY(2)
+        new_tape, fn = _make_probs(tape)
+        assert len(new_tape) == 1
+        assert np.isclose(fn([1]), 1)
+
+
+class TestComputeclassicalFisher:
+    """Testing that given p and dp, _compute_cfim() computes the correct outputs"""
+
+    @pytest.mark.parametrize("n_params", np.arange(1, 10))
+    @pytest.mark.parametrize("n_wires", np.arange(1, 5))
+    def test_construction_of_compute_cfim(self, n_params, n_wires):
+        """Test that the construction in _compute_cfim is correct"""
+        dp = np.arange(2**n_wires * n_params, dtype=float).reshape(2**n_wires, n_params)
+        p = np.ones(2**n_wires)
+
+        res = _compute_cfim(p, dp)
+
+        assert np.allclose(res, res.T)
+        assert all(
+            [
+                res[i, j] == np.sum(dp[:, i] * dp[:, j] / p)
+                for i in range(n_params)
+                for j in range(n_params)
+            ]
+        )
+
+    @pytest.mark.parametrize("n_params", np.arange(1, 10))
+    @pytest.mark.parametrize("n_wires", np.arange(1, 5))
+    def test_compute_cfim_trivial_distribution(self, n_params, n_wires):
+        """Test that the classical fisher information matrix (CFIM) is
+        constructed correctly for the trivial distirbution p=(1,0,0,...)
+        and some dummy gradient dp"""
+        # implicitly also does the unit test for the correct output shape
+
+        dp = np.ones(2**n_wires * n_params, dtype=float).reshape(2**n_wires, n_params)
+        p = np.zeros(2**n_wires)
+        p[0] = 1
+        res = _compute_cfim(p, dp)
+        assert np.allclose(res, np.ones((n_params, n_params)))
+
+
+class TestIntegration:
+    """Integration test of classical and quantum fisher information matrices"""
+
+    @pytest.mark.parametrize("n_wires", np.arange(1, 5))
+    @pytest.mark.parametrize("n_params", np.arange(1, 5))
+    def test_different_sizes(self, n_wires, n_params):
+        """Testing that for any number of wires and parameters, the correct size and values are computed"""
+        dev = qml.device("default.qubit", wires=n_wires)
+
+        @qml.qnode(dev, interface="autograd")
+        def circ(params):
+            for i in range(n_wires):
+                qml.Hadamard(wires=i)
+
+            for x in params:
+                for j in range(n_wires):
+                    qml.RX(x, wires=j)
+                    qml.RY(x, wires=j)
+                    qml.RZ(x, wires=j)
+
+            return qml.probs(wires=range(n_wires))
+
+        params = pnp.zeros(n_params, requires_grad=True)
+        res = classical_fisher(circ)(params)
+        assert np.allclose(res, n_wires * np.ones((n_params, n_params)))
+
+    def test_hardware_compatibility_classical_fisher(self):
+        """Testing that classical_fisher can be computed with finite shots"""
+        n_wires = 3
+        n_params = 3
+
+        dev = qml.device("default.qubit", wires=n_wires, shots=10000)
+
+        @qml.qnode(dev)
+        def circ(params):
+            for i in range(n_wires):
+                qml.Hadamard(wires=i)
+
+            for x in params:
+                for j in range(n_wires):
+                    qml.RX(x, wires=j)
+                    qml.RY(x, wires=j)
+                    qml.RZ(x, wires=j)
+
+            return qml.probs(wires=range(n_wires))
+
+        params = pnp.zeros(n_params, requires_grad=True)
+        res = qml.qinfo.classical_fisher(circ)(params)
+        assert np.allclose(res, n_wires * np.ones((n_params, n_params)), atol=1)
+
+    def test_quantum_fisher_info(self):
+        """Integration test of quantum fisher information matrix CFIM. This is just calling ``qml.metric_tensor`` or ``qml.adjoint_metric_tensor`` and multiplying by a factor of 4"""
+
+        n_wires = 2
+
+        dev = qml.device("default.qubit", wires=n_wires)
+        dev_hard = qml.device("default.qubit", wires=n_wires + 1, shots=1000)
+
+        def qfunc(params):
+            qml.RX(params[0], wires=0)
+            qml.RX(params[1], wires=0)
+            qml.CNOT(wires=(0, 1))
+            return qml.state()
+
+        params = pnp.random.random(2)
+
+        circ_hard = qml.QNode(qfunc, dev_hard)
+        QFIM_hard = quantum_fisher(circ_hard)(params)
+        QFIM1_hard = 4.0 * qml.metric_tensor(circ_hard)(params)
+
+        circ = qml.QNode(qfunc, dev)
+        QFIM = quantum_fisher(circ)(params)
+        QFIM1 = 4.0 * qml.adjoint_metric_tensor(circ)(params)
+        assert np.allclose(QFIM, QFIM1)
+        assert np.allclose(QFIM_hard, QFIM1_hard, atol=1e-1)
+
+
+class TestInterfacesClassicalFisher:
+    """Integration tests for the classical fisher information matrix CFIM"""
+
+    @pytest.mark.autograd
+    @pytest.mark.parametrize("n_wires", np.arange(1, 5))
+    def test_cfim_allnonzero_autograd(self, n_wires):
+        """Integration test of classical_fisher() with autograd for examples where all probabilities are all nonzero"""
+
+        dev = qml.device("default.qubit", wires=n_wires)
+
+        @qml.qnode(dev, interface="autograd")
+        def circ(params):
+            for i in range(n_wires):
+                qml.RX(params[0], wires=i)
+            for i in range(n_wires):
+                qml.RY(params[1], wires=i)
+            return qml.probs(wires=range(n_wires))
+
+        params = np.pi / 4 * pnp.ones(2, requires_grad=True)
+        cfim = classical_fisher(circ)(params)
+        assert np.allclose(cfim, (n_wires / 3.0) * np.ones((2, 2)))
+
+    @pytest.mark.autograd
+    @pytest.mark.parametrize("n_wires", np.arange(2, 5))
+    def test_cfim_contains_zeros_autograd(self, n_wires):
+        """Integration test of classical_fisher() with autograd for examples that have 0s in the probabilities and non-zero gradient"""
+        dev = qml.device("default.qubit", wires=n_wires)
+
+        @qml.qnode(dev, interface="autograd")
+        def circ(params):
+            qml.RZ(params[0], wires=0)
+            qml.RX(params[1], wires=1)
+            qml.RX(params[0], wires=1)
+            return qml.probs(wires=range(n_wires))
+
+        params = np.pi / 4 * pnp.ones(2, requires_grad=True)
+        cfim = classical_fisher(circ)(params)
+        assert np.allclose(cfim, np.ones((2, 2)))
+
+    @pytest.mark.jax
+    @pytest.mark.parametrize("n_wires", np.arange(1, 5))
+    def test_cfim_allnonzero_jax(self, n_wires):
+        """Integration test of classical_fisher() with jax for examples where all probabilities are all nonzero"""
+        import jax.numpy as jnp
+
+        dev = qml.device("default.qubit", wires=n_wires)
+
+        @qml.qnode(dev, interface="jax")
+        def circ(params):
+            for i in range(n_wires):
+                qml.RX(params[0], wires=i)
+            for i in range(n_wires):
+                qml.RY(params[1], wires=i)
+            return qml.probs(wires=range(n_wires))
+
+        params = np.pi / 4 * jnp.ones(2)
+        cfim = classical_fisher(circ)(params)
+        assert np.allclose(cfim, (n_wires / 3.0) * np.ones((2, 2)))
+
+    @pytest.mark.jax
+    @pytest.mark.parametrize("n_wires", np.arange(2, 5))
+    def test_cfim_contains_zeros_jax(self, n_wires):
+        """Integration test of classical_fisher() with jax for examples that have 0s in the probabilities and non-zero gradient"""
+        import jax.numpy as jnp
+
+        dev = qml.device("default.qubit", wires=n_wires)
+
+        @qml.qnode(dev, interface="jax")
+        def circ(params):
+            qml.RZ(params[0], wires=0)
+            qml.RX(params[1], wires=1)
+            qml.RX(params[0], wires=1)
+            return qml.probs(wires=range(n_wires))
+
+        params = np.pi / 4 * jnp.ones(2)
+        cfim = classical_fisher(circ)(params)
+        assert np.allclose(cfim, np.ones((2, 2)))
+
+    @pytest.mark.jax
+    def test_cfim_multiple_args_jax(self):
+        """Testing multiple args to be differentiated using jax"""
+        import jax.numpy as jnp
+
+        n_wires = 3
+
+        dev = qml.device("default.qubit", wires=n_wires)
+
+        @qml.qnode(dev, interface="jax")
+        def circ(x, y, z):
+            for xi in x:
+                qml.RX(xi, wires=0)
+                qml.RX(xi, wires=1)
+            for yi in y:
+                qml.RY(yi, wires=0)
+                qml.RY(yi, wires=1)
+            for zi in z:
+                qml.RZ(zi, wires=0)
+                qml.RZ(zi, wires=1)
+            return qml.probs(wires=range(n_wires))
+
+        x = jnp.pi / 8 * jnp.ones(2)
+        y = jnp.pi / 8 * jnp.ones(10)
+        z = jnp.ones(1)
+        cfim = qml.qinfo.classical_fisher(circ, argnums=(0, 1, 2))(x, y, z)
+        assert qml.math.allclose(cfim[0], 2.0 / 3.0 * np.ones((2, 2)))
+        assert qml.math.allclose(cfim[1], 2.0 / 3.0 * np.ones((10, 10)))
+        assert qml.math.allclose(cfim[2], np.zeros((1, 1)))
+
+    @pytest.mark.torch
+    @pytest.mark.parametrize("n_wires", np.arange(1, 5))
+    def test_cfim_allnonzero_torch(self, n_wires):
+        """Integration test of classical_fisher() with torch for examples where all probabilities are all nonzero"""
+        import torch
+
+        dev = qml.device("default.qubit", wires=n_wires)
+
+        @qml.qnode(dev, interface="torch")
+        def circ(params):
+            for i in range(n_wires):
+                qml.RX(params[0], wires=i)
+            for i in range(n_wires):
+                qml.RY(params[1], wires=i)
+            return qml.probs(wires=range(n_wires))
+
+        params = np.pi / 4 * torch.tensor([1.0, 1.0], requires_grad=True)
+        cfim = classical_fisher(circ)(params)
+        assert np.allclose(cfim.detach().numpy(), (n_wires / 3.0) * np.ones((2, 2)))
+
+    @pytest.mark.torch
+    @pytest.mark.parametrize("n_wires", np.arange(2, 5))
+    def test_cfim_contains_zeros_torch(self, n_wires):
+        """Integration test of classical_fisher() with torch for examples that have 0s in the probabilities and non-zero gradient"""
+        import torch
+
+        dev = qml.device("default.qubit", wires=n_wires)
+
+        @qml.qnode(dev, interface="torch")
+        def circ(params):
+            qml.RZ(params[0], wires=0)
+            qml.RX(params[1], wires=1)
+            qml.RX(params[0], wires=1)
+            return qml.probs(wires=range(n_wires))
+
+        params = np.pi / 4 * torch.tensor([1.0, 1.0], requires_grad=True)
+        cfim = classical_fisher(circ)(params)
+        assert np.allclose(cfim.detach().numpy(), np.ones((2, 2)))
+
+    @pytest.mark.torch
+    def test_cfim_multiple_args_torch(self):
+        """Testing multiple args to be differentiated using torch"""
+        import torch
+
+        n_wires = 3
+
+        dev = qml.device("default.qubit", wires=n_wires)
+
+        @qml.qnode(dev, interface="torch")
+        def circ(x, y, z):
+            for xi in x:
+                qml.RX(xi, wires=0)
+                qml.RX(xi, wires=1)
+            for yi in y:
+                qml.RY(yi, wires=0)
+                qml.RY(yi, wires=1)
+            for zi in z:
+                qml.RZ(zi, wires=0)
+                qml.RZ(zi, wires=1)
+            return qml.probs(wires=range(n_wires))
+
+        x = np.pi / 8 * torch.ones(2, requires_grad=True)
+        y = np.pi / 8 * torch.ones(10, requires_grad=True)
+        z = torch.ones(1, requires_grad=True)
+        cfim = qml.qinfo.classical_fisher(circ)(x, y, z)
+        assert np.allclose(cfim[0].detach().numpy(), 2.0 / 3.0 * np.ones((2, 2)))
+        assert np.allclose(cfim[1].detach().numpy(), 2.0 / 3.0 * np.ones((10, 10)))
+        assert np.allclose(cfim[2].detach().numpy(), np.zeros((1, 1)))
+
+    @pytest.mark.tf
+    @pytest.mark.parametrize("n_wires", np.arange(1, 5))
+    def test_cfim_allnonzero_tf(self, n_wires):
+        """Integration test of classical_fisher() with tf for examples where all probabilities are all nonzero"""
+        import tensorflow as tf
+
+        dev = qml.device("default.qubit", wires=n_wires)
+
+        @qml.qnode(dev, interface="tf")
+        def circ(params):
+            for i in range(n_wires):
+                qml.RX(params[0], wires=i)
+            for i in range(n_wires):
+                qml.RY(params[1], wires=i)
+            return qml.probs(wires=range(n_wires))
+
+        params = tf.Variable([np.pi / 4, np.pi / 4])
+        cfim = classical_fisher(circ)(params)
+        assert np.allclose(cfim, (n_wires / 3.0) * np.ones((2, 2)))
+
+    @pytest.mark.tf
+    @pytest.mark.parametrize("n_wires", np.arange(2, 5))
+    def test_cfim_contains_zeros_tf(self, n_wires):
+        """Integration test of classical_fisher() with tf for examples that have 0s in the probabilities and non-zero gradient"""
+        import tensorflow as tf
+
+        dev = qml.device("default.qubit", wires=n_wires)
+
+        @qml.qnode(dev, interface="tf")
+        def circ(params):
+            qml.RZ(params[0], wires=0)
+            qml.RX(params[1], wires=1)
+            qml.RX(params[0], wires=1)
+            return qml.probs(wires=range(n_wires))
+
+        params = tf.Variable([np.pi / 4, np.pi / 4])
+        cfim = classical_fisher(circ)(params)
+        assert np.allclose(cfim, np.ones((2, 2)))
+
+    @pytest.mark.tf
+    def test_cfim_multiple_args_tf(self):
+        """Testing multiple args to be differentiated using tf"""
+        import tensorflow as tf
+
+        n_wires = 3
+
+        dev = qml.device("default.qubit", wires=n_wires)
+
+        @qml.qnode(dev, interface="tf")
+        def circ(x, y, z):
+            for xi in tf.unstack(x):
+                qml.RX(xi, wires=0)
+                qml.RX(xi, wires=1)
+            for yi in tf.unstack(y):
+                qml.RY(yi, wires=0)
+                qml.RY(yi, wires=1)
+            for zi in tf.unstack(z):
+                qml.RZ(zi, wires=0)
+                qml.RZ(zi, wires=1)
+            return qml.probs(wires=range(n_wires))
+
+        x = tf.Variable(np.pi / 8 * np.ones(2), trainable=True)
+        y = tf.Variable(np.pi / 8 * np.ones(10), trainable=True)
+        z = tf.Variable([1.0], trainable=True)
+        cfim = qml.qinfo.classical_fisher(circ)(x, y, z)
+        assert np.allclose(cfim[0], 2.0 / 3.0 * np.ones((2, 2)))
+        assert np.allclose(cfim[1], 2.0 / 3.0 * np.ones((10, 10)))
+        assert np.allclose(cfim[2], np.zeros((1, 1)))
+
+
+class TestDiffCFIM:
+    """Testing differentiability of classical fisher info matrix (CFIM)"""
+
+    @pytest.mark.autograd
+    def test_diffability_autograd(self):
+        """Testing diffability with an analytic example for autograd. The CFIM of this single qubit is constant, so the gradient should be zero."""
+        dev = qml.device("default.qubit", wires=1)
+
+        # qml.disable_return()
+
+        @qml.qnode(dev, interface="autograd")
+        def circ(params):
+            qml.RY(params, wires=0)
+            return qml.probs(wires=range(1))
+
+        params = pnp.array(np.pi / 4, requires_grad=True)
+
+        assert np.allclose(qml.qinfo.classical_fisher(circ)(params), 1)
+
+        result = np.zeros((1, 1, 1), dtype="float64")
+        result_calc = qml.jacobian(qml.qinfo.classical_fisher(circ))(params)
+
+        assert np.allclose(result, result_calc, atol=1e-6)
+
+    @pytest.mark.jax
+    def test_diffability_jax(self):
+        """Testing diffability with an analytic example for jax. The CFIM of this single qubit is constant, so the gradient should be zero."""
+        import jax.numpy as jnp
+        import jax
+
+        dev = qml.device("default.qubit", wires=1)
+
+        @qml.qnode(dev, interface="jax")
+        def circ(params):
+            qml.RY(params, wires=0)
+            return qml.probs(wires=range(1))
+
+        # no matter what the input the CFIM here is always constant = 1
+        # so the derivative should be 0
+
+        params = jnp.array(np.pi / 4)
+
+        assert qml.math.allclose(qml.qinfo.classical_fisher(circ)(params), 1.0)
+
+        result = np.zeros((1, 1, 1), dtype="float64")
+        result_calc = jax.jacobian(qml.qinfo.classical_fisher(circ))(params)
+
+        assert np.allclose(result, result_calc, atol=1e-6)
+
+    @pytest.mark.tf
+    def test_diffability_tf(self):
+        """Testing diffability with an analytic example for tf. The CFIM of this single
+        qubit is constant, so the gradient should be zero."""
+        import tensorflow as tf
+
+        dev = qml.device("default.qubit", wires=1)
+
+        @qml.qnode(dev, interface="tf")
+        def circ(params):
+            qml.RY(params, wires=0)
+            return qml.probs(wires=range(1))
+
+        params = tf.Variable(np.pi / 4)
+
+        assert np.allclose(qml.qinfo.classical_fisher(circ)(params), 1)
+
+        result = np.zeros((1, 1, 1), dtype="float64")
+        with tf.GradientTape() as tape:
+            loss = qml.qinfo.classical_fisher(circ)(params)
+
+        result_calc = tape.jacobian(loss, params)
+
+        assert np.allclose(result, result_calc, atol=1e-6)
+
+    @pytest.mark.torch
+    def test_diffability_torch(self):
+        """Testing diffability with an analytic example for torch. The CFIM of this single qubit is constant, so the gradient should be zero."""
+        import torch
+
+        dev = qml.device("default.qubit", wires=1)
+
+        @qml.qnode(dev, interface="torch")
+        def circ(params):
+            qml.RY(params, wires=0)
+            return qml.probs(wires=range(1))
+
+        params = torch.tensor(np.pi / 4, requires_grad=True)
+
+        assert np.allclose(qml.qinfo.classical_fisher(circ)(params).detach().numpy(), 1)
+
+        result = np.zeros((1, 1, 1), dtype="float64")
+        result_calc = torch.autograd.functional.jacobian(qml.qinfo.classical_fisher(circ), params)
+
+        assert np.allclose(result, result_calc, atol=1e-6)
+
+    @pytest.mark.all_interfaces
+    def test_consistency(self):
+        """Testing that the derivative of the cfim is giving consistently the same results for all interfaces.
+        Currently failing as (jax and autograd) and (torch and tf) are giving two different results."""
+        import tensorflow as tf
+        import torch
+        import jax.numpy as jnp
+        import jax
+
+        dev = qml.device("default.qubit", wires=3)
+
+        def qfunc(weights):
+            qml.RX(weights[0], wires=0)
+            qml.RX(weights[1], wires=1)
+            qml.RX(weights[2], wires=2)
+            qml.CNOT(wires=[0, 1])
+            qml.RY(weights[4], wires=0)
+            qml.RY(weights[5], wires=1)
+            qml.RY(weights[6], wires=2)
+            qml.CNOT(wires=[1, 2])
+            qml.CRX(weights[7], wires=[0, 1])
+            qml.CRY(weights[8], wires=[1, 2])
+            qml.CRZ(weights[9], wires=[2, 0])
+            return qml.expval(qml.PauliZ(0) @ qml.PauliZ(1))
+
+        # Compute gradients of CFIM for different interfaces
+        circuit = qml.QNode(qfunc, dev, interface="torch")
+        weights = torch.tensor(
+            [0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0], requires_grad=True
+        )
+        grad_torch = torch.autograd.functional.jacobian(
+            qml.qinfo.classical_fisher(circuit), weights
+        )
+
+        circuit = qml.QNode(qfunc, dev, interface="autograd")
+        weights = pnp.array([0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0], requires_grad=True)
+        grad_autograd = qml.jacobian(qml.qinfo.classical_fisher(circuit))(weights)
+
+        circuit = qml.QNode(qfunc, dev, interface="jax")
+        weights = jnp.array([0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0])
+        grad_jax = jax.jacobian(qml.qinfo.classical_fisher(circuit))(weights)
+
+        circuit = qml.QNode(qfunc, dev, interface="tf")
+        weights = tf.Variable([0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0])
+        with tf.GradientTape() as tape:
+            loss = qml.qinfo.classical_fisher(circuit)(weights)
+        grad_tf = tape.jacobian(loss, weights)
+
+        # Evaluate and compare
+        grads = [grad_autograd, grad_tf, grad_torch, grad_jax]
+
+        is_same = np.zeros((4, 4))
+
+        for i, g1 in enumerate(grads):
+            for j, g2 in enumerate(grads):
+                is_same[i, j] = np.allclose(g1, g2, atol=1e-6)
+
+        assert np.allclose(is_same, np.ones((4, 4)))
diff --git a/tests/returntypes/transforms/test_metric_tensor_new.py b/tests/returntypes/transforms/test_metric_tensor_new.py
new file mode 100644
index 000000000..1587f9975
--- /dev/null
+++ b/tests/returntypes/transforms/test_metric_tensor_new.py
@@ -0,0 +1,1472 @@
+# 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
+
+#     http://www.apache.org/licenses/LICENSE-2.0
+
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+"""
+Unit tests for the metric tensor transform.
+"""
+import pytest
+from pennylane import numpy as np
+from scipy.linalg import block_diag
+
+import pennylane as qml
+from gate_data import Y, Z
+from pennylane.transforms.metric_tensor import _get_aux_wire
+
+
+class TestMetricTensor:
+    """Tests for metric tensor subcircuit construction and evaluation"""
+
+    @pytest.mark.parametrize("diff_method", ["parameter-shift", "backprop"])
+    def test_rot_decomposition(self, diff_method):
+        """Test that the rotation gate is correctly decomposed"""
+        dev = qml.device("default.qubit", wires=1)
+        params = np.array([1.0, 2.0, 3.0], requires_grad=True)
+
+        with qml.tape.QuantumTape() as circuit:
+            qml.Rot(params[0], params[1], params[2], wires=0)
+            qml.expval(qml.PauliX(0))
+
+        tapes, _ = qml.metric_tensor(circuit, approx="block-diag")
+        assert len(tapes) == 3
+
+        # first parameter subcircuit
+        assert len(tapes[0].operations) == 0
+
+        # Second parameter subcircuit
+        assert len(tapes[1].operations) == 4
+        assert isinstance(tapes[1].operations[0], qml.RZ)
+        assert tapes[1].operations[0].data == [1]
+        # PauliY decomp
+        assert isinstance(tapes[1].operations[1], qml.PauliZ)
+        assert isinstance(tapes[1].operations[2], qml.S)
+        assert isinstance(tapes[1].operations[3], qml.Hadamard)
+
+        # Third parameter subcircuit
+        assert len(tapes[2].operations) == 2
+        assert isinstance(tapes[2].operations[0], qml.RZ)
+        assert isinstance(tapes[2].operations[1], qml.RY)
+        assert tapes[2].operations[0].data == [1]
+        assert tapes[2].operations[1].data == [2]
+
+    @pytest.mark.parametrize("diff_method", ["parameter-shift", "backprop"])
+    def test_multirz_decomposition(self, diff_method):
+        """Test that the MultiRZ gate is correctly decomposed"""
+        dev = qml.device("default.qubit", wires=3)
+
+        @qml.qnode(dev, diff_method=diff_method)
+        def circuit(a, b):
+            qml.RX(a, wires=0)
+            qml.MultiRZ(b, wires=[0, 1, 2])
+            return qml.expval(qml.PauliX(0))
+
+        params = np.array([0.1, 0.2], requires_grad=True)
+        result = qml.metric_tensor(circuit, approx="block-diag")(*params)
+        assert isinstance(result, tuple) and len(result) == 2
+        assert qml.math.shape(result[0]) == ()
+        assert qml.math.shape(result[1]) == ()
+
+    @pytest.mark.parametrize("diff_method", ["parameter-shift", "backprop"])
+    def test_parameter_fan_out(self, diff_method):
+        """The metric tensor is with respect to the quantum circuit and ignores
+        classical processing if ``hybrid=False``. As a result, if there is
+        parameter fan-out, the returned metric tensor will be larger than
+        ``(len(args), len(args))`` if hybrid computation is deactivated.
+        """
+        dev = qml.device("default.qubit", wires=2)
+
+        def circuit(a):
+            qml.RX(a, wires=0)
+            qml.RX(a, wires=0)
+            return qml.expval(qml.PauliX(0))
+
+        circuit = qml.QNode(circuit, dev, diff_method=diff_method)
+        params = np.array([0.1], requires_grad=True)
+        result = qml.metric_tensor(circuit, hybrid=False, approx="block-diag")(*params)
+        assert result.shape == (2, 2)
+
+    def test_construct_subcircuit(self):
+        """Test correct subcircuits constructed"""
+        dev = qml.device("default.qubit", wires=2)
+
+        with qml.tape.QuantumTape() as tape:
+            qml.RX(np.array(1.0, requires_grad=True), wires=0)
+            qml.RY(np.array(1.0, requires_grad=True), wires=0)
+            qml.CNOT(wires=[0, 1])
+            qml.PhaseShift(np.array(1.0, requires_grad=True), wires=1)
+            return qml.expval(qml.PauliX(0)), qml.expval(qml.PauliX(1))
+
+        tapes, _ = qml.metric_tensor(tape, approx="block-diag")
+        assert len(tapes) == 3
+
+        # first parameter subcircuit
+        assert len(tapes[0].operations) == 1
+        assert isinstance(tapes[0].operations[0], qml.Hadamard)  # PauliX decomp
+
+        # second parameter subcircuit
+        assert len(tapes[1].operations) == 4
+        assert isinstance(tapes[1].operations[0], qml.RX)
+        # PauliY decomp
+        assert isinstance(tapes[1].operations[1], qml.PauliZ)
+        assert isinstance(tapes[1].operations[2], qml.S)
+        assert isinstance(tapes[1].operations[3], qml.Hadamard)
+
+        # third parameter subcircuit
+        assert len(tapes[2].operations) == 4
+        assert isinstance(tapes[2].operations[0], qml.RX)
+        assert isinstance(tapes[2].operations[1], qml.RY)
+        assert isinstance(tapes[2].operations[2], qml.CNOT)
+        # Phase shift generator
+        assert isinstance(tapes[2].operations[3], qml.QubitUnitary)
+
+    def test_construct_subcircuit_layers(self):
+        """Test correct subcircuits constructed
+        when a layer structure exists"""
+        dev = qml.device("default.qubit", wires=3)
+        params = np.ones([8])
+
+        with qml.tape.QuantumTape() as tape:
+            # section 1
+            qml.RX(params[0], wires=0)
+            # section 2
+            qml.RY(params[1], wires=0)
+            qml.CNOT(wires=[0, 1])
+            qml.CNOT(wires=[1, 2])
+            # section 3
+            qml.RX(params[2], wires=0)
+            qml.RY(params[3], wires=1)
+            qml.RZ(params[4], wires=2)
+            qml.CNOT(wires=[0, 1])
+            qml.CNOT(wires=[1, 2])
+            # section 4
+            qml.RX(params[5], wires=0)
+            qml.RY(params[6], wires=1)
+            qml.RZ(params[7], wires=2)
+            qml.CNOT(wires=[0, 1])
+            qml.CNOT(wires=[1, 2])
+            return qml.expval(qml.PauliX(0)), qml.expval(qml.PauliX(1)), qml.expval(qml.PauliX(2))
+
+        tapes, _ = qml.metric_tensor(tape, approx="block-diag")
+
+        # this circuit should split into 4 independent
+        # sections or layers when constructing subcircuits
+        assert len(tapes) == 4
+
+        # first layer subcircuit
+        assert len(tapes[0].operations) == 1
+        assert isinstance(tapes[0].operations[0], qml.Hadamard)  # PauliX decomp
+
+        # second layer subcircuit
+        assert len(tapes[1].operations) == 4
+        assert isinstance(tapes[1].operations[0], qml.RX)
+        # PauliY decomp
+        assert isinstance(tapes[1].operations[1], qml.PauliZ)
+        assert isinstance(tapes[1].operations[2], qml.S)
+        assert isinstance(tapes[1].operations[3], qml.Hadamard)
+
+        # # third layer subcircuit
+        assert len(tapes[2].operations) == 8
+        assert isinstance(tapes[2].operations[0], qml.RX)
+        assert isinstance(tapes[2].operations[1], qml.RY)
+        assert isinstance(tapes[2].operations[2], qml.CNOT)
+        assert isinstance(tapes[2].operations[3], qml.CNOT)
+        # PauliX decomp
+        assert isinstance(tapes[2].operations[4], qml.Hadamard)
+        # PauliY decomp
+        assert isinstance(tapes[2].operations[5], qml.PauliZ)
+        assert isinstance(tapes[2].operations[6], qml.S)
+        assert isinstance(tapes[2].operations[7], qml.Hadamard)
+
+        # # fourth layer subcircuit
+        assert len(tapes[3].operations) == 13
+        assert isinstance(tapes[3].operations[0], qml.RX)
+        assert isinstance(tapes[3].operations[1], qml.RY)
+        assert isinstance(tapes[3].operations[2], qml.CNOT)
+        assert isinstance(tapes[3].operations[3], qml.CNOT)
+        assert isinstance(tapes[3].operations[4], qml.RX)
+        assert isinstance(tapes[3].operations[5], qml.RY)
+        assert isinstance(tapes[3].operations[6], qml.RZ)
+        assert isinstance(tapes[3].operations[7], qml.CNOT)
+        assert isinstance(tapes[3].operations[8], qml.CNOT)
+        # PauliX decomp
+        assert isinstance(tapes[3].operations[9], qml.Hadamard)
+        # PauliY decomp
+        assert isinstance(tapes[3].operations[10], qml.PauliZ)
+        assert isinstance(tapes[3].operations[11], qml.S)
+        assert isinstance(tapes[3].operations[12], qml.Hadamard)
+
+    def test_evaluate_diag_metric_tensor(self, tol):
+        """Test that a diagonal metric tensor evaluates correctly for
+        block-diagonal and diagonal setting."""
+        dev = qml.device("default.qubit", wires=2)
+
+        def circuit(abc):
+            a, b, c = abc
+            qml.RX(a, wires=0)
+            qml.RY(b, wires=0)
+            qml.CNOT(wires=[0, 1])
+            qml.PhaseShift(c, wires=1)
+            return qml.expval(qml.PauliX(0)), qml.expval(qml.PauliX(1))
+
+        circuit = qml.QNode(circuit, dev)
+
+        abc = np.array([0.432, 0.12, -0.432], requires_grad=True)
+        a, b, c = abc
+
+        # evaluate metric tensor
+        g_diag = qml.metric_tensor(circuit, approx="diag")(abc)
+        g_blockdiag = qml.metric_tensor(circuit, approx="block-diag")(abc)
+
+        # check that the metric tensor is correct
+        expected = (
+            np.array(
+                [1, np.cos(a) ** 2, (3 - 2 * np.cos(a) ** 2 * np.cos(2 * b) - np.cos(2 * a)) / 4]
+            )
+            / 4
+        )
+        assert qml.math.allclose(g_diag, np.diag(expected), atol=tol, rtol=0)
+        assert qml.math.allclose(g_blockdiag, np.diag(expected), atol=tol, rtol=0)
+
+    @pytest.mark.parametrize("strategy", ["gradient", "device"])
+    def test_template_integration(self, strategy, tol):
+        """Test that the metric tensor transform acts on QNodes
+        correctly when the QNode contains a template"""
+        dev = qml.device("default.qubit", wires=3)
+
+        @qml.qnode(dev, expansion_strategy=strategy)
+        def circuit(weights):
+            qml.templates.StronglyEntanglingLayers(weights, wires=[0, 1, 2])
+            return qml.probs(wires=[0, 1])
+
+        weights = np.ones([2, 3, 3], dtype=np.float64, requires_grad=True)
+        res = qml.metric_tensor(circuit, approx="block-diag")(weights)
+        assert res.shape == (2, 3, 3, 2, 3, 3)
+
+    def test_evaluate_diag_metric_tensor_classical_processing(self, tol):
+        """Test that a diagonal metric tensor evaluates correctly
+        when the QNode includes classical processing."""
+        dev = qml.device("default.qubit", wires=2)
+
+        def circuit(a, b):
+            # The classical processing function is
+            #     f: ([a0, a1], b) -> (a1, a0, b)
+            # So the classical Jacobians will be a permutation matrix and an identity matrix:
+            #     classical_jacobian(circuit)(a, b) == ([[0, 1], [1, 0]], [[1]])
+            qml.RX(a[1], wires=0)
+            qml.RY(a[0], wires=0)
+            qml.CNOT(wires=[0, 1])
+            qml.U1(b, wires=1)
+            return qml.expval(qml.PauliX(0)), qml.expval(qml.PauliX(1))
+
+        circuit = qml.QNode(circuit, dev)
+
+        a = np.array([0.432, 0.1], requires_grad=True)
+        b = np.array(0.12, requires_grad=True)
+
+        # evaluate metric tensor
+        g = qml.metric_tensor(circuit, approx="block-diag")(a, b)
+        assert isinstance(g, tuple)
+        assert len(g) == 2
+        assert g[0].shape == (len(a), len(a))
+        assert g[1].shape == tuple()
+
+        # check that the metric tensor is correct
+        expected = np.array([np.cos(a[1]) ** 2, 1]) / 4
+        assert qml.math.allclose(g[0], np.diag(expected), atol=tol, rtol=0)
+
+        expected = (3 - 2 * np.cos(a[1]) ** 2 * np.cos(2 * a[0]) - np.cos(2 * a[1])) / 16
+        assert qml.math.allclose(g[1], expected, atol=tol, rtol=0)
+
+    @pytest.fixture(params=["parameter-shift", "backprop"])
+    def sample_circuit(self, request):
+        """Sample variational circuit fixture used in the
+        next couple of tests"""
+        dev = qml.device("default.qubit", wires=3)
+
+        def non_parametrized_layer(a, b, c):
+            qml.RX(a, wires=0)
+            qml.RX(b, wires=1)
+            qml.RX(c, wires=1)
+            qml.CNOT(wires=[0, 1])
+            qml.CNOT(wires=[1, 2])
+            qml.RZ(a, wires=0)
+            qml.Hadamard(wires=1)
+            qml.CNOT(wires=[0, 1])
+            qml.RZ(b, wires=1)
+            qml.Hadamard(wires=0)
+
+        a = 0.5
+        b = 0.1
+        c = 0.5
+
+        def final(params):
+            x, y, z, h, g, f = params
+            non_parametrized_layer(a, b, c)
+            qml.RX(x, wires=0)
+            qml.RY(-y, wires=1).inv()
+            qml.RZ(z, wires=2)
+            non_parametrized_layer(a, b, c)
+            qml.RY(f, wires=1)
+            qml.RZ(g, wires=2)
+            qml.RX(h, wires=1)
+            return qml.expval(qml.PauliX(0)), qml.expval(qml.PauliX(1)), qml.expval(qml.PauliX(2))
+
+        final = qml.QNode(final, dev, diff_method=request.param)
+
+        return dev, final, non_parametrized_layer, a, b, c
+
+    def test_evaluate_block_diag_metric_tensor(self, sample_circuit, tol):
+        """Test that a block-diagonal metric tensor evaluates correctly,
+        by comparing it to a known analytic result as well as numerical
+        computation."""
+        dev, circuit, non_parametrized_layer, a, b, c = sample_circuit
+
+        params = np.array(
+            [-0.282203, 0.145554, 0.331624, -0.163907, 0.57662, 0.081272],
+            requires_grad=True,
+        )
+
+        G = qml.metric_tensor(circuit, approx="block-diag")(params)
+
+        # ============================================
+        # Test block-diag metric tensor of first layer is correct.
+        # We do this by comparing against the known analytic result.
+        # First layer includes the non_parametrized_layer,
+        # followed by observables corresponding to generators of:
+        #   qml.RX(x, wires=0)
+        #   qml.RY(y, wires=1)
+        #   qml.RZ(z, wires=2)
+
+        G1 = np.zeros([3, 3])
+
+        # diag elements
+        G1[0, 0] = np.sin(a) ** 2 / 4
+        G1[1, 1] = (
+            16 * np.cos(a) ** 2 * np.sin(b) ** 3 * np.cos(b) * np.sin(2 * c)
+            + np.cos(2 * b) * (2 - 8 * np.cos(a) ** 2 * np.sin(b) ** 2 * np.cos(2 * c))
+            + np.cos(2 * (a - b))
+            + np.cos(2 * (a + b))
+            - 2 * np.cos(2 * a)
+            + 14
+        ) / 64
+        G1[2, 2] = (3 - np.cos(2 * a) - 2 * np.cos(a) ** 2 * np.cos(2 * (b + c))) / 16
+
+        # off diag elements
+        G1[0, 1] = np.sin(a) ** 2 * np.sin(b) * np.cos(b + c) / 4
+        G1[0, 2] = np.sin(a) ** 2 * np.cos(b + c) / 4
+        G1[1, 2] = (
+            -np.sin(b)
+            * (
+                np.cos(2 * (a - b - c))
+                + np.cos(2 * (a + b + c))
+                + 2 * np.cos(2 * a)
+                + 2 * np.cos(2 * (b + c))
+                - 6
+            )
+            / 32
+        )
+
+        G1[1, 0] = G1[0, 1]
+        G1[2, 0] = G1[0, 2]
+        G1[2, 1] = G1[1, 2]
+
+        assert qml.math.allclose(G[:3, :3], G1, atol=tol, rtol=0)
+
+        # =============================================
+        # Test block-diag metric tensor of second layer is correct.
+        # We do this by computing the required expectation values
+        # numerically using multiple circuits.
+        # The second layer includes the non_parametrized_layer,
+        # RX, RY, RZ gates (x, y, z params), and a 2nd non_parametrized_layer.
+        #
+        # Observables are the generators of:
+        #   qml.RY(f, wires=1)
+        #   qml.RZ(g, wires=2)
+        G2 = np.zeros([2, 2])
+
+        def layer2_diag(params):
+            x, y, z, h, g, f = params
+            non_parametrized_layer(a, b, c)
+            qml.RX(x, wires=0)
+            qml.RY(y, wires=1)
+            qml.RZ(z, wires=2)
+            non_parametrized_layer(a, b, c)
+            return qml.var(qml.PauliZ(2)), qml.var(qml.PauliY(1))
+
+        layer2_diag = qml.QNode(layer2_diag, dev)
+
+        def layer2_off_diag_first_order(params):
+            x, y, z, h, g, f = params
+            non_parametrized_layer(a, b, c)
+            qml.RX(x, wires=0)
+            qml.RY(y, wires=1)
+            qml.RZ(z, wires=2)
+            non_parametrized_layer(a, b, c)
+            return qml.expval(qml.PauliZ(2)), qml.expval(qml.PauliY(1))
+
+        layer2_off_diag_first_order = qml.QNode(layer2_off_diag_first_order, dev)
+
+        def layer2_off_diag_second_order(params):
+            x, y, z, h, g, f = params
+            non_parametrized_layer(a, b, c)
+            qml.RX(x, wires=0)
+            qml.RY(y, wires=1)
+            qml.RZ(z, wires=2)
+            non_parametrized_layer(a, b, c)
+            return qml.expval(qml.PauliY(1) @ qml.PauliZ(2))
+
+        layer2_off_diag_second_order = qml.QNode(layer2_off_diag_second_order, dev)
+
+        # calculate the diagonal terms
+        varK0, varK1 = layer2_diag(params)
+        G2[0, 0] = varK0 / 4
+        G2[1, 1] = varK1 / 4
+
+        # calculate the off-diagonal terms
+        exK0, exK1 = layer2_off_diag_first_order(params)
+        exK01 = layer2_off_diag_second_order(params)
+
+        G2[0, 1] = (exK01 - exK0 * exK1) / 4
+        G2[1, 0] = (exK01 - exK0 * exK1) / 4
+
+        assert qml.math.allclose(G[4:6, 4:6], G2, atol=tol, rtol=0)
+
+        # =============================================
+        # Test block-diag metric tensor of third layer is correct.
+        # We do this by computing the required expectation values
+        # numerically.
+        # The third layer includes the non_parametrized_layer,
+        # RX, RY, RZ gates (x, y, z params), a 2nd non_parametrized_layer,
+        # followed by the qml.RY(f, wires=2) operation.
+        #
+        # Observable is simply generator of:
+        #   qml.RY(f, wires=2)
+        #
+        # Note: since this layer only consists of a single parameter,
+        # only need to compute a single diagonal element.
+
+        def layer3_diag(params):
+            x, y, z, h, g, f = params
+            non_parametrized_layer(a, b, c)
+            qml.RX(x, wires=0)
+            qml.RY(y, wires=1)
+            qml.RZ(z, wires=2)
+            non_parametrized_layer(a, b, c)
+            qml.RY(f, wires=2)
+            return qml.var(qml.PauliX(1))
+
+        layer3_diag = qml.QNode(layer3_diag, dev)
+        G3 = layer3_diag(params) / 4
+        assert qml.math.allclose(G[3:4, 3:4], G3, atol=tol, rtol=0)
+
+        # ============================================
+        # Finally, double check that the entire metric
+        # tensor is as computed.
+
+        G_expected = block_diag(G1, G3, G2)
+        assert qml.math.allclose(G, G_expected, atol=tol, rtol=0)
+
+    def test_evaluate_diag_approx_metric_tensor(self, sample_circuit, tol):
+        """Test that a metric tensor under the diagonal approximation evaluates
+        correctly."""
+        dev, circuit, non_parametrized_layer, a, b, c = sample_circuit
+        params = np.array(
+            [-0.282203, 0.145554, 0.331624, -0.163907, 0.57662, 0.081272],
+            requires_grad=True,
+        )
+
+        G = qml.metric_tensor(circuit, approx="diag")(params)
+
+        # ============================================
+        # Test block-diag metric tensor of first layer is correct.
+        # We do this by comparing against the known analytic result.
+        # First layer includes the non_parametrized_layer,
+        # followed by observables corresponding to generators of:
+        #   qml.RX(x, wires=0)
+        #   qml.RY(y, wires=1)
+        #   qml.RZ(z, wires=2)
+
+        G1 = np.zeros([3, 3])
+
+        # diag elements
+        G1[0, 0] = np.sin(a) ** 2 / 4
+        G1[1, 1] = (
+            16 * np.cos(a) ** 2 * np.sin(b) ** 3 * np.cos(b) * np.sin(2 * c)
+            + np.cos(2 * b) * (2 - 8 * np.cos(a) ** 2 * np.sin(b) ** 2 * np.cos(2 * c))
+            + np.cos(2 * (a - b))
+            + np.cos(2 * (a + b))
+            - 2 * np.cos(2 * a)
+            + 14
+        ) / 64
+        G1[2, 2] = (3 - np.cos(2 * a) - 2 * np.cos(a) ** 2 * np.cos(2 * (b + c))) / 16
+
+        assert qml.math.allclose(G[:3, :3], G1, atol=tol, rtol=0)
+
+        # =============================================
+        # Test block-diag metric tensor of second layer is correct.
+        # We do this by computing the required expectation values
+        # numerically using multiple circuits.
+        # The second layer includes the non_parametrized_layer,
+        # RX, RY, RZ gates (x, y, z params), and a 2nd non_parametrized_layer.
+        #
+        # Observables are the generators of:
+        #   qml.RY(f, wires=1)
+        #   qml.RZ(g, wires=2)
+        G2 = np.zeros([2, 2])
+
+        def layer2_diag(params):
+            x, y, z, h, g, f = params
+            non_parametrized_layer(a, b, c)
+            qml.RX(x, wires=0)
+            qml.RY(y, wires=1)
+            qml.RZ(z, wires=2)
+            non_parametrized_layer(a, b, c)
+            return qml.var(qml.PauliZ(2)), qml.var(qml.PauliY(1))
+
+        layer2_diag = qml.QNode(layer2_diag, dev)
+
+        # calculate the diagonal terms
+        varK0, varK1 = layer2_diag(params)
+        G2[0, 0] = varK0 / 4
+        G2[1, 1] = varK1 / 4
+
+        assert qml.math.allclose(G[4:6, 4:6], G2, atol=tol, rtol=0)
+
+        # =============================================
+        # Test metric tensor of third layer is correct.
+        # We do this by computing the required expectation values
+        # numerically.
+        # The third layer includes the non_parametrized_layer,
+        # RX, RY, RZ gates (x, y, z params), a 2nd non_parametrized_layer,
+        # followed by the qml.RY(f, wires=2) operation.
+        #
+        # Observable is simply generator of:
+        #   qml.RY(f, wires=2)
+        #
+        # Note: since this layer only consists of a single parameter,
+        # only need to compute a single diagonal element.
+
+        def layer3_diag(params):
+            x, y, z, h, g, f = params
+            non_parametrized_layer(a, b, c)
+            qml.RX(x, wires=0)
+            qml.RY(y, wires=1)
+            qml.RZ(z, wires=2)
+            non_parametrized_layer(a, b, c)
+            qml.RY(f, wires=2)
+            return qml.var(qml.PauliX(1))
+
+        layer3_diag = qml.QNode(layer3_diag, dev)
+        G3 = layer3_diag(params) / 4
+        assert qml.math.allclose(G[3:4, 3:4], G3, atol=tol, rtol=0)
+
+        # ============================================
+        # Finally, double check that the entire metric
+        # tensor is as computed.
+
+        G_expected = block_diag(G1, G3, G2)
+        assert qml.math.allclose(G, G_expected, atol=tol, rtol=0)
+
+    def test_multi_qubit_gates(self):
+        """Test that a tape with Ising gates has the correct metric tensor tapes."""
+
+        dev = qml.device("default.qubit", wires=3)
+        with qml.tape.QuantumTape() as tape:
+            qml.Hadamard(0)
+            qml.Hadamard(2)
+            qml.IsingXX(0.2, wires=[0, 1])
+            qml.IsingXX(-0.6, wires=[1, 2])
+            qml.IsingZZ(1.02, wires=[0, 1])
+            qml.IsingZZ(-4.2, wires=[1, 2])
+
+        tapes, proc_fn = qml.metric_tensor(tape, approx="block-diag")
+        assert len(tapes) == 4
+        assert [len(tape.operations) for tape in tapes] == [3, 5, 4, 5]
+        assert [len(tape.measurements) for tape in tapes] == [1] * 4
+        expected_ops = [
+            [qml.Hadamard, qml.Hadamard, qml.Hadamard],
+            [qml.Hadamard, qml.Hadamard, qml.IsingXX, qml.Hadamard, qml.Hadamard],
+            [qml.Hadamard, qml.Hadamard, qml.IsingXX, qml.IsingXX],
+            [qml.Hadamard, qml.Hadamard, qml.IsingXX, qml.IsingXX, qml.IsingZZ],
+        ]
+        assert [[type(op) for op in tape.operations] for tape in tapes] == expected_ops
+
+    @pytest.mark.autograd
+    def test_no_trainable_params_qnode_autograd(self):
+        """Test that the correct ouput and warning is generated in the absence of any trainable
+        parameters"""
+
+        dev = qml.device("default.qubit", wires=3)
+
+        @qml.qnode(dev, interface="autograd")
+        def circuit(weights):
+            qml.RX(weights[0], wires=0)
+            qml.RY(weights[1], wires=0)
+            return qml.expval(qml.PauliZ(0) @ qml.PauliZ(1))
+
+        weights = [0.1, 0.2]
+        with pytest.warns(UserWarning, match="tensor of a QNode with no trainable parameters"):
+            res = qml.metric_tensor(circuit)(weights)
+
+        assert res == ()
+
+    @pytest.mark.torch
+    def test_no_trainable_params_qnode_torch(self):
+        """Test that the correct ouput and warning is generated in the absence of any trainable
+        parameters"""
+
+        dev = qml.device("default.qubit", wires=3)
+
+        @qml.qnode(dev, interface="torch")
+        def circuit(weights):
+            qml.RX(weights[0], wires=0)
+            qml.RY(weights[1], wires=0)
+            return qml.expval(qml.PauliZ(0) @ qml.PauliZ(1))
+
+        weights = [0.1, 0.2]
+        with pytest.warns(UserWarning, match="tensor of a QNode with no trainable parameters"):
+            res = qml.metric_tensor(circuit)(weights)
+
+        assert res == ()
+
+    @pytest.mark.tf
+    def test_no_trainable_params_qnode_tf(self):
+        """Test that the correct ouput and warning is generated in the absence of any trainable
+        parameters"""
+
+        dev = qml.device("default.qubit", wires=3)
+
+        @qml.qnode(dev, interface="tf")
+        def circuit(weights):
+            qml.RX(weights[0], wires=0)
+            qml.RY(weights[1], wires=0)
+            return qml.expval(qml.PauliZ(0) @ qml.PauliZ(1))
+
+        weights = [0.1, 0.2]
+        with pytest.warns(UserWarning, match="tensor of a QNode with no trainable parameters"):
+            res = qml.metric_tensor(circuit)(weights)
+
+        assert res == ()
+
+    @pytest.mark.jax
+    def test_no_trainable_params_qnode_jax(self):
+        """Test that the correct ouput and warning is generated in the absence of any trainable
+        parameters"""
+
+        dev = qml.device("default.qubit", wires=3)
+
+        @qml.qnode(dev, interface="jax")
+        def circuit(weights):
+            qml.RX(weights[0], wires=0)
+            qml.RY(weights[1], wires=0)
+            return qml.expval(qml.PauliZ(0) @ qml.PauliZ(1))
+
+        weights = [0.1, 0.2]
+        with pytest.warns(UserWarning, match="tensor of a QNode with no trainable parameters"):
+            res = qml.metric_tensor(circuit)(weights)
+
+        assert res == ()
+
+    def test_no_trainable_params_tape(self):
+        """Test that the correct ouput and warning is generated in the absence of any trainable
+        parameters"""
+        dev = qml.device("default.qubit", wires=3)
+
+        weights = [0.1, 0.2]
+        with qml.tape.QuantumTape() as tape:
+            qml.RX(weights[0], wires=0)
+            qml.RY(weights[1], wires=0)
+            qml.expval(qml.PauliZ(0) @ qml.PauliZ(1))
+
+        # TODO: remove once #2155 is resolved
+        tape.trainable_params = []
+
+        with pytest.warns(UserWarning, match="tensor of a tape with no trainable parameters"):
+            mt_tapes, post_processing = qml.metric_tensor(tape)
+        res = post_processing(qml.execute(mt_tapes, dev, None))
+
+        assert mt_tapes == []
+        assert res == ()
+
+
+fixed_pars = np.array([-0.2, 0.2, 0.5, 0.3, 0.7], requires_grad=False)
+
+
+def fubini_ansatz0(params, wires=None):
+    qml.RX(params[0], wires=0)
+    qml.RY(fixed_pars[0], wires=0)
+    qml.CNOT(wires=[wires[0], wires[1]])
+    qml.RZ(params[1], wires=0)
+    qml.CNOT(wires=[wires[0], wires[1]])
+
+
+def fubini_ansatz1(params, wires=None):
+    qml.RX(fixed_pars[1], wires=0)
+    for wire in wires:
+        qml.Rot(*params[0][wire], wires=wire)
+    qml.CNOT(wires=[0, 1])
+    qml.RY(fixed_pars[1], wires=0)
+    qml.CNOT(wires=[1, 2])
+    for wire in wires:
+        qml.Rot(*params[1][wire], wires=wire)
+    qml.CNOT(wires=[1, 2])
+    qml.RX(fixed_pars[2], wires=1)
+
+
+def fubini_ansatz2(params, wires=None):
+    params0 = params[0]
+    params1 = params[1]
+    qml.RX(fixed_pars[1], wires=0)
+    qml.Rot(*fixed_pars[2:5], wires=1)
+    qml.CNOT(wires=[0, 1])
+    qml.RY(params0, wires=0)
+    qml.RY(params0, wires=1)
+    qml.CNOT(wires=[0, 1])
+    qml.RX(params1, wires=0).inv()
+    qml.RX(params1, wires=1).inv()
+
+
+def fubini_ansatz3(params, wires=None):
+    params0 = params[0]
+    params1 = params[1]
+    params2 = params[2]
+    qml.RX(fixed_pars[1], wires=0)
+    qml.RX(fixed_pars[3], wires=1)
+    qml.CNOT(wires=[0, 1])
+    qml.CNOT(wires=[1, 2])
+    qml.RX(params0, wires=0).inv()
+    qml.RX(params0, wires=1)
+    qml.CNOT(wires=[0, 1])
+    qml.CNOT(wires=[1, 2])
+    qml.CNOT(wires=[2, 0])
+    qml.RY(params1, wires=0)
+    qml.RY(params1, wires=1)
+    qml.RY(params1, wires=2)
+    qml.RZ(params2, wires=0)
+    qml.RZ(params2, wires=1)
+    qml.RZ(params2, wires=2)
+
+
+def fubini_ansatz4(params00, params_rest, wires=None):
+    params01 = params_rest[0]
+    params10 = params_rest[1]
+    params11 = params_rest[2]
+    qml.RY(fixed_pars[3], wires=0)
+    qml.RY(fixed_pars[2], wires=1)
+    qml.CNOT(wires=[0, 1])
+    qml.CNOT(wires=[1, 2])
+    qml.RY(fixed_pars[4], wires=0)
+    qml.RX(params00, wires=0)
+    qml.CNOT(wires=[0, 1])
+    qml.RX(params01, wires=1)
+    qml.RZ(params10, wires=1)
+    qml.CNOT(wires=[0, 1])
+    qml.RZ(params11, wires=1)
+
+
+def fubini_ansatz5(params, wires=None):
+    fubini_ansatz4(params[0], [params[0], params[1], params[1]], wires=wires)
+
+
+def fubini_ansatz6(params, wires=None):
+    fubini_ansatz4(params[0], [params[0], params[1], -params[1]], wires=wires)
+
+
+def fubini_ansatz7(x, wires=None):
+    qml.RX(fixed_pars[0], wires=0)
+    qml.RX(x, wires=0)
+
+
+def fubini_ansatz8(params, wires=None):
+    params0 = params[0]
+    params1 = params[1]
+    qml.RX(fixed_pars[1], wires=[0])
+    qml.RY(fixed_pars[3], wires=[0])
+    qml.RZ(fixed_pars[2], wires=[0])
+    qml.RX(fixed_pars[2], wires=[1])
+    qml.RY(fixed_pars[2], wires=[1])
+    qml.RZ(fixed_pars[4], wires=[1])
+    qml.CNOT(wires=[0, 1])
+    qml.RX(fixed_pars[0], wires=[0])
+    qml.RY(fixed_pars[1], wires=[0])
+    qml.RZ(fixed_pars[3], wires=[0])
+    qml.RX(fixed_pars[1], wires=[1])
+    qml.RY(fixed_pars[2], wires=[1])
+    qml.RZ(fixed_pars[0], wires=[1])
+    qml.CNOT(wires=[0, 1])
+    qml.RX(params0, wires=[0])
+    qml.RX(params0, wires=[1])
+    qml.CNOT(wires=[0, 1])
+    qml.RY(fixed_pars[4], wires=[1])
+    qml.RY(params1, wires=[0])
+    qml.RY(params1, wires=[1])
+    qml.CNOT(wires=[0, 1])
+    qml.RX(fixed_pars[2], wires=[1])
+
+
+fubini_ansatze = [
+    fubini_ansatz0,
+    fubini_ansatz1,
+    fubini_ansatz2,
+    fubini_ansatz3,
+    fubini_ansatz4,
+    fubini_ansatz5,
+    fubini_ansatz6,
+    fubini_ansatz7,
+    fubini_ansatz8,
+]
+
+B = np.array(
+    [
+        [
+            [0.73, 0.49, 0.04],
+            [0.29, 0.45, 0.59],
+            [0.64, 0.06, 0.26],
+        ],
+        [
+            [0.93, 0.14, 0.46],
+            [0.31, 0.83, 0.79],
+            [0.25, 0.40, 0.16],
+        ],
+    ],
+    requires_grad=True,
+)
+fubini_params = [
+    (np.array([0.3434, -0.7245345], requires_grad=True),),
+    (B,),
+    (np.array([-0.1111, -0.2222], requires_grad=True),),
+    (np.array([-0.1111, -0.2222, 0.4554], requires_grad=True),),
+    (
+        np.array(-0.1735, requires_grad=True),
+        np.array([-0.1735, -0.2846, -0.2846], requires_grad=True),
+    ),
+    (np.array([-0.1735, -0.2846], requires_grad=True),),
+    (np.array([-0.1735, -0.2846], requires_grad=True),),
+    (np.array(-0.1735, requires_grad=True),),
+    (np.array([-0.1111, 0.3333], requires_grad=True),),
+]
+
+
+def autodiff_metric_tensor(ansatz, num_wires):
+    """Compute the metric tensor by full state vector
+    differentiation via autograd."""
+    dev = qml.device("default.qubit", wires=num_wires)
+
+    @qml.qnode(dev)
+    def qnode(*params):
+        ansatz(*params, wires=dev.wires)
+        return qml.state()
+
+    def mt(*params):
+        state = qnode(*params)
+        rqnode = lambda *params: np.real(qnode(*params))
+        iqnode = lambda *params: np.imag(qnode(*params))
+        rjac = qml.jacobian(rqnode)(*params)
+        ijac = qml.jacobian(iqnode)(*params)
+
+        if isinstance(rjac, tuple):
+            out = []
+            for rc, ic in zip(rjac, ijac):
+                c = rc + 1j * ic
+                psidpsi = np.tensordot(np.conj(state), c, axes=([0], [0]))
+                out.append(
+                    np.real(
+                        np.tensordot(np.conj(c), c, axes=([0], [0]))
+                        - np.tensordot(np.conj(psidpsi), psidpsi, axes=0)
+                    )
+                )
+            return tuple(out)
+
+        jac = rjac + 1j * ijac
+        psidpsi = np.tensordot(np.conj(state), jac, axes=([0], [0]))
+        return np.real(
+            np.tensordot(np.conj(jac), jac, axes=([0], [0]))
+            - np.tensordot(np.conj(psidpsi), psidpsi, axes=0)
+        )
+
+    return mt
+
+
+class TestFullMetricTensor:
+
+    num_wires = 3
+
+    @pytest.mark.autograd
+    @pytest.mark.parametrize("ansatz, params", zip(fubini_ansatze, fubini_params))
+    def test_correct_output_autograd(self, ansatz, params):
+        expected = autodiff_metric_tensor(ansatz, self.num_wires)(*params)
+        dev = qml.device("default.qubit.autograd", wires=self.num_wires + 1)
+
+        @qml.qnode(dev, interface="autograd")
+        def circuit(*params):
+            """Circuit with dummy output to create a QNode."""
+            ansatz(*params, dev.wires[:-1])
+            return qml.expval(qml.PauliZ(0))
+
+        mt = qml.metric_tensor(circuit, approx=None)(*params)
+
+        if isinstance(mt, tuple):
+            assert all(qml.math.allclose(_mt, _exp) for _mt, _exp in zip(mt, expected))
+        else:
+            assert qml.math.allclose(mt, expected)
+
+    @pytest.mark.jax
+    @pytest.mark.skip(reason="JAX does not support the forward pass metric tensor.")
+    @pytest.mark.parametrize("ansatz, params", zip(fubini_ansatze, fubini_params))
+    def test_correct_output_jax(self, ansatz, params):
+        from jax import numpy as jnp
+
+        expected = autodiff_metric_tensor(ansatz, self.num_wires)(*params)
+        dev = qml.device("default.qubit.jax", wires=self.num_wires + 1)
+
+        params = tuple(jnp.array(p) for p in params)
+
+        @qml.qnode(dev, interface="jax")
+        def circuit(*params):
+            """Circuit with dummy output to create a QNode."""
+            ansatz(*params, dev.wires[:-1])
+            return qml.expval(qml.PauliZ(0))
+
+        mt = qml.metric_tensor(circuit, approx=None)(*params)
+
+        if isinstance(mt, tuple):
+            assert all(qml.math.allclose(_mt, _exp) for _mt, _exp in zip(mt, expected))
+        else:
+            assert qml.math.allclose(mt, expected)
+
+    @pytest.mark.torch
+    @pytest.mark.parametrize("ansatz, params", zip(fubini_ansatze, fubini_params))
+    def test_correct_output_torch(self, ansatz, params):
+        import torch
+
+        expected = autodiff_metric_tensor(ansatz, self.num_wires)(*params)
+        dev = qml.device("default.qubit.torch", wires=self.num_wires + 1)
+
+        params = tuple(torch.tensor(p, dtype=torch.float64, requires_grad=True) for p in params)
+
+        @qml.qnode(dev, interface="torch")
+        def circuit(*params):
+            """Circuit with dummy output to create a QNode."""
+            ansatz(*params, dev.wires[:-1])
+            return qml.expval(qml.PauliZ(0))
+
+        mt = qml.metric_tensor(circuit, approx=None)(*params)
+
+        if isinstance(mt, tuple):
+            assert all(qml.math.allclose(_mt, _exp) for _mt, _exp in zip(mt, expected))
+        else:
+            assert qml.math.allclose(mt, expected)
+
+    @pytest.mark.tf
+    @pytest.mark.parametrize("ansatz, params", zip(fubini_ansatze, fubini_params))
+    def test_correct_output_tf(self, ansatz, params):
+        import tensorflow as tf
+
+        expected = autodiff_metric_tensor(ansatz, self.num_wires)(*params)
+        dev = qml.device("default.qubit.tf", wires=self.num_wires + 1)
+
+        params = tuple(tf.Variable(p, dtype=tf.float64) for p in params)
+
+        @qml.qnode(dev, interface="tf")
+        def circuit(*params):
+            """Circuit with dummy output to create a QNode."""
+            ansatz(*params, dev.wires[:-1])
+            return qml.expval(qml.PauliZ(0))
+
+        with tf.GradientTape() as t:
+            qml.metric_tensor(circuit, approx="block-diag")(*params)
+            mt = qml.metric_tensor(circuit, approx=None)(*params)
+
+        if isinstance(mt, tuple):
+            assert all(qml.math.allclose(_mt, _exp) for _mt, _exp in zip(mt, expected))
+        else:
+            assert qml.math.allclose(mt, expected)
+
+
+def diffability_ansatz_0(weights, wires=None):
+    qml.RX(weights[0], wires=0)
+    qml.RX(weights[1], wires=0)
+    qml.CNOT(wires=[0, 1])
+    qml.RZ(weights[2], wires=1)
+
+
+expected_diag_jac_0 = lambda weights: np.array(
+    [
+        [0, 0, 0],
+        [0, 0, 0],
+        [
+            np.cos(weights[0] + weights[1]) * np.sin(weights[0] + weights[1]) / 2,
+            np.cos(weights[0] + weights[1]) * np.sin(weights[0] + weights[1]) / 2,
+            0,
+        ],
+    ]
+)
+
+
+def diffability_ansatz_1(weights, wires=None):
+    qml.RX(weights[0], wires=0)
+    qml.RY(weights[1], wires=0)
+    qml.CNOT(wires=[0, 1])
+    qml.RZ(weights[2], wires=1)
+
+
+expected_diag_jac_1 = lambda weights: np.array(
+    [
+        [0, 0, 0],
+        [-np.sin(2 * weights[0]) / 4, 0, 0],
+        [
+            np.cos(weights[0]) * np.cos(weights[1]) ** 2 * np.sin(weights[0]) / 2,
+            np.cos(weights[0]) ** 2 * np.sin(2 * weights[1]) / 4,
+            0,
+        ],
+    ]
+)
+
+
+def diffability_ansatz_2(weights, wires=None):
+    qml.RX(weights[0], wires=0)
+    qml.RY(weights[1], wires=1)
+    qml.CNOT(wires=[0, 1])
+    qml.RZ(weights[2], wires=1)
+
+
+expected_diag_jac_2 = lambda weights: np.array(
+    [
+        [0, 0, 0],
+        [0, 0, 0],
+        [
+            np.cos(weights[1]) ** 2 * np.sin(2 * weights[0]) / 4,
+            np.cos(weights[0]) ** 2 * np.sin(2 * weights[1]) / 4,
+            0,
+        ],
+    ]
+)
+
+weights = np.array([0.432, 0.12, -0.292], requires_grad=True)
+
+
+@pytest.mark.parametrize("diff_method", ["backprop", "parameter-shift"])
+@pytest.mark.parametrize(
+    "ansatz, weights, expected_diag_jac",
+    [
+        (diffability_ansatz_0, (weights,), expected_diag_jac_0),
+        (diffability_ansatz_1, (weights,), expected_diag_jac_1),
+        (diffability_ansatz_2, (weights,), expected_diag_jac_2),
+    ],
+)
+class TestDifferentiability:
+    """Test for metric tensor differentiability"""
+
+    def get_circuit(self, ansatz):
+        def circuit(*weights):
+            ansatz(*weights)
+            return qml.expval(qml.PauliX(0)), qml.expval(qml.PauliX(1))
+
+        return circuit
+
+    dev = qml.device("default.qubit", wires=3)
+
+    @pytest.mark.autograd
+    def test_autograd_diag(self, diff_method, tol, ansatz, weights, expected_diag_jac):
+        """Test metric tensor differentiability in the autograd interface"""
+        circuit = self.get_circuit(ansatz)
+        qnode = qml.QNode(circuit, self.dev, interface="autograd", diff_method=diff_method)
+        qnode(*weights)
+
+        def cost_diag(*weights):
+            mt = qml.metric_tensor(qnode, approx="block-diag")(*weights)
+            if isinstance(mt, tuple):
+                diag = qml.math.hstack(
+                    [qml.math.diag(_mt) if len(qml.math.shape(_mt)) == 2 else _mt for _mt in mt]
+                )
+
+            else:
+                diag = qml.math.diag(mt)
+            return diag
+
+        jac = qml.jacobian(cost_diag)(*weights)
+        if isinstance(jac, tuple):
+            assert all(
+                qml.math.allclose(j, e, atol=tol, rtol=0)
+                for j, e in zip(jac, expected_diag_jac(*weights))
+            )
+        else:
+            assert qml.math.allclose(jac, expected_diag_jac(*weights), atol=tol, rtol=0)
+
+    @pytest.mark.autograd
+    def test_autograd(self, diff_method, tol, ansatz, weights, expected_diag_jac):
+        """Test metric tensor differentiability in the autograd interface"""
+        circuit = self.get_circuit(ansatz)
+        qnode = qml.QNode(circuit, self.dev, interface="autograd", diff_method=diff_method)
+
+        def cost_full(*weights):
+            return np.array(qml.metric_tensor(qnode, approx=None)(*weights))
+
+        _cost_full = lambda *weights: np.array(autodiff_metric_tensor(ansatz, 3)(*weights))
+        _c = _cost_full(*weights)
+        c = cost_full(*weights)
+        assert all(
+            qml.math.allclose(_sub_c, sub_c, atol=tol, rtol=0) for _sub_c, sub_c in zip(_c, c)
+        )
+        for argnum in range(len(weights)):
+            expected_full = qml.jacobian(_cost_full, argnum=argnum)(*weights)
+            jac = qml.jacobian(cost_full, argnum=argnum)(*weights)
+            assert qml.math.allclose(expected_full, jac, atol=tol, rtol=0)
+
+    @pytest.mark.jax
+    def test_jax_diag(self, diff_method, tol, ansatz, weights, expected_diag_jac):
+        """Test metric tensor differentiability in the JAX interface"""
+        if diff_method == "parameter-shift":
+            pytest.skip("Does not support parameter-shift")
+
+        import jax
+        from jax import numpy as jnp
+
+        circuit = self.get_circuit(ansatz)
+        qnode = qml.QNode(circuit, self.dev, interface="jax", diff_method=diff_method)
+
+        def cost_diag(*weights):
+            return jnp.diag(qml.metric_tensor(qnode, approx="block-diag")(*weights))
+
+        jac = jax.jacobian(cost_diag)(jnp.array(*weights))
+        assert qml.math.allclose(jac, expected_diag_jac(*weights), atol=tol, rtol=0)
+
+    @pytest.mark.jax
+    def test_jax(self, diff_method, tol, ansatz, weights, expected_diag_jac):
+        """Test metric tensor differentiability in the JAX interface"""
+        if diff_method == "parameter-shift":
+            pytest.skip("Does not support parameter-shift")
+
+        import jax
+
+        circuit = self.get_circuit(ansatz)
+        qnode = qml.QNode(circuit, self.dev, interface="jax", diff_method=diff_method)
+
+        def cost_full(*weights):
+            return qml.metric_tensor(qnode, approx=None)(*weights)
+
+        _cost_full = lambda *weights: autodiff_metric_tensor(ansatz, num_wires=3)(*weights)
+        assert qml.math.allclose(_cost_full(*weights), cost_full(*weights), atol=tol, rtol=0)
+        jac = jax.jacobian(cost_full)(*weights)
+        expected_full = qml.jacobian(_cost_full)(*weights)
+        assert qml.math.allclose(expected_full, jac, atol=tol, rtol=0)
+
+    @pytest.mark.tf
+    def test_tf_diag(self, diff_method, tol, ansatz, weights, expected_diag_jac):
+        """Test metric tensor differentiability in the TF interface"""
+        import tensorflow as tf
+
+        circuit = self.get_circuit(ansatz)
+        qnode = qml.QNode(circuit, self.dev, interface="tf", diff_method=diff_method)
+
+        weights_t = tuple(tf.Variable(w) for w in weights)
+        with tf.GradientTape() as tape:
+            loss_diag = tf.linalg.diag_part(
+                qml.metric_tensor(qnode, approx="block-diag")(*weights_t)
+            )
+        jac = tape.jacobian(loss_diag, weights_t)
+        assert qml.math.allclose(jac, expected_diag_jac(*weights), atol=tol, rtol=0)
+
+    @pytest.mark.tf
+    def test_tf(self, diff_method, tol, ansatz, weights, expected_diag_jac):
+        """Test metric tensor differentiability in the TF interface"""
+        import tensorflow as tf
+
+        circuit = self.get_circuit(ansatz)
+        qnode = qml.QNode(circuit, self.dev, interface="tf", diff_method=diff_method)
+
+        weights_t = tuple(tf.Variable(w) for w in weights)
+        with tf.GradientTape() as tape:
+            loss_full = qml.metric_tensor(qnode, approx=None)(*weights_t)
+        jac = tape.jacobian(loss_full, weights_t)
+        _cost_full = lambda *weights: autodiff_metric_tensor(ansatz, num_wires=3)(*weights)
+        assert qml.math.allclose(_cost_full(*weights), loss_full, atol=tol, rtol=0)
+        expected_full = qml.jacobian(_cost_full)(*weights)
+        assert qml.math.allclose(expected_full, jac, atol=tol, rtol=0)
+
+    @pytest.mark.torch
+    def test_torch_diag(self, diff_method, tol, ansatz, weights, expected_diag_jac):
+        """Test metric tensor differentiability in the torch interface"""
+        import torch
+
+        circuit = self.get_circuit(ansatz)
+        qnode = qml.QNode(circuit, self.dev, interface="torch", diff_method=diff_method)
+
+        weights_t = tuple(torch.tensor(w, requires_grad=True) for w in weights)
+
+        def cost_diag(*weights):
+            mt = qml.metric_tensor(qnode, approx="block-diag")(*weights)
+            if isinstance(mt, tuple):
+                diag = qml.math.hstack(
+                    [qml.math.diag(_mt) if len(qml.math.shape(_mt)) == 2 else _mt for _mt in mt]
+                )
+
+            else:
+                diag = qml.math.diag(mt)
+            return diag
+
+        jac = torch.autograd.functional.jacobian(cost_diag, weights_t)
+
+        if isinstance(jac, tuple) and len(jac) != 1:
+            assert all(
+                qml.math.allclose(j.detach().numpy(), e, atol=tol, rtol=0)
+                for j, e in zip(jac, expected_diag_jac(*weights))
+            )
+        else:
+            if isinstance(jac, tuple) and len(jac) == 1:
+                jac = jac[0]
+            assert qml.math.allclose(
+                jac.detach().numpy(), expected_diag_jac(*weights), atol=tol, rtol=0
+            )
+
+    @pytest.mark.torch
+    def test_torch(self, diff_method, tol, ansatz, weights, expected_diag_jac):
+        """Test metric tensor differentiability in the torch interface"""
+        import torch
+
+        circuit = self.get_circuit(ansatz)
+        qnode = qml.QNode(circuit, self.dev, interface="torch", diff_method=diff_method)
+        weights_t = tuple(torch.tensor(w, requires_grad=True) for w in weights)
+        qnode(*weights_t)
+        cost_full = qml.metric_tensor(qnode, approx=None)
+        _cost_full = autodiff_metric_tensor(ansatz, num_wires=3)
+        jac = torch.autograd.functional.jacobian(cost_full, weights_t)
+        expected_full = qml.jacobian(_cost_full)(*weights)
+        assert qml.math.allclose(
+            _cost_full(*weights), cost_full(*weights_t).detach().numpy(), atol=tol, rtol=0
+        )
+        if isinstance(jac, tuple) and len(jac) != 1:
+            assert all(
+                qml.math.allclose(j.detach().numpy(), e, atol=tol, rtol=0)
+                for j, e in zip(jac, expected_full)
+            )
+        else:
+            if isinstance(jac, tuple) and len(jac) == 1:
+                jac = jac[0]
+            assert qml.math.allclose(jac.detach().numpy(), expected_full, atol=tol, rtol=0)
+
+
+@pytest.mark.parametrize("approx", [True, False, "Invalid", 2])
+def test_invalid_value_for_approx(approx):
+    """Test exception is raised if ``approx`` is invalid."""
+    with qml.tape.QuantumTape() as tape:
+        qml.RX(np.array(0.5, requires_grad=True), wires=0)
+        qml.expval(qml.PauliX(0))
+
+    with pytest.raises(ValueError, match="keyword argument approx"):
+        qml.metric_tensor(tape, approx=approx)
+
+
+def test_generator_no_expval(monkeypatch):
+    """Test exception is raised if subcircuit contains an
+    operation with generator object that is not an observable"""
+    with monkeypatch.context() as m:
+        m.setattr("pennylane.RX.generator", lambda self: qml.RX(0.1, wires=0))
+
+        with qml.tape.QuantumTape() as tape:
+            qml.RX(np.array(0.5, requires_grad=True), wires=0)
+            qml.expval(qml.PauliX(0))
+
+        with pytest.raises(qml.QuantumFunctionError, match="is not an observable"):
+            qml.metric_tensor(tape, approx="block-diag")
+
+
+def test_error_missing_aux_wire():
+    """Tests that a special error is raised if the requested (or default, if not given)
+    auxiliary wire for the Hadamard test is missing."""
+    dev = qml.device("default.qubit", wires=qml.wires.Wires(["wire1", "wire2"]))
+
+    @qml.qnode(dev)
+    def circuit(x, z):
+        qml.RX(x, wires="wire1")
+        qml.RZ(z, wires="wire2")
+        qml.CNOT(wires=["wire1", "wire2"])
+        qml.RX(x, wires="wire1")
+        qml.RZ(z, wires="wire2")
+        return qml.expval(qml.PauliZ("wire2"))
+
+    x = np.array(0.5, requires_grad=True)
+    z = np.array(0.1, requires_grad=True)
+
+    with pytest.warns(UserWarning, match="The device does not have a wire that is not used"):
+        qml.metric_tensor(circuit, approx=None)(x, z)
+
+
+def test_error_not_available_aux_wire():
+    """Tests that a special error is raised if aux wires is not available."""
+
+    dev = qml.device("default.qubit", wires=1)
+    x = np.array(0.5, requires_grad=True)
+
+    @qml.qnode(dev)
+    def circuit(x):
+        qml.RX(x, wires=0)
+        qml.RY(x, wires=0)
+        return qml.expval(qml.PauliZ(0))
+
+    with pytest.warns(UserWarning, match="An auxiliary wire is not available."):
+        qml.metric_tensor(circuit, aux_wire=404)(x)
+
+
+def test_error_aux_wire_replaced():
+    """Tests that even if an aux_wire is provided, it is superseded by a device
+    wire if it does not exist itself on the device, so that the metric_tensor is
+    successfully computed."""
+    dev = qml.device("default.qubit", wires=qml.wires.Wires(["wire1", "wire2", "hidden_wire"]))
+
+    @qml.qnode(dev)
+    def circuit(x, z):
+        qml.RX(x, wires="wire1")
+        qml.RZ(z, wires="wire2")
+        qml.CNOT(wires=["wire1", "wire2"])
+        qml.RX(x, wires="wire1")
+        qml.RZ(z, wires="wire2")
+        return qml.expval(qml.PauliZ("wire2"))
+
+    x = np.array(0.5, requires_grad=True)
+    z = np.array(0.1, requires_grad=True)
+
+    qml.metric_tensor(circuit, approx=None, aux_wire="wire3")(x, z)
+
+
+@pytest.mark.parametrize("allow_nonunitary", [True, False])
+def test_error_generator_not_registered(allow_nonunitary, monkeypatch):
+    """Tests that an error is raised if an operation doe not have a
+    controlled-generator operation registered."""
+    dev = qml.device("default.qubit", wires=qml.wires.Wires(["wire1", "wire2", "wire3"]))
+
+    x = np.array(0.5, requires_grad=True)
+    z = np.array(0.1, requires_grad=True)
+
+    @qml.qnode(dev)
+    def circuit(x, z):
+        qml.CRX(x, wires=["wire1", "wire2"])
+        qml.RZ(z, wires="wire2")
+        return qml.expval(qml.PauliZ("wire2"))
+
+    with monkeypatch.context() as m:
+        exp_fn = lambda tape, *args, **kwargs: tape
+        m.setattr("pennylane.transforms.metric_tensor.expand_fn", exp_fn)
+
+        if allow_nonunitary:
+            qml.metric_tensor(circuit, approx=None, allow_nonunitary=allow_nonunitary)(x, z)
+        else:
+            with pytest.raises(ValueError, match="Generator for operation"):
+                qml.metric_tensor(circuit, approx=None, allow_nonunitary=allow_nonunitary)(x, z)
+
+    class RX(qml.RX):
+        def generator(self):
+            return qml.Hadamard(self.wires)
+
+    @qml.qnode(dev)
+    def circuit(x, z):
+        RX(x, wires="wire1")
+        qml.RZ(z, wires="wire1")
+        return qml.expval(qml.PauliZ("wire2"))
+
+    with monkeypatch.context() as m:
+        exp_fn = lambda tape, *args, **kwargs: tape
+        m.setattr("pennylane.transforms.metric_tensor.expand_fn", exp_fn)
+
+        if allow_nonunitary:
+            qml.metric_tensor(circuit, approx=None, allow_nonunitary=allow_nonunitary)(x, z)
+        else:
+            with pytest.raises(ValueError, match="Generator for operation"):
+                qml.metric_tensor(circuit, approx=None, allow_nonunitary=allow_nonunitary)(x, z)
+
+
+def test_no_error_missing_aux_wire_not_used(recwarn):
+    """Tests that a no error is raised if the requested (or default, if not given)
+    auxiliary wire for the Hadamard test is missing but it is not used, either
+    because ``approx`` is used or because there only is a diagonal contribution."""
+    dev = qml.device("default.qubit", wires=qml.wires.Wires(["wire1", "wire2"]))
+
+    @qml.qnode(dev)
+    def circuit_single_block(x, z):
+        """This circuit has a metric tensor that consists
+        of a single block in the block diagonal "approximation"."""
+        qml.RX(x, wires="wire1")
+        qml.RZ(z, wires="wire2")
+        qml.CNOT(wires=["wire1", "wire2"])
+        return qml.expval(qml.PauliZ("wire2"))
+
+    @qml.qnode(dev)
+    def circuit_multi_block(x, z):
+        """This circuit has a metric tensor that consists
+        of multiple blocks and thus is approximated when only
+        computing the block diagonal."""
+        qml.RX(x, wires="wire1")
+        qml.RZ(z, wires="wire2")
+        qml.CNOT(wires=["wire1", "wire2"])
+        qml.RX(x, wires="wire1")
+        qml.RZ(z, wires="wire2")
+        return qml.expval(qml.PauliZ("wire2"))
+
+    x = np.array(0.5, requires_grad=True)
+    z = np.array(0.1, requires_grad=True)
+
+    qml.metric_tensor(circuit_single_block, approx=None)(x, z)
+    qml.metric_tensor(circuit_single_block, approx=None, aux_wire="aux_wire")(x, z)
+    qml.metric_tensor(circuit_multi_block, approx="block-diag")(x, z)
+    qml.metric_tensor(circuit_multi_block, approx="block-diag", aux_wire="aux_wire")(x, z)
+
+    assert len(recwarn) == 0
+
+
+def aux_wire_ansatz_0(x, y):
+    qml.RX(x, wires=0)
+    qml.RY(x, wires=2)
+
+
+def aux_wire_ansatz_1(x, y):
+    qml.RX(x, wires=0)
+    qml.RY(x, wires=1)
+
+
+@pytest.mark.parametrize("aux_wire", [None, "aux", 3])
+@pytest.mark.parametrize("ansatz", [aux_wire_ansatz_0, aux_wire_ansatz_1])
+def test_get_aux_wire(aux_wire, ansatz):
+    """Test ``_get_aux_wire`` without device_wires."""
+    x, y = np.array([0.2, 0.1], requires_grad=True)
+    with qml.tape.QuantumTape() as tape:
+        ansatz(x, y)
+    out = _get_aux_wire(aux_wire, tape, None)
+
+    if aux_wire is not None:
+        assert out == aux_wire
+    else:
+        assert out == (1 if 1 not in tape.wires else 2)
+
+
+def test_get_aux_wire_with_device_wires():
+    """Test ``_get_aux_wire`` with device_wires."""
+    x, y = np.array([0.2, 0.1], requires_grad=True)
+    with qml.tape.QuantumTape() as tape:
+        qml.RX(x, wires=0)
+        qml.RX(x, wires="one")
+
+    device_wires = qml.wires.Wires([0, "aux", "one"])
+
+    assert _get_aux_wire(0, tape, device_wires) == 0
+    assert _get_aux_wire("one", tape, device_wires) == "one"
+    assert _get_aux_wire(None, tape, device_wires) == "aux"
+
+
+def test_get_aux_wire_with_unavailable_aux():
+    """Test ``_get_aux_wire`` with device_wires and a requested ``aux_wire`` that is missing."""
+    x, y = np.array([0.2, 0.1], requires_grad=True)
+    with qml.tape.QuantumTape() as tape:
+        qml.RX(x, wires=0)
+        qml.RX(x, wires="one")
+    device_wires = qml.wires.Wires([0, "one"])
+    with pytest.raises(qml.wires.WireError, match="The requested aux_wire does not exist"):
+        _get_aux_wire("two", tape, device_wires)