Repository: quantumlib/qsim (language: python)

Issue:
QSimSimulator.simulate_expectation_values for identity operators is wrong when the coefficient is other than one
The QSimSimulator.simulate_expectation_values simply just returns one for the identity operators, which is correct as long as there is no coefficient. When the observable has a coefficient the expectation value should be equal to the coefficient and not one.
The following example calculates the expectation value of an observable in three different ways.
- Using the qsimcirq simulator
- Using cirq simulator
- Calculating the expectation value from the Hamiltonian

```python
import cirq
import qsimcirq


pauli = 'IIIIIII'
num_qubits = 7
moments = 9
density = 0.95
qubits = cirq.LineQubit.range(num_qubits)
cirq_circuit = cirq.testing.random_circuit(qubits=qubits,
                                      n_moments=moments,
                                      op_density=density)

coeff = 0.25
hamiltonian = float(coeff) * cirq.PauliString(
                cirq.I(cirq.LineQubit(i)) for i, p in enumerate(pauli))

qubits = cirq.LineQubit.range(num_qubits)
energy_ham = hamiltonian.expectation_from_state_vector(
    state_vector= cirq.final_state_vector(cirq_circuit), 
    qubit_map={q: i for i, q in enumerate(qubits)},
    atol=0.01)

qsimSim = qsimcirq.QSimSimulator()
qsimcirq_result = qsimSim.simulate_expectation_values(cirq_circuit, hamiltonian) 

cirqSim = cirq.Simulator()
cirq_result = cirqSim.simulate_expectation_values(cirq_circuit, hamiltonian)  

#print(cirq_circuit)
print(f'energy_ham = {energy_ham}')
print(f'qsimcirq_result = {qsimcirq_result}')
print(f'cirq_result = {cirq_result}')
```
which returns:
```
energy_ham = (0.24999994039535522+0j)
qsimcirq_result = [(1+0j)]
cirq_result = [(0.24999994039535522+0j)]
```
As can be seen qsimcirq returns 1 regardless of the coefficient which I think is related to [#565](https://github.com/quantumlib/qsim/issues/565) and [#566](https://github.com/quantumlib/qsim/issues/565). I'm using `0.14.1.dev20220804`.


Discussion / hints:

Fix patch:
diff --git a/pybind_interface/pybind_main.cpp b/pybind_interface/pybind_main.cpp
index b4fe4859..64df5d2a 100644
--- a/pybind_interface/pybind_main.cpp
+++ b/pybind_interface/pybind_main.cpp
@@ -804,8 +804,12 @@ class SimulatorHelper {
       const auto& opsum = std::get<0>(opsum_qubit_count_pair);
       const auto& opsum_qubits = std::get<1>(opsum_qubit_count_pair);
       if (opsum_qubits == 0) {
-        // This represents an identity, which always has EV=1.
-        results.push_back(std::complex<double>(1.0, 0.0));
+        // This represents an identity, which always has EV = weight.
+        auto result = std::complex<double>(0.0, 0.0);
+        for (const auto& str : opsum) {
+          result += str.weight;
+        }
+        results.push_back(result);
       } else if (opsum_qubits <= 6) {
         results.push_back(ExpectationValue<IO, Fuser>(opsum, simulator, state));
       } else {


Test patch:
diff --git a/qsimcirq_tests/qsimcirq_test.py b/qsimcirq_tests/qsimcirq_test.py
index 26e4f7d8..322f9ca2 100644
--- a/qsimcirq_tests/qsimcirq_test.py
+++ b/qsimcirq_tests/qsimcirq_test.py
@@ -523,7 +523,7 @@ def test_expectation_values(mode: str):
     ]
     psum1 = cirq.Z(a) + 3 * cirq.X(b)
     psum2 = cirq.X(a) - 3 * cirq.Z(b)
-    psum3 = cirq.I(a)
+    psum3 = cirq.I(a) * 3
 
     if mode == "noisy":
         circuit.append(NoiseTrigger().on(a))