Noisy simulators in Qiskit Runtime

This notebook shows how to set up ibmq_qasm_simulator and map a basic noise model for an IBM Quantum hardware device in Qiskit Runtime, and use this noise model to perform noisy simulations of QuantumCircuits using Sampler and Estimator to study the effects of errors which occur on real devices.

Set up your local development environment

This tutorial requires a Qiskit Runtime service instance to be setup. If you haven’t done so already, follow these steps to set one up.

# load necessary Runtime libraries
from qiskit_ibm_runtime import QiskitRuntimeService, Sampler, Estimator, Session, Options

service = QiskitRuntimeService(channel="ibm_quantum")

Preparing the environment

To demonstrate the routine, we shall proceed with running an example routine. One of the major benefits of using primitives is simplification of binding multiple parameters in parameterized circuits. To check this, here is an example circuit with a controlled P-gate as implemented in the following code. Here, we parametrise the P-gate with a rotation parameter theta. To learn how to create circuits and bind parameters to them by using Qiskit, see the Circuit Basics and Advanced Circuits in Qiskit documentation.

from qiskit.circuit import Parameter
from qiskit import QuantumCircuit

theta = Parameter('theta')

qc = QuantumCircuit(2,1)
qc.x(1)
qc.h(0)
qc.cp(theta,0,1)
qc.h(0)
qc.measure(0,0)

qc.draw('mpl')

The circuit shown by the previous cell is parameterized with the eigenvalue being kicked back into qubit 0 to be measured. The amount of kickback will be determined by the parameter theta. Now in the following cell, we shall define our parameters for our circuit as a list. The parameters here will be from \(0\) to \(2\pi\) divided over 50 evenly spaced points.

import numpy as np

phases = np.linspace(0, 2*np.pi, 50)

# phases need to be expressed as a list of lists in order to work
individual_phases = [[phase] for phase in phases]

Running on the ideal simulator

Set the backend and options to use

First we shall demonstrate a run using an ideal case without any noise_model, optimization_level or resilience_level for both Sampler and Estimator. We shall proceed to setup the options in the following code:

backend = "ibmq_qasm_simulator" # use the simulator

options = Options()
options.simulator.seed_simulator = 42
options.execution.shots = 1000
options.optimization_level = 0 # no optimization
options.resilience_level = 0 # no error mitigation

Run the circuits on Sampler

We shall now sample the circuit to get the result probability distribution using the Sampler primitive to do the same. To learn how to use the Sampler primitive and how to get started using Qiskit Runtime Sessions, you can check this tutorial: Get started with the Sampler primitive.

with Session(service=service, backend=backend):
    sampler = Sampler(options=options)
    job = sampler.run(
        circuits=[qc]*len(phases),
        parameter_values=individual_phases
    )
    result = job.result()

import matplotlib.pyplot as plt

# the probablity of being in the 1 state for each of these values
prob_values = [dist.get(1, 0) for dist in result.quasi_dists]

plt.plot(phases, prob_values, 'o', label='Simulator')
plt.plot(phases, np.sin(phases/2,)**2, label='Theory')
plt.xlabel('Phase')
plt.ylabel('Probability')
plt.legend()

<matplotlib.legend.Legend at 0x7f7fd233b6d0>

Run the circuits on Estimator

To learn how to start a session for Estimator, you may check this tutorial: Get started with the Estimator primitive.

The Estimator will bind single-qubit rotations to get Hamiltonians before it returns expectation values of quantum operators. Therefore, the circuit doesn’t require any measurements. Currently the circuit qc has measurements so we will remove these with remove_final_measurements.

qc_no_meas = qc.remove_final_measurements(inplace=False)
qc_no_meas.draw('mpl')

from qiskit.quantum_info import SparsePauliOp

ZZ = SparsePauliOp.from_list([("ZZ", 1)])
print(f"  > Observable: {ZZ.paulis}")

> Observable: ['ZZ']

With this observable, the expectation value is calculated by the following equation.

\[\langle ZZ\rangle =\langle \psi | ZZ | \psi\rangle=\langle \psi|(|0\rangle\langle 0| -|1\rangle\langle 1|)\otimes(|0\rangle\langle 0| - |1\rangle\langle 1|) |\psi\rangle =|\langle 00|\psi\rangle|^2 - |\langle 01 | \psi\rangle|^2 - |\langle 10 | \psi\rangle|^2 + |\langle 11|\psi\rangle|^2\]

The next cell will implement this as shown.

with Session(service=service, backend=backend):
    estimator = Estimator(options=options)
    job = estimator.run(
        circuits=[qc_no_meas]*len(phases),
        parameter_values=individual_phases,
        observables=[ZZ]*len(phases)
    )
    result = job.result()

exp_values = result.values

plt.plot(phases, exp_values, 'o', label='Simulator')
plt.plot(phases, 2*np.sin(phases/2)**2-1, label='Theory')
plt.xlabel('Phase')
plt.ylabel('Expectation')
plt.legend()

<matplotlib.legend.Legend at 0x7f7fd0ed8820>

Running a noisy simulation

Now we’ll setup our simulator to run a noisy simulation rather than the ideal one. We can pass a custom noise_model to the simulator on Runtime by specifying it in the Options parameter. Here we will try to mimic a real backend and map on the noise_model from a FakeBackend class. The noise model can be extracted from the FakeBackend and passed as a simulator parameter in options. If you want to know more about fake_provider, check Fake Provider in Qiskit documentation.

