Improving the results of quantum phase estimation
Using Fire Opal's error suppression to enhance the performance of quantum phase estimation on real hardware
Fire Opal reduces errors and increases algorithmic performance when running quantum algorithms on real quantum hardware. This is done by optimizing the quantum circuit using techniques from quantum control theory to mitigate the noise on NISQ devices and beyond. Fire Opal applies error suppression and error mitigation across various supported backends without requiring any user configuration—everything happens automatically.
The benefits of Fire Opal are broadly applicable and can be used to improve the performance of any type of quantum algorithm. In this application note, we demonstrate the algorithmic performance improvements of Fire Opal when applied to the quantum algorithm of quantum phase estimation, which is a commonly used subroutine that serves as a building block for various algorithms.
This application note covers the following:
- An introduction to quantum phase estimation and how it is incorporated into quantum simulation algorithms.
- Definition of a simple instance of quantum phase estimation.
- Execution of the circuits on an ideal simulator and on IBM hardware with and without Fire Opal.
- Analysis of the results and conclusions.
The example in this notebook is inspired by the Qiskit documentation, and we compare the performance of Fire Opal and IBM default.
Some cells in this notebook require an account with IBM-Q to execute correctly. If you want to run them, please go to the IBM-Q experience to set up an account.
1. Introduction to quantum phase estimation
In this section, we introduce the problem of quantum phase estimation and explain its relevance in quantum chemistry algorithms.
1.1 Definition of the problem of quantum phase estimation
The goal of quantum phase estimation is to provide an estimate on the eigenvalues of a unitary operation $U$. More specifically, quantum phase estimation is defined by the following task. For an $n$ qubit unitary $U$, we are given:
- An $n$ qubit quantum circuit that implements $U$. Furthermore, we need to be able to apply controlled $U$ operations,
- An eigenstate $\vert \varphi \rangle$ of $U$ according to the eigenvalue equation $U \vert \varphi \rangle = e^{2 \pi i \varphi} \vert \varphi \rangle$, and
- $m$ additional qubits, which we refer to as the "counting qubits".
The algorithm of quantum phase estimation then computes the phase $\varphi$ up to an error $\epsilon$ by encoding it into the counting qubit register. Through this procedure, the error $\epsilon$ becomes exponentially small in the number of counting qubits, $\epsilon = O(1 / 2^m)$. The quantum circuit of quantum phase estimation looks as follows:
From this circuit, we can see that we need to call the quantum circuit implementing $U$ $O(1 / \epsilon)$ many times. We will discuss later how to encode the phase $\varphi$ in the final quantum state and how we can retrieve it from measuring the counting qubits. For the mathematical details on how this quantum circuit works, we refer to the Qiskit documentation.
If we only have access to an input state $\vert \psi \rangle$, which is not an eigenstate of $U$, we can still perform quantum phase estimation. In this case, the algorithm projects the input qubits to a specific eigenstate $\vert \varphi_j \rangle$ of $U$ with a probability $p_j = \vert \langle \psi \vert \varphi_j \rangle \vert^2$ and computes the corresponding phase $\varphi_j$. Therefore, in general quantum phase estimation samples the eigenvalues of $U$ from a probability distribution defined by the input state. By repeating the algorithm many times, we thus obtain an estimate of all eigenvalues of $U$ with nonzero $p_j$. Notice that with a large number of repititions it becomes exponentially likely to sample a given eigenvalue at least once. For example, even for a small overlap of say $p_j = 0.01$, a shot count of $1,000$ suffices to recover $\varphi_j$ with a high probability.
1.2 Relevance of quantum phase estimation
Quantum phase estimation is a subroutine in many important quantum algorithms including the simulation of chemical materials such as drugs or fertilizers and Shor's algorithm for RSA decryption. Here, we focus on its incorporation in quantum chemistry simulation, which is proposed as a promising application of quantum computers. Using quantum chemistry algorithms, we could potentially find novel medicine or improve the reaction outcomes in industrial production of chemicals.
Many physical properties of chemical materials are derived from the eigenvalues of their corresponding Hamiltonian $H$. Specifically, the ground state energy (smallest eigenvalue) is an important figure of merit. By encoding $H$ into a unitary $U$, we can use quantum phase estimation to find its eigenvalues. There are various ways to perform this encoding, for example:
- Hamiltonian simulation through Trotterisation,
- Linear combination of unitaries, and
- Qubitisation through a block encoding.
Decoding the results from quantum phase estimation for the corresponding unitary $U$ allows us to conclude what are the eigenvalues of $H$.
1.3 Encoding of the eigenvalues
By using phase kickback, quantum phase estimation provides a binary encoding of the phase $\varphi$ in the final state of the counting qubit register. For example, consider the case of $\varphi = 0.25$ on a quantum computer with $m = 3$ counting qubits. This quantum computer has $2^m = 2^3 = 8$ different states, which are labeled by all possible 3-digit bitstrings $\vert 000 \rangle, \vert 001 \rangle, \dots,\vert 111 \rangle$. Each of these bitstrings encodes a number $q$ between $0$ and $7$ – in the exact same way as on a classical computer. That is, $\vert 000 \rangle \equiv 0, \vert 001 \rangle \equiv 1, \dots$. Dividing each of these numbers by $8$, that is computing $q / 2^m$, results in numbers between $0$ and $1$, which encode the phase (up to a multiple of $2 \pi$). In our case, $\varphi = 0.25 = 2 / 8$, and hence we expect to measure the state $\vert 010 \rangle \equiv 2$ after quantum phase estimation. The above discussion is summarized in the following table:
| Quantum state (bitstring) | Encoded number $q$ | Represented phase $\varphi = q / 2^m$ |
|---|---|---|
| $\vert 000 \rangle$ | $0$ | $0$ |
| $\vert 001 \rangle$ | $1$ | $0.125$ |
| $\mathbf{\color{#680CE9}{\vert 010 \rangle}}$ | $\mathbf{\color{#680CE9}{2}}$ | $\mathbf{\color{#680CE9}{0.25}}$ |
| $\vert 011 \rangle$ | $3$ | $0.375$ |
| $\vert 100 \rangle$ | $4$ | $0.5$ |
| $\vert 101 \rangle$ | $5$ | $0.625$ |
| $\vert 110 \rangle$ | $6$ | $0.75$ |
| $\vert 111 \rangle$ | $7$ | $0.875$ |
But what happens if the actual phase was not an inverse multiple of $8$? For instance, imagine a scenario where $\varphi = 0.3$. In this case, there is no $q / 2^m$, which represents $\varphi$ exactly. However, we can still perform quantum phase estimation to obtain an estimate of the phase. The phase $\varphi = 0.3$ lies between $0.25$ and $0.375$ but is closer to $0.25$. For this reason, after performing the quantum phase estimation algorithm, we will obtain a quantum state which has its maximum probability amplitude in $\vert 010 \rangle$ encoding a phase of $0.25$. Nevertheless, as opposed to the previous example of $\varphi = 0.25$, we also receive a nonzero probability amplitude in other quantum states. In general, quantum phase estimation will give us an estimate for the actual phase by providing a lower and upper bound: One bound is always given by the highest amplitude bitstring. The other bound is given by the higher amplitude state next to it. In the case of $\varphi = 0.3$, the two bitstrings next to the highest amplitude one are given by the state $\vert 001 \rangle \equiv 0.125$ (left) and $\vert 011 \rangle \equiv 0.375$ (right). Thus, we will see that the state $\vert 011 \rangle$ has a higher amplitude than the state $\vert 001 \rangle$ as the desired phase is between $0.25$ and $0.375$. If we want to increase the accuracy of the estimate, we have to increase the number of counting qubits $m$. See the following figures for an illustration of the ideal probability distributions:
2. Imports and initialization
The following section sets up the necessary imports and helper functions, which will be used in creating and evaluating the circuits.
# Import packages.
import math
import fireopal
import numpy as np
import matplotlib.pyplot as plt
import json
import qiskit
from qiskit_aer import AerSimulator
from qiskit_ibm_runtime import QiskitRuntimeService
from qiskit_ibm_runtime import SamplerV2 as Sampler
from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager
from qctrlvisualizer import get_qctrl_style
plt.style.use(get_qctrl_style())
# Define the functions to create the quantum phase estimation circuits.
def power_of_two_q(number):
"""
Ask whether a number is a power of 2.
"""
return math.ceil(np.log2(number)) == math.floor(np.log2(number))
def inverse_quantum_fourier_transform(quantum_circuit, number_of_qubits):
"""
Apply an inverse Quantum Fourier Transform the first `number_of_qubits` qubits in the
`quantum_circuit`.
"""
for qubit in range(number_of_qubits // 2):
quantum_circuit.swap(qubit, number_of_qubits - qubit - 1)
for j in range(number_of_qubits):
for m in range(j):
quantum_circuit.cp(-math.pi / float(2 ** (j - m)), m, j)
quantum_circuit.h(j)
def create_circuit(number_of_counting_qubits, phase):
"""
Create a quantum circuit for quantum phase estimation. The unitary, whose phase is evaluated
is a single qubit phase gate. Its phase can be specified by the input parameter `phase`.
"""
# Initialize the quantum circuit.
quantum_circuit = qiskit.QuantumCircuit(
number_of_counting_qubits + 1, number_of_counting_qubits
)
quantum_circuit.x(number_of_counting_qubits)
for qubit in range(number_of_counting_qubits):
quantum_circuit.h(qubit)
# Add the controlled phase gates to the circuit.
repetitions = 1
for counting_qubit in range(number_of_counting_qubits):
for _ in range(repetitions):
quantum_circuit.cp(phase, counting_qubit, number_of_counting_qubits)
repetitions *= 2
quantum_circuit.barrier()
# Apply the inverse Quantum Fourier Transform.
inverse_quantum_fourier_transform(quantum_circuit, number_of_counting_qubits)
# Measure the ancilla register.
quantum_circuit.barrier()
for qubit in range(number_of_counting_qubits):
quantum_circuit.measure(qubit, qubit)
return quantum_circuit
def bitstring_count_to_probabilities(data, shot_count, number_of_counting_qubits):
"""
Convenience function to manipulate the measured and simulated results.
This function turns an unsorted dictionary of bitstring counts into a sorted dictionary
of probabilities. It also turns bitstrings from base 0 to base 2 and adds bitstrings,
with zero count that do not appear in the original dictionary.
"""
# Check if the base of the bitstrings is 0 and turn it into base 2 in this case.
if any("x" in key for key in data.keys()):
base = 0
else:
base = 2
# Turn the bitstring counts into probabilities.
data = {
format(int(bitstring, base), f"0{number_of_counting_qubits}b"): bitstring_count
/ shot_count
for bitstring, bitstring_count in data.items()
}
# Add bitstrings with zero count that do not appear in the dictionary.
all_bitstrings = [
format(i, f"0{number_of_counting_qubits}b")
for i in range(2**number_of_counting_qubits)
]
for bitstring in all_bitstrings:
if bitstring not in data:
data[bitstring] = 0.0
# Sort the dictionary.
sorted_data = dict(sorted(data.items()))
return sorted_data3. Setting up the problem
The following section involves the problem setup and creation of our circuits.
3.1 Setting up our IBM account
We will need to set up our IBM account credentials and session in order to run our circuits on real hardware.
# These are the properties for the publicly available provider for IBM backends.
# If you have access to a private provider and wish to use it, replace these values.
hub = "ibm-q"
group = "open"
project = "main"
token = "YOUR_TOKEN_HERE"
credentials = fireopal.credentials.make_credentials_for_ibmq(
token=token, hub=hub, group=group, project=project
)
instance = hub + "/" + group + "/" + project
service = QiskitRuntimeService(token=token, instance=instance, channel="ibm_quantum")
# Replace "DESIRED_BACKEND_NAME_HERE" with the name of the device you wish to use.
# Run fireopal.show_supported_devices(credentials) to get a list of suppported backend names.
backend_name = "DESIRED_BACKEND_NAME_HERE"3.2 Creating the circuits
Here, we choose 10 different phases and run the circuit for each phase twice. In this way, we obtain 20 data points per number of counting qubits. We exclude phases that are an inverse power of two to avoid exact results. This is a more realistic scenario closer to a real quantum chemistry application.
# Set parameters.
shot_count = 4096
phases = [1 / number for number in np.arange(1, 15) if power_of_two_q(number) is False]
experiments_per_phase = 2
# Set up the quantum circuits.
# 3 counting qubits.
circuits_quantum_phase_estimation_3 = [
create_circuit(3, phase) for phase in phases
] * experiments_per_phase
print(
"quantum phase estimation circuit with m=3 counting qubits and a phase of 𝜑=0.2/2π:"
)
display(circuits_quantum_phase_estimation_3[1].draw(fold=-1))
# 4 counting qubits.
circuits_quantum_phase_estimation_4 = [
create_circuit(4, phase) for phase in phases
] * experiments_per_phase
print(
"quantum phase estimation circuit with m=4 counting qubits and a phase of 𝜑=0.2/2π:"
)
display(circuits_quantum_phase_estimation_4[1].draw(fold=-1))
# 5 counting qubits.
circuits_quantum_phase_estimation_5 = [
create_circuit(5, phase) for phase in phases
] * experiments_per_phase
print(
"quantum phase estimation circuit with m=5 counting qubits and a phase of 𝜑=0.2/2π:"
)
display(circuits_quantum_phase_estimation_5[1].draw(fold=-1))
# 6 counting qubits.
circuits_quantum_phase_estimation_6 = [
create_circuit(6, phase) for phase in phases
] * experiments_per_phase
print(
"quantum phase estimation circuit with m=6 counting qubits and a phase of 𝜑=0.2/2π:"
)
display(circuits_quantum_phase_estimation_6[1].draw(fold=-1))quantum phase estimation circuit with m=3 counting qubits and a phase of 𝜑=0.2/2π:
┌───┐ ░ ┌───┐ ░ ┌─┐
q_0: ┤ H ├─■──────────────────────────────────────────────────────────────░──X─┤ H ├─■──────────────■────────────────────────░─┤M├──────
├───┤ │ ░ │ └───┘ │P(-π/2) ┌───┐ │ ░ └╥┘┌─┐
q_1: ┤ H ├─┼────────■────────■────────────────────────────────────────────░──┼───────■────────┤ H ├─┼─────────■──────────────░──╫─┤M├───
├───┤ │ │ │ ░ │ └───┘ │P(-π/4) │P(-π/2) ┌───┐ ░ ║ └╥┘┌─┐
q_2: ┤ H ├─┼────────┼────────┼────────■────────■────────■────────■────────░──X──────────────────────■─────────■────────┤ H ├─░──╫──╫─┤M├
├───┤ │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) ░ └───┘ ░ ║ ║ └╥┘
q_3: ┤ X ├─■────────■────────■────────■────────■────────■────────■────────░──────────────────────────────────────────────────░──╫──╫──╫─
└───┘ ░ ░ ║ ║ ║
c: 3/═══════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════╩══╩══╩═
0 1 2 quantum phase estimation circuit with m=4 counting qubits and a phase of 𝜑=0.2/2π:
┌───┐ ░ ┌───┐ ░ ┌─┐
q_0: ┤ H ├─■──────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────░──X─┤ H ├─■──────────────■───────────────────■───────────────────────────────────────░─┤M├─────────
├───┤ │ ░ │ └───┘ │P(-π/2) ┌───┐ │ │ ░ └╥┘┌─┐
q_1: ┤ H ├─┼────────■────────■────────────────────────────────────────────────────────────────────────────────────────────────────────────────────░──┼───X───■────────┤ H ├─┼─────────■─────────┼──────────────■────────────────────────░──╫─┤M├──────
├───┤ │ │ │ ░ │ │ └───┘ │P(-π/4) │P(-π/2) │ ┌───┐ │ ░ ║ └╥┘┌─┐
q_2: ┤ H ├─┼────────┼────────┼────────■────────■────────■────────■────────────────────────────────────────────────────────────────────────────────░──┼───X──────────────────■─────────■─────────┼────────┤ H ├─┼─────────■──────────────░──╫──╫─┤M├───
├───┤ │ │ │ │ │ │ │ ░ │ │P(-π/8) └───┘ │P(-π/4) │P(-π/2) ┌───┐ ░ ║ ║ └╥┘┌─┐
q_3: ┤ H ├─┼────────┼────────┼────────┼────────┼────────┼────────┼────────■────────■────────■────────■────────■────────■────────■────────■────────░──X──────────────────────────────────────────■──────────────■─────────■────────┤ H ├─░──╫──╫──╫─┤M├
├───┤ │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) ░ └───┘ ░ ║ ║ ║ └╥┘
q_4: ┤ X ├─■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────░─────────────────────────────────────────────────────────────────────────────────────░──╫──╫──╫──╫─
└───┘ ░ ░ ║ ║ ║ ║
c: 4/══════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════╩══╩══╩══╩═
0 1 2 3 quantum phase estimation circuit with m=5 counting qubits and a phase of 𝜑=0.2/2π:
┌───┐ ░ ┌───┐ ░ ┌─┐
q_0: ┤ H ├─■──────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────░──X─┤ H ├─■──────────────■───────────────────■────────────────────────■────────────────────────────────────────────────────────────░─┤M├────────────
├───┤ │ ░ │ └───┘ │P(-π/2) ┌───┐ │ │ │ ░ └╥┘┌─┐
q_1: ┤ H ├─┼────────■────────■────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────░──┼───X───■────────┤ H ├─┼─────────■─────────┼──────────────■─────────┼────────────────────■───────────────────────────────────────░──╫─┤M├─────────
├───┤ │ │ │ ░ │ │ └───┘ │P(-π/4) │P(-π/2) │ ┌───┐ │ │ │ ░ ║ └╥┘┌─┐
q_2: ┤ H ├─┼────────┼────────┼────────■────────■────────■────────■────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────░──┼───┼──────────────────■─────────■─────────┼────────┤ H ├─┼─────────┼──────────■─────────┼──────────────■────────────────────────░──╫──╫─┤M├──────
├───┤ │ │ │ │ │ │ │ ░ │ │ │P(-π/8) └───┘ │P(-π/4) │ │P(-π/2) │ ┌───┐ │ ░ ║ ║ └╥┘┌─┐
q_3: ┤ H ├─┼────────┼────────┼────────┼────────┼────────┼────────┼────────■────────■────────■────────■────────■────────■────────■────────■────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────░──┼───X──────────────────────────────────────■──────────────■─────────┼──────────■─────────┼────────┤ H ├─┼─────────■──────────────░──╫──╫──╫─┤M├───
├───┤ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ ░ │ │P(-π/16) │P(-π/8) └───┘ │P(-π/4) │P(-π/2) ┌───┐ ░ ║ ║ ║ └╥┘┌─┐
q_4: ┤ H ├─┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────░──X───────────────────────────────────────────────────────────────────■────────────────────■──────────────■─────────■────────┤ H ├─░──╫──╫──╫──╫─┤M├
├───┤ │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) ░ └───┘ ░ ║ ║ ║ ║ └╥┘
q_5: ┤ X ├─■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────░───────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────░──╫──╫──╫──╫──╫─
└───┘ ░ ░ ║ ║ ║ ║ ║
c: 5/════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════╩══╩══╩══╩══╩═
0 1 2 3 4 quantum phase estimation circuit with m=6 counting qubits and a phase of 𝜑=0.2/2π:
┌───┐ ░ ┌───┐ ░ ┌─┐
q_0: ┤ H ├─■──────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────░──X─┤ H ├─■──────────────■───────────────────■────────────────────────■──────────────────────────────■──────────────────────────────────────────────────────────────────────────────────────░─┤M├───────────────
├───┤ │ ░ │ └───┘ │P(-π/2) ┌───┐ │ │ │ │ ░ └╥┘┌─┐
q_1: ┤ H ├─┼────────■────────■────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────░──┼───X───■────────┤ H ├─┼─────────■─────────┼──────────────■─────────┼────────────────────■─────────┼─────────────────────────■────────────────────────────────────────────────────────────░──╫─┤M├────────────
├───┤ │ │ │ ░ │ │ └───┘ │P(-π/4) │P(-π/2) │ ┌───┐ │ │ │ │ │ ░ ║ └╥┘┌─┐
q_2: ┤ H ├─┼────────┼────────┼────────■────────■────────■────────■────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────░──┼───┼──────X───────────■─────────■─────────┼────────┤ H ├─┼─────────┼──────────■─────────┼─────────┼───────────────■─────────┼────────────────────■───────────────────────────────────────░──╫──╫─┤M├─────────
├───┤ │ │ │ │ │ │ │ ░ │ │ │ │P(-π/8) └───┘ │P(-π/4) │ │P(-π/2) │ │ ┌───┐ │ │ │ ░ ║ ║ └╥┘┌─┐
q_3: ┤ H ├─┼────────┼────────┼────────┼────────┼────────┼────────┼────────■────────■────────■────────■────────■────────■────────■────────■────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────░──┼───┼──────X───────────────────────────────■──────────────■─────────┼──────────■─────────┼─────────┼─────────┤ H ├─┼─────────┼──────────■─────────┼──────────────■────────────────────────░──╫──╫──╫─┤M├──────
├───┤ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ ░ │ │ │P(-π/16) │P(-π/8) │ └───┘ │P(-π/4) │ │P(-π/2) │ ┌───┐ │ ░ ║ ║ ║ └╥┘┌─┐
q_4: ┤ H ├─┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────░──┼───X───────────────────────────────────────────────────────────────■────────────────────■─────────┼───────────────■─────────┼──────────■─────────┼────────┤ H ├─┼─────────■──────────────░──╫──╫──╫──╫─┤M├───
├───┤ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ │ ░ │ │P(-π/32) │P(-π/16) │P(-π/8) └───┘ │P(-π/4) │P(-π/2) ┌───┐ ░ ║ ║ ║ ║ └╥┘┌─┐
q_5: ┤ H ├─┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────┼────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────░──X──────────────────────────────────────────────────────────────────────────────────────────────────■─────────────────────────■────────────────────■──────────────■─────────■────────┤ H ├─░──╫──╫──╫──╫──╫─┤M├
├───┤ │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) │P(0.2) ░ └───┘ ░ ║ ║ ║ ║ ║ └╥┘
q_6: ┤ X ├─■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────■────────░────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────░──╫──╫──╫──╫──╫──╫─
└───┘ ░ ░ ║ ║ ║ ║ ║ ║
c: 6/═════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════════╩══╩══╩══╩══╩══╩═
0 1 2 3 4 5 4. Running the algorithm
In this section, we run the results on an ideal simulator and on real hardware both with and without Fire Opal for comparison.
4.1 Running on an ideal simulator
In order to assess the quality of the measured results, we first run the quantum circuits on an ideal, noiseless simulator. This will serve as our benchmark for comparison.
# Run the algorithm on the IBM Aer simulator.
aer_simulator = AerSimulator(method="statevector")
# 3 counting qubits.
transpiled_circuits_3 = qiskit.transpile(
circuits_quantum_phase_estimation_3, aer_simulator
)
simulated_result_3 = (
aer_simulator.run(transpiled_circuits_3, shots=shot_count).result().get_counts()
)
# 4 counting qubits.
transpiled_circuits_4 = qiskit.transpile(
circuits_quantum_phase_estimation_4, aer_simulator
)
simulated_result_4 = (
aer_simulator.run(transpiled_circuits_4, shots=shot_count).result().get_counts()
)
# 5 counting qubits.
transpiled_circuits_5 = qiskit.transpile(
circuits_quantum_phase_estimation_5, aer_simulator
)
simulated_result_5 = (
aer_simulator.run(transpiled_circuits_5, shots=shot_count).result().get_counts()
)
# 6 counting qubits.
transpiled_circuits_6 = qiskit.transpile(
circuits_quantum_phase_estimation_6, aer_simulator
)
simulated_result_6 = (
aer_simulator.run(transpiled_circuits_6, shots=shot_count).result().get_counts()
)
# Turn bitstring counts into probabilities.
simulated_result_3 = [
bitstring_count_to_probabilities(result, shot_count, 3)
for result in simulated_result_3
]
simulated_result_4 = [
bitstring_count_to_probabilities(result, shot_count, 4)
for result in simulated_result_4
]
simulated_result_5 = [
bitstring_count_to_probabilities(result, shot_count, 5)
for result in simulated_result_5
]
simulated_result_6 = [
bitstring_count_to_probabilities(result, shot_count, 6)
for result in simulated_result_6
]4.2 Running the algorithm on quantum hardware without Fire Opal
We run the quantum circuits for quantum phase estimation on a real IBM device both with and without Fire Opal. Note that there might be delays in obtaining the results due to IBM device queues. If you do not want to wait for the circuits to be executed, you can import data we collected previously by skipping this and the next subsection and proceed to the section afterwards.
# Run the algorithm with IBM default.
ibm_backend = service.backend(backend_name, instance=instance)
pm = generate_preset_pass_manager(backend=ibm_backend, optimization_level=1)
sampler = Sampler(ibm_backend)
# 3 counting qubits.
isa_circuits_3 = pm.run(circuits_quantum_phase_estimation_3)
ibm_default_result_3 = sampler.run(isa_circuits_3, shots=shot_count).result()
# 4 counting qubits.
isa_circuits_4 = pm.run(circuits_quantum_phase_estimation_4)
ibm_default_result_4 = sampler.run(isa_circuits_4, shots=shot_count).result()
# 5 counting qubits.
isa_circuits_5 = pm.run(circuits_quantum_phase_estimation_5)
ibm_default_result_5 = sampler.run(isa_circuits_5, shots=shot_count).result()
# 6 counting qubits.
isa_circuits_6 = pm.run(circuits_quantum_phase_estimation_6)
ibm_default_result_6 = sampler.run(isa_circuits_6, shots=shot_count).result()4.3 Running the algorithm on quantum hardware with Fire Opal
Again, it's worth noting that there will be delays in receiving results due to IBM queue times. To see where you are in the queue, you can go to https://quantum-computing.ibm.com/jobs.
# 3 counting qubits.
fire_opal_result_3 = fireopal.execute(
circuits=[circuit.qasm() for circuit in circuits_quantum_phase_estimation_3],
shot_count=shot_count,
credentials=credentials,
backend_name=backend_name,
)
# 4 counting qubits.
fire_opal_result_4 = fireopal.execute(
circuits=[circuit.qasm() for circuit in circuits_quantum_phase_estimation_4],
shot_count=shot_count,
credentials=credentials,
backend_name=backend_name,
)
# 5 counting qubits.
fire_opal_result_5 = fireopal.execute(
circuits=[circuit.qasm() for circuit in circuits_quantum_phase_estimation_5],
shot_count=shot_count,
credentials=credentials,
backend_name=backend_name,
)
# 6 counting qubits.
fire_opal_result_6 = fireopal.execute(
circuits=[circuit.qasm() for circuit in circuits_quantum_phase_estimation_6],
shot_count=shot_count,
credentials=credentials,
backend_name=backend_name,
)4.4 Alternative: Importing the data from a previous run
In case you prefer not to wait in the IBM queue until the above computations are performed, you can also import their results from a run we performed previosuly. Data is available for the following backends:
- "ibm_lagos" for 3-6 counting qubits,
- "ibm_perth" for 3-6 counting qubits,
- "ibmq_guadalupe" for 3-9 counting qubits, and
- "ibmq_jakarta" for 3-6 counting qubits.
def import_data(backend_name, number_of_counting_qubits):
"""
Import the data from a previous run of both IBM default and Fire Opal.
"""
# Import IBM default data.
with open(
f"resources/phase_estimation_{number_of_counting_qubits}_{backend_name}_ibm_default.json",
"r+",
) as file:
results_ibm_default = json.load(file)["results"]
# Import Fire Opal data.
with open(
f"resources/phase_estimation_{number_of_counting_qubits}_{backend_name}_fire_opal.json",
"r+",
) as file:
results_fire_opal = json.load(file)["results"]
# Manipulate the data dictionaries so that keys are sorted bitstrings and values are probabilities.
results_ibm_default = np.array(
[
bitstring_count_to_probabilities(
result["data"]["counts"], shot_count, number_of_counting_qubits
)
for result in results_ibm_default
]
)
results_fire_opal = np.array(
[
bitstring_count_to_probabilities(
result["data"]["counts"], shot_count, number_of_counting_qubits
)
for result in results_fire_opal
]
)
# Place each experiment into a subarray.
results_ibm_default = np.split(
results_ibm_default[0 : experiments_per_phase * len(phases)],
experiments_per_phase,
)
results_fire_opal = np.split(
results_fire_opal[0 : experiments_per_phase * len(phases)],
experiments_per_phase,
)
return results_ibm_default, results_fire_opal
# Import the data for 3 counting qubits.
ibm_default_result_3, fire_opal_result_3 = import_data(backend_name, 3)
data_3 = {
"simulation": simulated_result_3,
"ibm_default": ibm_default_result_3,
"fire_opal": fire_opal_result_3,
}
# Import the data for 4 counting qubits.
ibm_default_result_4, fire_opal_result_4 = import_data(backend_name, 4)
data_4 = {
"simulation": simulated_result_4,
"ibm_default": ibm_default_result_4,
"fire_opal": fire_opal_result_4,
}
# Import the data for 5 counting qubits.
ibm_default_result_5, fire_opal_result_5 = import_data(backend_name, 5)
data_5 = {
"simulation": simulated_result_5,
"ibm_default": ibm_default_result_5,
"fire_opal": fire_opal_result_5,
}
# Import the data for 6 counting qubits.
ibm_default_result_6, fire_opal_result_6 = import_data(backend_name, 6)
data_6 = {
"simulation": simulated_result_6,
"ibm_default": ibm_default_result_6,
"fire_opal": fire_opal_result_6,
}5. Analyzing the data
Let's assess the data we obtained. For this, we first plot the ideal and measured distributions.
# Define the functions to plot the distributions.
def highest_count_bitstring(data):
"""
Return the bitstring with the highest count for the measurement results of a single shot.
"""
return max(data, key=data.get)
def counting_qubits(data):
"""
Check whether all passed dictionaries with single shot measurements belong to the same
number of counting qubits and return this number.
"""
first_keys = iter(list(result.keys())[0] for result in data.values())
number_of_counting_qubits = len(next(first_keys))
if not all(len(key) == number_of_counting_qubits for key in first_keys):
raise ValueError(
"The input data must contain dictionaries corresponding to the same number of counting qubits."
)
return number_of_counting_qubits
def plot_distributions(data, phase_index):
"""
Plot the distributions of all single shot measurements for all three experiments
"simulation", "ibm_default" and "fire_opal" in one plot.
"""
number_of_counting_qubits = counting_qubits(data)
phase = f"1/{int(1/phases[phase_index])}"
simulation = data["simulation"]
ibm_default = data["ibm_default"]
fire_opal = data["fire_opal"]
correct_result = highest_count_bitstring(simulation)
max_probability = max(
max(simulation.values()), max(ibm_default.values()), max(fire_opal.values())
)
figure, axes = plt.subplots(1, 3, layout="constrained")
# Ideal simulation.
axes[0].bar(
list(simulation.keys()),
list(simulation.values()),
color="#32A4A8",
label="Ideal simulation",
)
axes[0].legend()
# IBM default.
bars_ibm_default = axes[1].bar(
list(ibm_default.keys()),
list(ibm_default.values()),
color="#6C5C71",
label="IBM default",
)
axes[1].legend()
bars_ibm_default[list(ibm_default.keys()).index(correct_result)].set_color(
"#D84144"
)
# Fire Opal.
bars_fire_opal = axes[2].bar(
list(fire_opal.keys()),
list(fire_opal.values()),
color="#680CE9",
label="Fire Opal",
)
axes[2].legend()
bars_fire_opal[list(fire_opal.keys()).index(correct_result)].set_color("#D84144")
for axis in axes.flat:
axis.set(xlabel="Bitstrings", ylabel="Probabilities")
axis.set_ylim([0, max_probability])
axis.tick_params(axis="x", labelrotation=90)
# Hide labels and tick labels for inner plots.
for axis in axes.flat:
axis.label_outer()
figure.suptitle(
rf"{number_of_counting_qubits} counting qubits, $2 \pi \varphi$={phase}"
)
figure.set_size_inches(20, 5)
plt.show()
# Plot the distributions.
# Specify, which of the data sets to plot by specifying the index of the
# experiment (that is, run 0 or run 1) and the index of the phase in the list
# `phases`. Remember that we performed the algoritm twice for each phase,
# so the `experiment_index` specifies, which run is plotted.
experiment_index = 0
phase_index = 0
# Plot the distributions for 3 counting qubits.
data_3_first_phase = {
"simulation": data_3["simulation"][phase_index],
"ibm_default": data_3["ibm_default"][experiment_index][phase_index],
"fire_opal": data_3["fire_opal"][experiment_index][phase_index],
}
plot_distributions(data_3_first_phase, phase_index)
# Plot the distributions for 4 counting qubits.
distributions_4 = {
"simulation": data_4["simulation"][phase_index],
"ibm_default": data_4["ibm_default"][experiment_index][phase_index],
"fire_opal": data_4["fire_opal"][experiment_index][phase_index],
}
plot_distributions(distributions_4, phase_index)
# Plot the distributions for 5 counting qubits.
distributions_5 = {
"simulation": data_5["simulation"][phase_index],
"ibm_default": data_5["ibm_default"][experiment_index][phase_index],
"fire_opal": data_5["fire_opal"][experiment_index][phase_index],
}
plot_distributions(distributions_5, phase_index)
# Plot the distributions for 6 counting qubits.
distributions_6 = {
"simulation": data_6["simulation"][phase_index],
"ibm_default": data_6["ibm_default"][experiment_index][phase_index],
"fire_opal": data_6["fire_opal"][experiment_index][phase_index],
}
plot_distributions(distributions_6, phase_index)
As you can see, the results obtained without Fire Opal are very different to the ideal, correct distribution. Due to high noise, the qubits get corrupted while the quantum phase estimation circuit is running. This leads to measurement distributions, which look pretty much random and meaningless. On the contrary, Fire Opal produces the correct results and the measurement distributions look very similar to the ideal ones. This observation is independent of the number of counting qubits. Even for $m = 6$ counting qubits, Fire Opal recovers the same phase as the ideal simulation. Due to a decent number of counting qubits, we obtain a good estimate for the phase in this case. The actual phase here is $\varphi = \frac{1}{6 \pi} \approx 0.0531$. Fire Opal predicts a phase between $3 / 2^6 \approx 0.0469$ and $4 / 2^6 = 0.0625$, but closer to $0.0469$, which is a reasonable approximation. Furthermore, due to its scalability Fire Opal allows for obtaining even more accurate estimates by increasing the number of counting qubits.
Let's examine the success rates, that is how often we predict the correct highest count bitstring.
def success_rate(data):
"""
Compare the highest count bitstring of the ideal simulation with the ones obtained
from the measurements. Return the percentage of correct predictions.
"""
simulation = data["simulation"]
ibm_default = data["ibm_default"]
fire_opal = data["fire_opal"]
simulation_results = [highest_count_bitstring(ideal) for ideal in simulation]
ibm_default_results = [
highest_count_bitstring(result) for results in ibm_default for result in results
]
fire_opal_results = [
highest_count_bitstring(result) for results in fire_opal for result in results
]
correct_results_ibm_default = 0
correct_results_fire_opal = 0
for phase_index, correct_result in enumerate(simulation_results):
if ibm_default_results[phase_index] == correct_result:
correct_results_ibm_default += 1
if fire_opal_results[phase_index] == correct_result:
correct_results_fire_opal += 1
return (
100 * correct_results_ibm_default / len(ibm_default_results),
100 * correct_results_fire_opal / len(fire_opal_results),
)
def plot_success_rates(data):
"""
Plot the success rates of a Fire Opal experiment and an IBM default experiment for a different number
of counting qubits.
"""
ibm_default, fire_opal = zip(*data.values())
bar_width = 0.3
x_ticks = np.array(list(data.keys()))
plt.xticks(ticks=x_ticks)
# Fire Opal.
plt.bar(
x_ticks - bar_width / 2,
fire_opal,
color="#680CE9",
width=bar_width,
label="Fire Opal",
)
# IBM default.
plt.bar(
x_ticks + bar_width / 2,
ibm_default,
color="#6C5C71",
width=bar_width,
label="IBM default",
)
plt.title("Highest Count Success Rates")
plt.xlabel(r"Number of counting qubits $m$")
plt.ylabel("Success rate in %")
plt.legend(loc="upper left", bbox_to_anchor=(1, 1))
plt.show()
success_rates = {
3: success_rate(data_3),
4: success_rate(data_4),
5: success_rate(data_5),
6: success_rate(data_6),
}
plot_success_rates(data=success_rates)
As you can see, Fire Opal consistently predicts the correct answer even as the number of counting qubits increases. Without Fire Opal, the algorithm was unable to scale past even 3 counting qubits.
6. Conclusion
We performed quantum phase estimation for 10 different phases and repeated each experiment twice. In this way, we obtained 20 data points to average over in order to compute the success rates. From the plot, we can see that the IBM default configuration generally does not recover the correct bitstring as its highest count for four or more counting qubits. In the case of five or six counting qubits, the success rate is even zero, that is all predictions are wrong. On the other hand, Fire Opal achieves a success rate of 100 % no matter how many counting qubits are used.
Fire Opal allows you to get correct answers from quantum phase estimation without additional overhead, whereas running on the default hardware typically results in pure noise. By allowing you to scale to higher numbers of counting qubits without any performance degradation, Fire Opal provides better accuracy for the estimated phase. This example shows how Fire Opal is an excellent tool for near term quantum chemistry applications, just as it can be applied to other algorithms to similarly improve the likelihood of achieving quality results.
The package versions below were used to produce this notebook.
from fireopal import print_package_versions
print_package_versions()| Package | Version |
| --------------------- | ------- |
| Python | 3.11.0 |
| networkx | 2.8.8 |
| numpy | 1.26.2 |
| sympy | 1.12 |
| fire-opal | 6.7.0 |
| qctrl-workflow-client | 2.2.0 |
