Update qml.matrix docs for wire ordering with qnodes (#4874)

**Context:**
`qml.matrix` wire ordering with the new `default.qubit` was causing
problems when `device.wires=None`. This is because when
`device.wires=None`, the tape wires are used, the order of which is
determined by which wires are used first, even with sequential integer
wire labels.

**Description of the Change:**
* Added documentation to `qml.matrix` to provide more details on how the
wire order is determined when using `qml.matrix` with a `QNode`.

**Benefits:**
Users may be less confused if `qml.matrix` gives a matrix with an
unexpected wire order.

**Possible Drawbacks:**

**Related GitHub Issues:**

---------

Co-authored-by: Matthew Silverman <[email protected]>
Co-authored-by: lillian542 <[email protected]>

doc/releases/changelog-dev.md

@@ -188,6 +188,10 @@
   to their top-level pages (rather than their module-level pages, eg. `qml.measurements.expval`).
   [(#4750)](https://github.com/PennyLaneAI/pennylane/pull/4750)
 
+* Added information to documentation of `qml.matrix` about wire ordering when using `qml.matrix` on a
+  `QNode` which uses a device with `device.wires=None`.
+  [(#4874)](https://github.com/PennyLaneAI/pennylane/pull/4874)
+
 <h3>Bug fixes 🐛</h3>
 
 * `qml.grad` and `qml.jacobian` now explicitly raise errors if trainable parameters are integers.

pennylane/ops/functions/matrix.py

@@ -30,7 +30,15 @@ def matrix(op: qml.operation.Operator, wire_order=None) -> TensorLike:
     Args:
         op (Operator or QNode or QuantumTape or Callable): A quantum operator or quantum circuit.
         wire_order (Sequence[Any], optional): Order of the wires in the quantum circuit.
-            Defaults to the order in which the wires appear in the quantum function.
+            The default wire order depends on the type of ``op``:
+
+            - If ``op`` is a :class:`~.QNode`, then the wire order is determined by the
+              associated device's wires, if provided.
+
+            - Otherwise, the wire order is determined by the order in which wires
+              appear in the circuit.
+
+            - See the usage details for more information.
 
     Returns:
         TensorLike or qnode (QNode) or quantum function (Callable) or tuple[List[QuantumTape], function]:
@@ -123,6 +131,41 @@ def cost(theta):
         1.9775421558720845
         >>> qml.grad(cost)(theta)
         -0.14943813247359922
+
+        .. note::
+
+            When using ``qml.matrix`` with a ``QNode``, unless specified, the device wire order will
+            be used. If the device wires are not set, the wire order will be inferred
+            from the quantum function used to create the ``QNode``. Consider the following example:
+
+            .. code-block:: python
+
+                def circuit():
+                    qml.Hadamard(wires=1)
+                    qml.CZ(wires=[0, 1])
+                    qml.Hadamard(wires=1)
+                    return qml.state()
+
+                dev_with_wires = qml.device("default.qubit", wires=[0, 1])
+                dev_without_wires = qml.device("default.qubit")
+
+                qnode_with_wires = qml.QNode(circuit, dev_with_wires)
+                qnode_without_wires = qml.QNode(circuit, dev_without_wires)
+
+            >>> qml.matrix(qnode_with_wires)().round(2)
+            array([[ 1.+0.j, -0.+0.j,  0.+0.j,  0.+0.j],
+                   [-0.+0.j,  1.+0.j,  0.+0.j,  0.+0.j],
+                   [ 0.+0.j,  0.+0.j, -0.+0.j,  1.+0.j],
+                   [ 0.+0.j,  0.+0.j,  1.+0.j, -0.+0.j]])
+            >>> qml.matrix(qnode_without_wires)().round(2)
+            array([[ 1.+0.j,  0.+0.j, -0.+0.j,  0.+0.j],
+                   [ 0.+0.j, -0.+0.j,  0.+0.j,  1.+0.j],
+                   [-0.+0.j,  0.+0.j,  1.+0.j,  0.+0.j],
+                   [ 0.+0.j,  1.+0.j,  0.+0.j, -0.+0.j]])
+
+            The second matrix above uses wire order ``[1, 0]`` because the device does not have
+            wires specified, and this is the order in which wires appear in ``circuit()``.
+
     """
     if not isinstance(op, qml.operation.Operator):
         if not isinstance(op, (qml.tape.QuantumScript, qml.QNode)) and not callable(op):