Since we are trying to mimic a real backend, we can also pass in the coupling_map that the backend topology has and the basis_gates that the backend supports to have a more realistic noisy simulation.

from qiskit.providers.fake_provider import FakeManila
from qiskit_aer.noise import NoiseModel

# Make a noise model
fake_backend = FakeManila()
noise_model = NoiseModel.from_backend(fake_backend)

# Set options to include the noise model
options = Options()
options.simulator = {
    "noise_model": noise_model,
    "basis_gates": fake_backend.configuration().basis_gates,
    "coupling_map": fake_backend.configuration().coupling_map,
    "seed_simulator": 42
}

# Set number of shots, optimization_level and resilience_level
options.execution.shots = 1000
options.optimization_level = 0
options.resilience_level = 0

set_backend() is the syntactic sugar for setting options. The following code is equivalent.

from qiskit.providers.fake_provider import FakeManila

# Make a noise model
fake_backend = FakeManila()

# Set options to include the noise model
options = Options()
options.simulator.set_backend(fake_backend)
options.simulator.seed_simulator = 42

# Set number of shots, optimization_level and resilience_level
options.execution.shots = 1000
options.optimization_level = 0
options.resilience_level = 0

The ibmq_qasm_simulator allows for the activation of the resilience_levels offered by the Qiskit Runtime Service, and use of these levels on simulators is best demonstrated using the noisy simulation as we have described previously.

To see the comparison, we shall define two set of Options. The ibmq_qasm_simulator allows for the activation of the resilience levels offered by Qiskit Runtime, and the use of these levels on simulators is best demonstrated using the noisy simulation that we have built. Here, options is set toresilience level = 0 to represent a normal run without error mitigation, and options with em is set to resilience level = 1 to represent a run with error mitigation enabled.

# Set options to include the noise model with error mitigation
options_with_em = Options()
options_with_em.simulator = {
    "noise_model": noise_model,
    "basis_gates": fake_backend.configuration().basis_gates,
    "coupling_map": fake_backend.configuration().coupling_map,
    "seed_simulator": 42
}

# Set number of shots, optimization_level and resilience_level
options_with_em.execution.shots = 1000
options_with_em.optimization_level = 0 # no optimization
options_with_em.resilience_level = 1 # M3 for Sampler and T-REx for Estimator

When you set the resilience_level to 1, M3 is activated in Sampler. All available resilience level configurations can be found here.

with Session(service=service, backend=backend):
    # include the noise model without M3
    sampler = Sampler(options=options)
    job = sampler.run(
        circuits=[qc]*len(phases),
        parameter_values=individual_phases
    )
    result = job.result()
    prob_values = [1-dist[0] for dist in result.quasi_dists]

    # include the noise model with M3
    sampler = Sampler(options=options_with_em)
    job = sampler.run(
        circuits=[qc]*len(phases),
        parameter_values=individual_phases
    )
    result = job.result()
    prob_values_with_em = [1-dist[0] for dist in result.quasi_dists]

plt.plot(phases, prob_values, 'o', label='Noisy')
plt.plot(phases, prob_values_with_em, 'o', label='Mitigated')
plt.plot(phases, np.sin(phases/2,)**2, label='Theory')
plt.xlabel('Phase')
plt.ylabel('Probability')
plt.legend()

<matplotlib.legend.Legend at 0x7f7fb4230700>

T-REx is triggered in Estimator when the resilience level is set to 1.

with Session(service=service, backend=backend):
    # include the noise model without T-REx
    estimator = Estimator(options=options)
    job = estimator.run(
        circuits=[qc_no_meas]*len(phases),
        parameter_values=individual_phases,
        observables=[ZZ]*len(phases)
    )
    result = job.result()
    exp_values = result.values

    # include the noise model with T-REx
    estimator = Estimator(options=options_with_em)
    job = estimator.run(
        circuits=[qc_no_meas]*len(phases),
        parameter_values=individual_phases,
        observables=[ZZ]*len(phases))
    result = job.result()
    exp_values_with_em = result.values

plt.plot(phases, exp_values, 'o', label='Noisy')
plt.plot(phases, exp_values_with_em, 'o', label='Mitigated')
plt.plot(phases, 2*np.sin(phases/2)**2-1, label='Theory')
plt.xlabel('Phase')
plt.ylabel('Expectation')
plt.legend()

<matplotlib.legend.Legend at 0x7f7f7006ca00>

Resilience levels are currently in beta so sampling overhead and solution quality will vary from circuit to circuit. New features, advanced options and management tools will be released on a rolling basis. You can also play around with higher levels of resilience and explore additional options offered by them. If you want to learn more about activating features like Digital-ZNE, PEC in addition to M3 and T-REx as shown in the previous examples, check out this tutorial: Error suppression and error mitigation with Qiskit Runtime.

import qiskit_ibm_runtime
qiskit_ibm_runtime.version.get_version_info()

'0.8.0'

from qiskit.tools.jupyter import *
%qiskit_version_table

Version Information

Qiskit Software	Version
qiskit-terra	0.22.2
qiskit-aer	0.11.1
qiskit-ibmq-provider	0.19.2
qiskit	0.39.2
qiskit-nature	0.5.0
qiskit-finance	0.3.4
qiskit-optimization	0.4.0
qiskit-machine-learning	0.5.0
System information
Python version	3.8.13
Python compiler	GCC 10.3.0
Python build	default, Mar 25 2022 06:04:10
OS	Linux
CPUs	8
Memory (Gb)	31.211326599121094
Wed Nov 30 02:43:41 2022 UTC