Note

This page was generated from tutorials/algorithms/03_vqe_simulation_with_noise.ipynb.

VQE with Qiskit Aer Primitives

This notebook demonstrates how to leverage the Qiskit Aer Primitives to run both noiseless and noisy simulations locally. Qiskit Aer not only allows you to define your own custom noise model, but also to easily create a noise model based on the properties of a real quantum device. This notebook will show an example of the latter, to illustrate the general workflow of running algorithms with local noisy simulators.

For further information on the Qiskit Aer noise model, you can consult the Qiskit Aer documentation, as well the tutorial for building noise models.

The algorithm of choice is once again VQE, where the task consists on finding the minimum (ground state) energy of a Hamiltonian. As shown in previous tutorials, VQE takes in a qubit operator as input. Here, you will take a set of Pauli operators that were originally computed by Qiskit Nature for the H2 molecule, using the SparsePauliOp class.

from qiskit.quantum_info import SparsePauliOp

H2_op = SparsePauliOp.from_list(
    [
        ("II", -1.052373245772859),
        ("IZ", 0.39793742484318045),
        ("ZI", -0.39793742484318045),
        ("ZZ", -0.01128010425623538),
        ("XX", 0.18093119978423156),
    ]
)

print(f"Number of qubits: {H2_op.num_qubits}")

Number of qubits: 2

As the above problem is still easily tractable classically, you can use NumPyMinimumEigensolver to compute a reference value to compare the results later.

from qiskit.algorithms import NumPyMinimumEigensolver
from qiskit.opflow import PauliSumOp

numpy_solver = NumPyMinimumEigensolver()
result = numpy_solver.compute_minimum_eigenvalue(operator=PauliSumOp(H2_op))
ref_value = result.eigenvalue.real
print(f"Reference value: {ref_value:.5f}")

Reference value: -1.85728

/tmp/ipykernel_3539/1936825166.py:4: DeprecationWarning: The class ``qiskit.algorithms.minimum_eigen_solvers.numpy_minimum_eigen_solver.NumPyMinimumEigensolver`` is deprecated as of qiskit-terra 0.24.0. It will be removed no earlier than 3 months after the release date. Instead, use the class ``qiskit.algorithms.minimum_eigensolvers.NumPyMinimumEigensolver``. See https://qisk.it/algo_migration for a migration guide.
  numpy_solver = NumPyMinimumEigensolver()
/tmp/ipykernel_3539/1936825166.py:5: DeprecationWarning: The class ``qiskit.opflow.primitive_ops.pauli_sum_op.PauliSumOp`` is deprecated as of qiskit-terra 0.24.0. It will be removed no earlier than 3 months after the release date. For code migration guidelines, visit https://qisk.it/opflow_migration.
  result = numpy_solver.compute_minimum_eigenvalue(operator=PauliSumOp(H2_op))

The following examples will all use the same ansatz and optimizer, defined as follows:

# define ansatz and optimizer
from qiskit.circuit.library import TwoLocal
from qiskit.algorithms.optimizers import SPSA

iterations = 125
ansatz = TwoLocal(rotation_blocks="ry", entanglement_blocks="cz")
spsa = SPSA(maxiter=iterations)

Performance without noise

Let’s first run the VQE on the default Aer simulator without adding noise, with a fixed seed for the run and transpilation to obtain reproducible results. This result should be relatively close to the reference value from the exact computation.

# define callback
# note: Re-run this cell to restart lists before training
counts = []
values = []


def store_intermediate_result(eval_count, parameters, mean, std):
    counts.append(eval_count)
    values.append(mean)

# define Aer Estimator for noiseless statevector simulation
from qiskit.utils import algorithm_globals
from qiskit_aer.primitives import Estimator as AerEstimator

seed = 170
algorithm_globals.random_seed = seed

noiseless_estimator = AerEstimator(
    run_options={"seed": seed, "shots": 1024},
    transpile_options={"seed_transpiler": seed},
)

# instantiate and run VQE
from qiskit.algorithms.minimum_eigensolvers import VQE

vqe = VQE(
    noiseless_estimator, ansatz, optimizer=spsa, callback=store_intermediate_result
)
result = vqe.compute_minimum_eigenvalue(operator=H2_op)

print(f"VQE on Aer qasm simulator (no noise): {result.eigenvalue.real:.5f}")
print(
    f"Delta from reference energy value is {(result.eigenvalue.real - ref_value):.5f}"
)

VQE on Aer qasm simulator (no noise): -1.85160
Delta from reference energy value is 0.00567

You captured the energy values above during the convergence, so you can track the process in the graph below.

import pylab

pylab.rcParams["figure.figsize"] = (12, 4)
pylab.plot(counts, values)
pylab.xlabel("Eval count")
pylab.ylabel("Energy")
pylab.title("Convergence with no noise")

Text(0.5, 1.0, 'Convergence with no noise')

Performance with noise

Now, let’s add noise to our simulation. In particular, you will extract a noise model from a (fake) device. As stated in the introduction, it is also possible to create custom noise models from scratch, but this task is beyond the scope of this notebook.

