Real Hardware: N=2 Case

Authors: Kyle Godbey, Alexandra Semposki

Maintainer: Kyle Godbey

Here we’ll try running \(N=2\) deuteron case from the VQE chapter on real quantum hardware. This will, of course, be noisy just as the last chapter, but now our noise model is a little more true to life!

%matplotlib inline

import matplotlib.pyplot as plt
from pennylane import numpy as np
import pennylane as qml
import warnings
warnings.filterwarnings('ignore')

import mitiq as mq
from mitiq.zne.scaling import fold_global
from mitiq.zne.inference import RichardsonFactory
from pennylane.transforms import mitigate_with_zne
from qiskit.providers.aer import AerSimulator
import qiskit.providers.aer.noise as noise
from qiskit import IBMQ
from qiskit.providers.ibmq import RunnerResult
from pennylane_qiskit import upload_vqe_runner, vqe_runner

VQE for N=2 (on IBM Nairobi)

Here we’ll pull from some of the VQE notebook to define the Hamiltonian and circuit again. If you would like to revisit the VQE process, go to this page.

#define the Hamiltonian we're using again, for N=2
coeffs = [5.906709, 0.218291, -6.125, -2.143304, -2.143304]
obs = [qml.Identity(0), qml.PauliZ(0), qml.PauliZ(1), qml.PauliX(0) @ qml.PauliX(1), qml.PauliY(0) @ qml.PauliY(1)]

H = qml.Hamiltonian(coeffs, obs)

print(H)

  (-6.125) [Z1]
+ (0.218291) [Z0]
+ (5.906709) [I0]
+ (-2.143304) [X0 X1]
+ (-2.143304) [Y0 Y1]

Now we’ll set up our qiskit runtime job. For this, be sure to add your own IBMQ API token in token to authenticate with the IBMQ system.

token = "XXX"

IBMQ.enable_account(token)

program_id = upload_vqe_runner(hub="ibm-q", group="open", project="main")

Now we define some initial parameter and convergence info, as before, but we also need to define our job that wil be submitted to the qiskit runtime.

#now we'll define the circuit we want to use with this Hamiltonian
def circuit(params):
    t0 = params
    qml.PauliX(wires=0)
    qml.RY(t0, wires=1)
    qml.CNOT(wires=[1,0])
   # return qml.expval(H)

#parameter array
init_params = np.array([2.5,])

#convergence information and step size
max_iterations = 120
conv_tol = 1e-06
step_size = 0.01
shots = 10000

job = vqe_runner(
    program_id=program_id,
    backend="ibm_nairobi",
    hamiltonian=H,
    ansatz=circuit,
    x0=init_params,
    shots=shots,
    optimizer="SPSA",
    optimizer_config={"maxiter": 50},
    kwargs={"hub": "ibm-q", "group": "open", "project": "main"},
)

# Get runtime job result.
result = job.result()

print(result)

     fun: -0.7380420742000001
 message: 'Optimization terminated successfully.'
    nfev: 150
     nit: 50
 success: True
       x: array([2.70905085])

Phew, that took a long time! But that’s the curse of running on real hardware with other real humans interested in running their circuits too. The value returned is also not great.. let’s see if we can improve that!

ZNE for N=2

#set up the extrapolation step and scale factors to use
extrapolate = RichardsonFactory.extrapolate
scale_factors = [1, 2, 3, 4]

params = result['x']
unmitig_energy = result['fun']

dev_ibm = qml.device('qiskit.ibmq', wires=2, backend='ibm_nairobi', ibmqx_token=token)

def circuit(params, wires):
    t0 = params[0]
    qml.PauliX(wires=0)
    qml.RY(t0, wires=1)
    qml.CNOT(wires=[1,0])
    return qml.expval(H)

ibm_qnode = qml.QNode(circuit, dev_ibm)

#use ZNE to mitigate error
mitigated_qnode = mitigate_with_zne(scale_factors, fold_global, extrapolate)(ibm_qnode)
zne_result = mitigated_qnode(params,wires=2, shots=2**14)

#print the result of the ground state energy
print('N=2 ground state energy result with ZNE: {} MeV'.format(zne_result))

N=2 ground state energy result with ZNE: -1.4298696529541 MeV

The result is still pretty bad… but again, that’s one side effect of the current generation of hardware and our implementation on it. To do this better you could enable error mitigation during the VQE iteration, at the cost of additional time. That can be left as an exercise for the reader ;)