Repository: PennyLaneAI/pennylane (python)

Issue:
[BUG] imaginary part of expectation values of observables with complex coefficients is silently discarded
### Expected behavior

The PL operator arithmetic allows for operators with imaginary/complex  coefficients. If such expectation values of such operators are measured I expect to either get a imaginary/complex value (the real expectation value of the Hermitian operator times the complex coefficient) or an exception that only Hermitian observables can be measured.

### Actual behavior

Instead PennyLane silently discards the imaginary part, as the attached example code shows. 

### Additional information

_No response_

### Source code

```shell
import pennylane as qml

obs = (0.2 + 0.7j) * qml.PauliZ(0)

@qml.qnode(qml.device("default.qubit", 1))
def qnode():
    return qml.expval(obs)

qnode()
```


### Tracebacks

_No response_

### System information

```shell
Name: PennyLane
Version: 0.35.1
Summary: PennyLane is a cross-platform Python library for quantum computing, quantum machine learning, and quantum chemistry. Train a quantum computer the same way as a neural network.
Home-page: https://github.com/PennyLaneAI/pennylane
Author:
Author-email:
License: Apache License 2.0
Location: /home/cvjjm/miniforge3/envs/[...]/lib/python3.11/site-packages
Requires: appdirs, autograd, autoray, cachetools, networkx, numpy, pennylane-lightning, requests, rustworkx, scipy, semantic-version, toml, typing-extensions
Required-by: matchgate_shadows, PennyLane_Lightning

Platform info:           Linux-5.15.153.1-microsoft-standard-WSL2-x86_64-with-glibc2.31
Python version:          3.11.9
Numpy version:           1.23.5
Scipy version:           1.13.0
Installed devices:
- lightning.qubit (PennyLane_Lightning-0.35.1)
- default.clifford (PennyLane-0.35.1)
- default.gaussian (PennyLane-0.35.1)
- default.mixed (PennyLane-0.35.1)
- default.qubit (PennyLane-0.35.1)
- default.qubit.autograd (PennyLane-0.35.1)
- default.qubit.jax (PennyLane-0.35.1)
- default.qubit.legacy (PennyLane-0.35.1)
- default.qubit.tf (PennyLane-0.35.1)
- default.qubit.torch (PennyLane-0.35.1)
- default.qutrit (PennyLane-0.35.1)
- null.qubit (PennyLane-0.35.1)
```


### Existing GitHub issues

- [X] I have searched existing GitHub issues to make sure the issue does not already exist.

Pull request:
Update `qml.expval` to correctly cast value type
**Context:**

Prior to this PR we had the expectation value being silently converted to a real type, 
```python
obs = (0.2 + 0.7j) * qml.PauliZ(0)

@qml.qnode(qml.device("default.qubit", 1))
def qnode():
    return qml.expval(obs)

>>> qnode()
0.2
```

**Description of the Change:**

This PR stops casting to a real number when calculating expval and instead casts to whatever the eigenvalue type is.

**Benefits:**

Code is now correct with the following behaviour,
```python
obs = (0.2 + 0.7j) * qml.PauliZ(0)

@qml.qnode(qml.device("default.qubit", 1))
def qnode():
    return qml.expval(obs)

>>> qnode()
(0.2+0.7j)
```

**Possible Drawbacks:** The math is now correct, but this result is non-physical. However, after some discussion we do not want to support silently type casting and instead leave this to the researchers.

**Related GitHub Issues:** Fixes #6076 

[sc-70629]


Patch:
diff --git a/doc/releases/changelog-dev.md b/doc/releases/changelog-dev.md
index 1e56914c5..b5eb82fb2 100644
--- a/doc/releases/changelog-dev.md
+++ b/doc/releases/changelog-dev.md
@@ -282,6 +282,9 @@
 
 <h3>Bug fixes 🐛</h3>
 