First, you need to get an actual device backend and from its configuration and properties you can setup a coupling map and a noise model to match the device. Note: You can also use this coupling map as the entanglement map for the variational form if you choose to.

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

# fake providers contain data from real IBM Quantum devices stored in Qiskit Terra,
# and are useful for extracting realistic noise models.
device = FakeVigo()

coupling_map = device.configuration().coupling_map
noise_model = NoiseModel.from_backend(device)

print(noise_model)

NoiseModel:
  Basis gates: ['cx', 'id', 'rz', 'sx', 'x']
  Instructions with noise: ['cx', 'measure', 'x', 'sx', 'id']
  Qubits with noise: [0, 1, 2, 3, 4]
  Specific qubit errors: [('id', (0,)), ('id', (1,)), ('id', (2,)), ('id', (3,)), ('id', (4,)), ('sx', (0,)), ('sx', (1,)), ('sx', (2,)), ('sx', (3,)), ('sx', (4,)), ('x', (0,)), ('x', (1,)), ('x', (2,)), ('x', (3,)), ('x', (4,)), ('cx', (3, 4)), ('cx', (4, 3)), ('cx', (3, 1)), ('cx', (1, 3)), ('cx', (1, 2)), ('cx', (2, 1)), ('cx', (0, 1)), ('cx', (1, 0)), ('measure', (0,)), ('measure', (1,)), ('measure', (2,)), ('measure', (3,)), ('measure', (4,))]

Once the noise model is defined, you can run VQE using an Aer Estimator, where you can pass the noise model to the underlying simulator using the backend_options dictionary. Please note that this simulation will take longer than the noiseless one.

noisy_estimator = AerEstimator(
    backend_options={
        "method": "density_matrix",
        "coupling_map": coupling_map,
        "noise_model": noise_model,
    },
    run_options={"seed": seed, "shots": 1024},
    transpile_options={"seed_transpiler": seed},
)

Instead of defining a new instance of the VQE class, you can now simply assign a new estimator to our previous VQE instance. As the callback method will be re-used, you will also need to re-start the counts and values variables to be able to plot the convergence graph later on.

# re-start callback variables
counts = []
values = []

vqe.estimator = noisy_estimator

result1 = vqe.compute_minimum_eigenvalue(operator=H2_op)

print(f"VQE on Aer qasm simulator (with noise): {result1.eigenvalue.real:.5f}")
print(
    f"Delta from reference energy value is {(result1.eigenvalue.real - ref_value):.5f}"
)

VQE on Aer qasm simulator (with noise): -1.75227
Delta from reference energy value is 0.10501

if counts or values:
    pylab.rcParams["figure.figsize"] = (12, 4)
    pylab.plot(counts, values)
    pylab.xlabel("Eval count")
    pylab.ylabel("Energy")
    pylab.title("Convergence with noise")

Summary

In this tutorial, you compared three calculations for the H2 molecule ground state. First, you produced a reference value using a classical minimum eigensolver. Then, you proceeded to run VQE using the Qiskit Aer Estimator with 1024 shots. Finally, you extracted a noise model from a backend and used it to define a new Estimator for noisy simulations. The results are:

print(f"Reference value: {ref_value:.5f}")
print(f"VQE on Aer qasm simulator (no noise): {result.eigenvalue.real:.5f}")
print(f"VQE on Aer qasm simulator (with noise): {result1.eigenvalue.real:.5f}")

Reference value: -1.85728
VQE on Aer qasm simulator (no noise): -1.85160
VQE on Aer qasm simulator (with noise): -1.75227

You can notice that, while the noiseless simulation’s result is closer to the exact reference value, there is still some difference. This is due to the sampling noise, introduced by limiting the number of shots to 1024. A larger number of shots would decrease this sampling error and close the gap between these two values.

As for the noise introduced by real devices (or simulated noise models), it could be tackled through a wide variety of error mitigation techniques. The Qiskit Runtime Primitives have enabled error mitigation through the resilience_level option. This option is currently available for remote simulators and real backends accessed via the Runtime Primitives, you can consult this tutorial for further information.

import qiskit.tools.jupyter

%qiskit_version_table
%qiskit_copyright

Version Information

Qiskit Software	Version
qiskit-terra	0.24.1
qiskit-aer	0.12.0
qiskit-ibmq-provider	0.20.2
qiskit	0.43.1
qiskit-nature	0.6.2
qiskit-finance	0.3.4
qiskit-optimization	0.5.0
qiskit-machine-learning	0.6.1
System information
Python version	3.8.16
Python compiler	GCC 11.3.0
Python build	default, Jan 11 2023 00:28:51
OS	Linux
CPUs	2
Memory (Gb)	6.7812042236328125
Fri Jun 09 12:46:59 2023 UTC

This code is a part of Qiskit

© Copyright IBM 2017, 2023.

This code is licensed under the Apache License, Version 2.0. You may
obtain a copy of this license in the LICENSE.txt file in the root directory
of this source tree or at http://www.apache.org/licenses/LICENSE-2.0.

Any modifications or derivative works of this code must retain this
copyright notice, and modified files need to carry a notice indicating
that they have been altered from the originals.