# Evaluating expectation values of operators in Qiskit

I'm wondering how in Qiskit one can calculate the expectation value of an operator given as a WeightedPauli (or, at least, of a single Pauli operator...) in a certain state (given as a QuantumCircuit object ⁠— meaning that the actual state is the result of the action of this circuit on the computational basis state). I would like the inputs of such a procedure to be floats, not Parameters (it is an essential requirement — I'm using an external library to form the circuit for each set of parameters, and then converting it gate-by-gate to Qiskit format).

This would be useful if, say, we wanted to manually implement VQE, and for that needed a function calculating the expectation value of the Hamiltonian on a quantum computer. More importantly, we would need this for implementing generalizations of VQE, such as subspace search.

I guess, PauliBasisChange may be involved...

• – Simon Crane May 22 at 16:17
• @Simon, to my understanding, this is an unphysical method. I'm looking for smth that would work on an anctual device, based on sampling. – mavzolej May 22 at 16:22
• The SnapshotExpectationValue is developed to be fast. It computes the expectation value via matrix multiplication, not using shots. – Cryoris May 27 at 16:48

The operators in Qiskit Aqua allow the evaluation of expectation values both exactly (via matrix multiplication) or on shot-based sampling (closer to real quantum computers). The basic principle is the same both times, it only differs in how the expectation value is evaluated in the end.

First, you need to define the operator $$O$$ you're interested in and the state $$|\psi\rangle$$ with respect to which you want to compute the expecation value. So we're looking for $$E = \langle\psi|O|\psi\rangle.$$ In the code below we have $$O$$ = op and $$|\psi\rangle$$ = psi. See also there for your use-case of a WeightedPauliOperator.

# you can define your operator as circuit
circuit = QuantumCircuit(2)
circuit.z(0)
circuit.z(1)
op = CircuitOp(circuit)  # and convert to an operator

# or if you have a WeightedPauliOperator, do
op = weighted_pauli_op.to_opflow()

# but here we'll use the H2-molecule Hamiltonian
from qiskit.aqua.operators import X, Y, Z, I
op =  (-1.0523732 * I^I) + (0.39793742 * I^Z) + (-0.3979374 * Z^I) \
+ (-0.0112801 * Z^Z) + (0.18093119 * X^X)

# define the state you w.r.t. which you want the expectation value
psi = QuantumCircuit(2)
psi.x(0)
psi.x(1)

# convert to a state
psi = CircuitStateFn(circuit)


There are now different ways to evaluate the expectation value. The straightforward, "mathematical", approach would be to take the adjoint of $$|\psi\rangle$$ (which is $$\langle\psi|$$) and multiply with $$O$$ and then $$|\psi\rangle$$ to get the expectation. You can actually do exactly this in Qiskit:

# easy expectation value, use for small systems only!


to get

Exact: -1.0636533199999998


This is only suitable for small systems though.

To use the simulators, and the also get the shot-based result, you can use the PauliExpectation (shots), AerPauliExpectation (exact) or MatrixExpectation (exact). Here's how to do it:

from qiskit import Aer
from qiskit.aqua import QuantumInstance
from qiskit.aqua.operators import PauliExpectation, CircuitSampler, StateFn

# define your backend or quantum instance (where you can add settings)
backend = Aer.get_backend('qasm_simulator')
q_instance = QuantumInstance(backend, shots=1024)

# define the state to sample
measurable_expression = StateFn(op, is_measurement=True).compose(psi)

# convert to expectation value
expectation = PauliExpectation().convert(measurable_expression)

# get state sampler (you can also pass the backend directly)
sampler = CircuitSampler(q_instance).convert(expectation)

# evaluate
print('Sampled:', sampler.eval().real)


which yields

Sampled: -1.0530518430859401


This result varies if you execute multiple times.

For comparison, here the other methods to evaluate the expecation value

expectation = AerPauliExpectation().convert(measurable_expression)
sampler = CircuitSampler(backend).convert(expectation)
print('Snapshot:', sampler.eval().real)

expectation = MatrixExpectation().convert(measurable_expression)
sampler = CircuitSampler(backend).convert(expectation)
print('Matrix:', sampler.eval().real)


which produces

Snapshot: -1.06365328
Matrix: -1.06365328


I hope that clarifies how to compute the expectation value!

• Thanks! A couple of questions regarding the sampler solution: 1) Will it work if I replace the qasm backend with a real device? 2) Where do I specify the number of samples? – mavzolej May 28 at 3:23
• 1) Yes it should, it can use the QuantumInstance like the other Aqua algorithms use too. 2) By passing a QuantumInstance, I updated the answer above. Could you accept the answer, if your question was answered? :) It helps so others can see that this is resolved and a solution was found. – Cryoris May 28 at 7:50
• What kinds of optimizations does Qikist do while using this method? For example, does it split the WeightedPauliSum into commuting sets, in order to reduce the number of measurements? – mavzolej May 28 at 12:54
• Following my previous comment - are there any circuit optimizations performed using the transpiler in this code, or should I first do that manually? – mavzolej May 28 at 13:33
• Also, the line q_instance = QuantumInstance(backend, shots = shots) generates the following warning, if the real hardware is used: WARNING - The skip Qobj validation does not work for IBMQ provider. Disable it.. Do you know why, and how to avoid it? – mavzolej May 28 at 13:39