+* `qml.expval` no longer silently casts to a real number when observable coefficients are imaginary.
+  [(#6939)](https://github.com/PennyLaneAI/pennylane/pull/6939)
+
 * Fixed `qml.wires.Wires` initialization to disallow `Wires` objects as wires labels.
   Now, `Wires` is idempotent, e.g. `Wires([Wires([0]), Wires([1])])==Wires([0, 1])`.
   [(#6933)](https://github.com/PennyLaneAI/pennylane/pull/6933)
diff --git a/pennylane/measurements/expval.py b/pennylane/measurements/expval.py
index b6ef68544..d21c58192 100644
--- a/pennylane/measurements/expval.py
+++ b/pennylane/measurements/expval.py
@@ -162,5 +162,4 @@ class ExpectationMP(SampleMeasurement, StateMeasurement):
         Args:
             probabilities (array): the probabilities of collapsing to eigen states
         """
-        eigvals = qml.math.cast_like(self.eigvals(), 1.0)
-        return qml.math.dot(probabilities, eigvals)
+        return qml.math.dot(probabilities, self.eigvals())


Test patch:
diff --git a/tests/measurements/test_expval.py b/tests/measurements/test_expval.py
index 4700fd82f..30f191888 100644
--- a/tests/measurements/test_expval.py
+++ b/tests/measurements/test_expval.py
@@ -45,9 +45,52 @@ def test_expval_identity_nowires_LQ():
         _ = qnode()
 
 
+# pylint: disable=too-many-public-methods
 class TestExpval:
     """Tests for the expval function"""
 
+    @pytest.mark.parametrize("coeffs", [1, 1j, 1 + 1j])
+    def test_process_counts_dtype(self, coeffs):
+        """Test that the return type of the process_counts function is correct"""
+        counts = {"000": 100, "100": 100}
+        wire_order = qml.wires.Wires((0, 1, 2))
+        res = qml.expval(coeffs * qml.Z(1)).process_counts(counts=counts, wire_order=wire_order)
+        assert np.allclose(res, coeffs)
+
+    @pytest.mark.parametrize("coeffs", [1, 1j, 1 + 1j])
+    def test_process_state_dtype(self, coeffs):
+        """Test that the return type of the process_state function is correct"""
+        res = qml.measurements.ExpectationMP(obs=coeffs * qml.Z(0)).process_state(
+            state=[1, 0], wire_order=qml.wires.Wires(0)
+        )
+        assert np.allclose(res, coeffs)
+
+    @pytest.mark.parametrize("coeffs", [1, 1j, 1 + 1j])
+    @pytest.mark.parametrize(
+        "state,expected",
+        [
+            ([[1.0, 0.0], [0.0, 0.0]], 1.0),  # Pure |0⟩ state
+            ([[0.0, 0.0], [0.0, 1.0]], -1.0),  # Pure |1⟩ state
+            ([[0.5, 0.5], [0.5, 0.5]], 0.0),  # Mixed state
+            ([[0.75, 0.0], [0.0, 0.25]], 0.5),  # Another mixed state
+        ],
+    )
+    def test_process_density_matrix_dtype(self, coeffs, state, expected):
+        mp = ExpectationMP(obs=coeffs * qml.PauliZ(0))
+        result = mp.process_density_matrix(state, wire_order=qml.wires.Wires([0]))
+        assert qml.math.allclose(result, expected * coeffs)
+
+    @pytest.mark.parametrize("coeffs", [1, 1j, 1 + 1j])
+    def test_process_samples_dtype(self, coeffs, seed):
+        """Test that the return type of the process_samples function is correct"""
+        shots = 100
+        rng = np.random.default_rng(seed)
+        samples = rng.choice([0, 1], size=(shots, 2)).astype(np.int64)
+        obs = coeffs * qml.PauliZ(0)
+        expected = qml.expval(obs).process_samples(samples, [0, 1])
+        result = ExpectationMP(obs=obs).process_samples(samples, [0, 1])
+        assert qml.math.allclose(result, expected)
+
     @pytest.mark.parametrize("shots", [None, 1111, [1111, 1111]])
     def test_value(self, tol, shots, seed):
         """Test that the expval interface works"""
@@ -358,3 +401,15 @@ class TestExpval:
         res = qml.expval(qml.Z(wire)).process_counts(counts=counts, wire_order=wire_order)
 
         assert np.allclose(res, expected)
+
+
+@pytest.mark.parametrize("coeffs", [1, 1j, 1 + 1j])
+def test_qnode_expval_dtype(coeffs):
+    """System level test to ensure dtype is correctly preserved."""
+
+    @qml.qnode(qml.device("default.qubit"))
+    def circuit(coeffs):
+        return qml.expval(coeffs * qml.PauliZ(0))
+
+    res = circuit(coeffs)
+    assert np.allclose(res, coeffs)