Hardware-efficient trial states for variational quantum eigensolvers
7 min readDec 11, 2020
Cahyati Supriyati Sangaji (My Note)
h2_hamiltonian = (-1.0523732) II +
(0.39793742) IZ +
(-0.3979374) ZI +
(-0.0112801) ZZ +
(0.18093119) XXNote that these co-efficients are functions of the interatomic distance between the atoms. A term like IZ is shorthand notation for a tensor product of two-qubit operators -- the identity operator (πΌ) on one qubit and pauli-Z operator (π) on the other qubit.
Installing necessary packages
Before we begin, you will need to install some prerequisites into your environment. Run the cell below to complete these installations. At the end, the cell outputs will be cleared.
!pip install -U qiskit==0.19
!pip install -U qiskit-ibmq-provider==0.7from IPython.display import clear_output
clear_output()
Hardware efficient trial states
from numpy import pi
from qiskit import QuantumCircuit, Aer, execute
from qiskit.visualization import plot_histogramdef prepare_hets_circuit(depth, angle1, angle2):hets_circ = QuantumCircuit(depth)
hets_circ.ry(angle1, 0)
hets_circ.rz(angle1, 0)
hets_circ.ry(angle1, 1)
hets_circ.rz(angle1, 1)for ii in range(depth):
hets_circ.cx(0,1)
hets_circ.ry(angle2,0)
hets_circ.rz(angle2,0)
hets_circ.ry(angle2,1)
hets_circ.rz(angle2,1)
return hets_circhets_circuit = prepare_hets_circuit(2, pi/2, pi/2)
hets_circuit.draw()
Measuring expectation values
Next, we measure expectation values. We will begin by measuring the ZZ expectation value, or β¨ππβ©. We will first create a copy of the hets_circ quantum circuit that we created above, and add measurements to it.
def measure_zz_circuit(given_circuit):
zz_meas = given_circuit.copy()
zz_meas.measure_all()
return zz_measzz_meas = measure_zz_circuit(hets_circuit)
zz_meas.draw()
Next, letβs execute this quantum circuit and see the measurement outcomes.
simulator = Aer.get_backend('qasm_simulator')result = execute(zz_meas, backend = simulator, shots=10000).result()
counts = result.get_counts(zz_meas)plot_histogram(counts)
We can analyze the counts and calculate the β¨ππβ© as follows:
def measure_zz(given_circuit, num_shots = 10000):zz_meas = measure_zz_circuit(given_circuit)
result = execute(zz_meas, backend = simulator, shots = num_shots).result()
counts = result.get_counts(zz_meas)if '00' not in counts:
counts['00'] = 0
if '01' not in counts:
counts['01'] = 0
if '10' not in counts:
counts['10'] = 0
if '11' not in counts:
counts['11'] = 0total_counts = counts['00'] + counts['11'] + counts['01'] + counts['10']
zz = counts['00'] + counts['11'] - counts['01'] - counts['10']
zz = zz / total_counts
return zzzz = measure_zz(hets_circuit)
print("<ZZ> =", str(zz))
Out:
<ZZ> = -0.0036What about β¨ππΌβ© and β¨πΌπβ©? Do these need new circuits?
The answer is no, and they can be computed from the results above.
def measure_zi(given_circuit, num_shots = 10000):
zz_meas = measure_zz_circuit(given_circuit)
result = execute(zz_meas, backend = simulator, shots = num_shots).result()
counts = result.get_counts(zz_meas)if '00' not in counts:
counts['00'] = 0
if '01' not in counts:
counts['01'] = 0
if '10' not in counts:
counts['10'] = 0
if '11' not in counts:
counts['11'] = 0total_counts = counts['00'] + counts['11'] + counts['01'] + counts['10']
zi = counts['00'] - counts['11'] + counts['01'] - counts['10']
zi = zi / total_counts
return zidef measure_iz(given_circuit, num_shots = 10000):
zz_meas = measure_zz_circuit(given_circuit)
result = execute(zz_meas, backend = simulator, shots = num_shots).result()
counts = result.get_counts(zz_meas)if '00' not in counts:
counts['00'] = 0
if '01' not in counts:
counts['01'] = 0
if '10' not in counts:
counts['10'] = 0
if '11' not in counts:
counts['11'] = 0total_counts = counts['00'] + counts['11'] + counts['01'] + counts['10']
iz = counts['00'] - counts['11'] - counts['01'] + counts['10']
iz = iz / total_counts
return izzi = measure_zi(hets_circuit)
print("<ZI> =", str(zi))iz = measure_iz(hets_circuit)
print("<IZ> =", str(iz))
Out:
<ZI> = -0.0072
<IZ> = -0.0046Next, we measure β¨ππβ©
def measure_xx_circuit(given_circuit):
xx_meas = given_circuit.copy()
### WRITE YOUR CODE BETWEEN THESE LINES - START
given_circuit.h(0)
given_circuit.h(1)
xx_meas.measure_all()
### WRITE YOUR CODE BETWEEN THESE LINES - END return xx_measxx_meas = measure_xx_circuit(hets_circuit)
xx_meas.draw()
def measure_xx(given_circuit, num_shots = 10000):
xx_meas = measure_xx_circuit(given_circuit)
result = execute(xx_meas, backend = simulator, shots = num_shots).result()
counts = result.get_counts(xx_meas)if '00' not in counts:
counts['00'] = 0
if '01' not in counts:
counts['01'] = 0
if '10' not in counts:
counts['10'] = 0
if '11' not in counts:
counts['11'] = 0total_counts = counts['00'] + counts['11'] + counts['01'] + counts['10']
xx = counts['00'] + counts['11'] - counts['01'] - counts['10']
xx = xx / total_counts
return xxxx = measure_xx(hets_circuit)
print("<XX> =", str(xx))
Out:
<XX> = -1.0Now we evaluate the energy of the trial state
def get_energy(given_circuit, num_shots = 10000):
zz = measure_zz(given_circuit, num_shots = num_shots)
iz = measure_iz(given_circuit, num_shots = num_shots)
zi = measure_zi(given_circuit, num_shots = num_shots)
xx = measure_xx(given_circuit, num_shots = num_shots)
energy = (-1.0523732)*1 + (0.39793742)*iz + (-0.3979374)*zi + (-0.0112801)*zz + (0.18093119)*xx
return energyenergy = get_energy(hets_circuit)
print("The energy of the trial state is", str(energy))
Out:
The energy of the trial state is -1.0541510196219999Computing gradients
hets_circuit_plus = None
hets_circuit_minus = None### WRITE YOUR CODE BETWEEN THESE LINES - STARTdef hets_circuit_plus(depth, angle1, angle2):hets_circ_plus = QuantumCircuit(depth)
hets_circ_plus.ry(angle1, 0)
hets_circ_plus.rz(angle1, 0)
hets_circ_plus.ry(angle1, 1)
hets_circ_plus.rz(angle1, 1)for ii in range(depth):
hets_circ_plus.cx(0,1)
hets_circ_plus.ry(angle2,0)
hets_circ_plus.rz(angle2,0)
hets_circ_plus.ry(angle2,1)
hets_circ_plus.rz(angle2,1)
return hets_circ_plusdef hets_circuit_minus(depth, angle1, angle2):hets_circ_minus = QuantumCircuit(depth)
hets_circ_minus.ry(angle1, 0)
hets_circ_minus.rz(angle1, 0)
hets_circ_minus.ry(angle1, 1)
hets_circ_minus.rz(angle1, 1)for ii in range(depth):
hets_circ_minus.cx(0,1)
hets_circ_minus.ry(angle2,0)
hets_circ_minus.rz(angle2,0)
hets_circ_minus.ry(angle2,1)
hets_circ_minus.rz(angle2,1)
return hets_circ_minushets_circuit_plus = hets_circuit_plus(depth=2,angle1 = pi/2 + 0.1*pi/2, angle2 = pi/2)
hets_circuit_minus = hets_circuit_minus(depth=2,angle1 = pi/2 - 0.1*pi/2, angle2 = pi/2)
print(hets_circuit_plus)
print(hets_circuit_minus)### WRITE YOUR CODE BETWEEN THESE LINES - END
energy_plus = get_energy(hets_circuit_plus, num_shots=100000)
energy_minus = get_energy(hets_circuit_minus, num_shots=100000)print(energy_plus, energy_minus)
Out:
-0.9600031637752001 -1.1304885633485997As you can see, one of these is certainly lower energy than the other, and is also lower energy than the case when angle1 = pi/2. This is a suitable next point for our iteration of a variational eigensolver.
Bonus 1
While this is not graded, explore whether the decision above would be easy if your execution ran different numbers of shots. In particular, measure energy_plus and energy_minus again with 100, 1000 and 10000 shots to explore how easy or difficult this decision gets with each one.
energy_plus_100, energy_plus_1000, energy_plus_10000 = 0, 0, 0
energy_minus_100, energy_minus_1000, energy_minus_10000 = 0, 0, 0### WRITE YOUR CODE BETWEEN THESE LINES - STARTdef hets_circuit_plus_100(depth, angle1, angle2):hets_circ_plus_100 = QuantumCircuit(depth)
hets_circ_plus_100.ry(angle1, 0)
hets_circ_plus_100.rz(angle1, 0)
hets_circ_plus_100.ry(angle1, 1)
hets_circ_plus_100.rz(angle1, 1)for ii in range(depth):
hets_circ_plus_100.cx(0,1)
hets_circ_plus_100.ry(angle2,0)
hets_circ_plus_100.rz(angle2,0)
hets_circ_plus_100.ry(angle2,1)
hets_circ_plus_100.rz(angle2,1)
return hets_circ_plus_100def hets_circuit_plus_1000(depth, angle1, angle2):hets_circ_plus_1000 = QuantumCircuit(depth)
hets_circ_plus_1000.ry(angle1, 0)
hets_circ_plus_1000.rz(angle1, 0)
hets_circ_plus_1000.ry(angle1, 1)
hets_circ_plus_1000.rz(angle1, 1)for ii in range(depth):
hets_circ_plus_1000.cx(0,1)
hets_circ_plus_1000.ry(angle2,0)
hets_circ_plus_1000.rz(angle2,0)
hets_circ_plus_1000.ry(angle2,1)
hets_circ_plus_1000.rz(angle2,1)
return hets_circ_plus_1000def hets_circuit_plus_10000(depth, angle1, angle2):hets_circ_plus_10000 = QuantumCircuit(depth)
hets_circ_plus_10000.ry(angle1, 0)
hets_circ_plus_10000.rz(angle1, 0)
hets_circ_plus_10000.ry(angle1, 1)
hets_circ_plus_10000.rz(angle1, 1)for ii in range(depth):
hets_circ_plus_10000.cx(0,1)
hets_circ_plus_10000.ry(angle2,0)
hets_circ_plus_10000.rz(angle2,0)
hets_circ_plus_10000.ry(angle2,1)
hets_circ_plus_10000.rz(angle2,1)
return hets_circ_plus_10000def hets_circuit_minus_100(depth, angle1, angle2):hets_circ_minus_100 = QuantumCircuit(depth)
hets_circ_minus_100.ry(angle1, 0)
hets_circ_minus_100.rz(angle1, 0)
hets_circ_minus_100.ry(angle1, 1)
hets_circ_minus_100.rz(angle1, 1)for ii in range(depth):
hets_circ_minus_100.cx(0,1)
hets_circ_minus_100.ry(angle2,0)
hets_circ_minus_100.rz(angle2,0)
hets_circ_minus_100.ry(angle2,1)
hets_circ_minus_100.rz(angle2,1)
return hets_circ_minus_100def hets_circuit_minus_1000(depth, angle1, angle2):hets_circ_minus_1000 = QuantumCircuit(depth)
hets_circ_minus_1000.ry(angle1, 0)
hets_circ_minus_1000.rz(angle1, 0)
hets_circ_minus_1000.ry(angle1, 1)
hets_circ_minus_1000.rz(angle1, 1)for ii in range(depth):
hets_circ_minus_1000.cx(0,1)
hets_circ_minus_1000.ry(angle2,0)
hets_circ_minus_1000.rz(angle2,0)
hets_circ_minus_1000.ry(angle2,1)
hets_circ_minus_1000.rz(angle2,1)
return hets_circ_minus_1000def hets_circuit_minus_10000(depth, angle1, angle2):hets_circ_minus_10000 = QuantumCircuit(depth)
hets_circ_minus_10000.ry(angle1, 0)
hets_circ_minus_10000.rz(angle1, 0)
hets_circ_minus_10000.ry(angle1, 1)
hets_circ_minus_10000.rz(angle1, 1)for ii in range(depth):
hets_circ_minus_10000.cx(0,1)
hets_circ_minus_10000.ry(angle2,0)
hets_circ_minus_10000.rz(angle2,0)
hets_circ_minus_10000.ry(angle2,1)
hets_circ_minus_10000.rz(angle2,1)
return hets_circ_minus_10000hets_circuit_plus_100 = hets_circuit_plus_100(depth=2,angle1 = pi/2 + 0.1*pi/2, angle2 = pi/2)
hets_circuit_plus_1000 = hets_circuit_plus_1000(depth=2,angle1 = pi/2 + 0.1*pi/2, angle2 = pi/2)
hets_circuit_plus_10000 = hets_circuit_plus_10000(depth=2,angle1 = pi/2 + 0.1*pi/2, angle2 = pi/2)hets_circuit_minus_100 = hets_circuit_minus_100(depth=2,angle1 = pi/2 - 0.1*pi/2, angle2 = pi/2)
hets_circuit_minus_1000 = hets_circuit_minus_1000(depth=2,angle1 = pi/2 - 0.1*pi/2, angle2 = pi/2)
hets_circuit_minus_10000 = hets_circuit_minus_10000(depth=2,angle1 = pi/2 - 0.1*pi/2, angle2 = pi/2)print(hets_circuit_plus_100)
print(hets_circuit_plus_1000)
print(hets_circuit_plus_10000)
print(hets_circuit_minus_100)
print(hets_circuit_minus_1000)
print(hets_circuit_minus_10000)### WRITE YOUR CODE BETWEEN THESE LINES - END
energy_plus_100 = get_energy(hets_circuit_plus_100, num_shots=100)
energy_plus_1000 = get_energy(hets_circuit_plus_100, num_shots=1000)
energy_plus_10000 = get_energy(hets_circuit_plus_100, num_shots=10000)energy_minus_100 = get_energy(hets_circuit_minus_100, num_shots=100)
energy_minus_1000 = get_energy(hets_circuit_minus_100, num_shots=1000)
energy_minus_10000 = get_energy(hets_circuit_minus_100, num_shots=10000)print(energy_plus_100, energy_minus_100, "difference = ", energy_minus_100 - energy_plus_100)
print(energy_plus_1000, energy_minus_1000, "difference = ", energy_minus_1000 - energy_plus_1000)
print(energy_plus_10000, energy_minus_10000, "difference = ", energy_minus_10000 - energy_plus_10000)
Out:
-0.9122630381999999 -1.1198517056 difference = -0.20758866740000004
-0.4401202818199999 -0.45621825961999984 difference = -0.01609797779999994
-0.96466067423 -1.1187116705179998 difference = -0.1540509962879998Bonus 2
While this is not graded, diagonalize the Hamiltonian by writing down the matrices for the Pauli operators I, X and Z, and find the exact ground state energy.
### WRITE YOUR CODE BETWEEN THESE LINES - START
import numpy as np
from qiskit.quantum_info.operators import Operator, Pauli
h = [[-1.0523732, Pauli(label='II')],
[0.39793742, Pauli(label='IZ')],
[-0.3979374, Pauli(label='ZI')],
[-0.0112801, Pauli(label='ZZ')],
[0.18093119, Pauli(label='XX')]
]
H_m = []
for hi in h:
H_m.append(hi[0]*Operator(hi[1]))
np.linalg.eigvalsh((H_m[0]+H_m[1]+H_m[2]+H_m[3]+H_m[4]).data)from qiskit.aqua.algorithms import NumPyEigensolver
from qiskit.aqua.operators import WeightedPauliOperatorH = WeightedPauliOperator(h)
print(H.print_details())ee = NumPyEigensolver(H)result = ee.run()
ref = result['eigenvalues']
print(ref)### WRITE YOUR CODE BETWEEN THESE LINES - END
Out:
II (-1.0523732+0j)
IZ (0.39793742+0j)
ZI (-0.3979374+0j)
ZZ (-0.0112801+0j)
XX (0.18093119+0j)
[-1.85727496-1.35566063e-16j]References:
Qiskit (Global Summer School), Introduction to Quantum Computing and Quantum Hardware β Lab 9.
