diff --git a/docs/tutorials/_toc.json b/docs/tutorials/_toc.json
index 9fbb6f4eae1..4e1bc219575 100644
--- a/docs/tutorials/_toc.json
+++ b/docs/tutorials/_toc.json
@@ -178,6 +178,10 @@
{
"title": "Circuit cutting for depth reduction",
"url": "/docs/tutorials/depth-reduction-with-circuit-cutting"
+ },
+ {
+ "title": "Readout error mitigation for the Sampler primitive using M3",
+ "url": "/docs/tutorials/readout-error-mitigation-sampler"
}
]
},
diff --git a/docs/tutorials/index.mdx b/docs/tutorials/index.mdx
index ba9110988ca..871a57a65d8 100644
--- a/docs/tutorials/index.mdx
+++ b/docs/tutorials/index.mdx
@@ -135,6 +135,8 @@ Addons enable advanced circuit manipulation, such as cutting, backpropagating ob
* [Circuit cutting for depth reduction](/docs/tutorials/depth-reduction-with-circuit-cutting)
+* [Readout error mitigation for the Sampler primitive using M3](/docs/tutorials/readout-error-mitigation-sampler)
+
diff --git a/docs/tutorials/readout-error-mitigation-sampler.ipynb b/docs/tutorials/readout-error-mitigation-sampler.ipynb
new file mode 100644
index 00000000000..a5bbfa0fcce
--- /dev/null
+++ b/docs/tutorials/readout-error-mitigation-sampler.ipynb
@@ -0,0 +1,1417 @@
+{
+ "cells": [
+ {
+ "cell_type": "markdown",
+ "id": "70810e0b-08f9-48ae-beb5-8c863819000b",
+ "metadata": {},
+ "source": [
+ "{/* cspell:ignore braket, Zgrm, newcommand, probs, quasis, topten */}\n",
+ "# Readout error mitigation for the Sampler primitive using M3\n",
+ "\n",
+ "*Estimated QPU usage: 1 minute (tested on IBM `ibm_kingston`)*\n",
+ "\n",
+ "## Background\n",
+ "\n",
+ "Unlike the Estimator primitive, the Sampler primitive does not have built-in support for error mitigation.\n",
+ "Several of the methods supported by the Estimator are specifically designed for expectation values, and hence are not applicable to the Sampler primitive.\n",
+ "However, readout error mitigation is a highly effective method that is also applicable to the Sampler primitive.\n",
+ "This tutorial explains how to use the [M3 Qiskit addon](https://qiskit.github.io/qiskit-addon-mthree/) to mitigate readout error for the Sampler primitive.\n",
+ "\n",
+ "The M3 Qiskit addon implements an efficient method for readout error mitigation.\n",
+ "\n",
+ "### What is readout error?\n",
+ "\n",
+ "Immediately before measurement, the state of a qubit register is\n",
+ "described by a superposition of computational basis states,\n",
+ "or by a density matrix.\n",
+ "Measurement of the qubit register into a classical bit register then proceeds in two steps.\n",
+ "First the quantum measurement proper is performed.\n",
+ "This means that the state of the qubit register\n",
+ "is projected onto a single basis state that is characterized\n",
+ "by a string of $1$s and $0$s.\n",
+ "The second step consists of reading the bitstring characterizing this basis state\n",
+ "and writing it into classical computer memory.\n",
+ "We call this step *readout*.\n",
+ "It turns out that readout error is more important than error due to to incorrectly collapsing the state.\n",
+ "This makes sense when you recall that readout requires detecting a microscopic\n",
+ "quantum state and amplifying it to the macroscopic realm. A readout resonator is coupled to\n",
+ "the (transmon) qubit, thereby experiencing a very small frequency shift. A microwave pulse\n",
+ "is then bounced off of the resonator, in turn experiencing small changes in its\n",
+ "characteristics. The reflected pulse is then amplified and analyzed. This is a delicate\n",
+ "process subject to a host of errors.\n",
+ "\n",
+ "The important point is that, while both quantum measurement and readout are subject to error, the\n",
+ "latter incurs the dominant error, called readout error.\n",
+ "In the following we be concerned only with readout error."
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "ba74fc81-30ac-4365-9ba6-a891ba082bd6",
+ "metadata": {
+ "jp-MarkdownHeadingCollapsed": true
+ },
+ "source": [
+ "### Theoretical background\n",
+ "\n",
+ "If the sampled bitstring (stored in classical memory) differs from the bitstring characterizing\n",
+ "the projected quantum state, we say that a readout error has occurred.\n",
+ "These errors are observed to be random and uncorrelated from sample to sample.\n",
+ "It has proven useful to model readout error as a _noisy classical channel_.\n",
+ "That is, for every pair of\n",
+ "bitstrings $i$ and $j$, there is fixed probability that a true value of $j$ will\n",
+ "be incorrectly read as $i$.\n",
+ "\n",
+ "More precisely, for every pair of bitstrings $(i, j)$, there is a (conditional) probability ${M}_{i,j}$\n",
+ "that $i$ is read, given that the true value is $j$.\n",
+ "That is,\n",
+ "$$\n",
+ " {M}_{i,j} = \\Pr(\\text{readout value is } i | \\text{true value is } j)\n",
+ " \\text{ for } i,j \\in (0,...,2^n - 1), \\tag{1}\n",
+ "$$\n",
+ "where $n$ is the number of bits in the readout register.\n",
+ "For concreteness, we assume that $i$ is a decimal integer whose binary representation is\n",
+ "the bitstring that labels the computational basis states.\n",
+ "We call the $2^n \\times 2^n$ matrix ${M}$ the _assignment matrix_.\n",
+ "For fixed true value $j$, summing the probability over all noisy outcomes $i$ must give one. That is\n",
+ "$$\n",
+ " \\sum_{i=0}^{2^n - 1} {M}_{i,j} = 1 \\text{ for all } j\n",
+ "$$\n",
+ "A matrix with no negative entries that satisfies (1) is called\n",
+ "_left-stochastic_.\n",
+ "A left-stochastic matrix is also called _column-stochastic_ because each of its columns sums to $1$.\n",
+ "We experimentally determine approximate values for each element ${M}_{i,j}$ by\n",
+ "repeatedly preparing each basis state $|j \\rangle$ and then computing the frequencies\n",
+ "of occurrence of sampled bitstrings.\n",
+ "\n",
+ "If an experiment involves estimating a probability distribution over output bitstrings by repeated sampling,\n",
+ "then we can use ${M}$ to mitigate readout error at the level of the distribution.\n",
+ "The first step is to repeat a fixed circuit of interest many times,\n",
+ "creating a histogram of sampled bitstrings.\n",
+ "The normalized histogram is the measured probability distribution over\n",
+ "the $2^n$ possible bitstrings, which we denote by ${\\tilde{p}} \\in \\mathbb{R}^{2^n}$.\n",
+ "The (estimated) probability ${{\\tilde{p}}}_i$ of sampling bitstring $i$\n",
+ "is equal to the sum over all true bitstrings $j$, each weighted by\n",
+ "the probability that it is mistaken for $i$.\n",
+ "This statement in matrix form is\n",
+ "$$\n",
+ " {\\tilde{p}} = {M} {\\vec{p}}, \\tag{2}\n",
+ "$$\n",
+ "where ${\\vec{p}}$ is the true distribution. In words, the readout error has the effect of multiplying\n",
+ "the ideal distribution over bitstrings ${\\vec{p}}$ by the assignment matrix ${M}$ to\n",
+ "produce the observed distribution ${\\tilde{p}}$.\n",
+ "We have measured ${\\tilde{p}}$ and ${M}$, but have no direct access to ${\\vec{p}}$. In principle, we will\n",
+ "obtain the true distribution of bitstrings for our circuit\n",
+ "by solving equation (2) for ${\\vec{p}}$ numerically.\n",
+ "\n",
+ "Before we move on, it's worth noting some important features of this naive approach.\n",
+ "\n",
+ "- In practice, equation (2) is not solved by inverting ${M}$. Linear algebra\n",
+ " routines in software libraries employ methods that are more stable, accurate, and efficient.\n",
+ "- When estimating ${M}$, we assumed that only readout errors occurred. In particular,\n",
+ " we assume there were no state preparation and quantum measurement errorsβ\n",
+ " or at least that they were otherwise mitigated.\n",
+ " To the extent that this is a good assumption, ${M}$ really represents\n",
+ " only readout error. But when we _use_ ${M}$ to correct a measured distribution\n",
+ " over bitstrings, we make no such assumption. In fact, we expect an interesting\n",
+ " circuit to introduce noise, for instance, gate errors. The \"true\" distribution\n",
+ " still includes effects from any errors that are not otherwise mitigated.\n",
+ "\n",
+ "This method, while useful in some circumstances, suffers from a few limitations.\n",
+ "\n",
+ "The space and time resources needed to estimate ${M}$ grow exponentially in $n$:\n",
+ "- The estimation of ${M}$ and ${\\tilde{p}}$ is subject to statistical error due to finite sampling.\n",
+ " This noise can be made as small as desired\n",
+ " at the cost of more shots (up to the timescale of drifting hardware parameters\n",
+ " that result in systematic errors in ${M}$.)\n",
+ " However, if no assumptions are made on the bitstrings observed\n",
+ " when performing mitigation, the number of shots required to estimate ${M}$ grows\n",
+ " at least exponentially in $n$.\n",
+ "- ${M}$ is a $2^n \\times 2^n$ matrix.\n",
+ " When $n>10$, the amount of memory required to store ${M}$ is\n",
+ " greater than the memory available in a powerful laptop.\n",
+ "\n",
+ "Further limitations are:\n",
+ "\n",
+ "- The recovered distribution ${\\vec{p}}$ may have one\n",
+ " or more negative probabilities (while still summing to one). One solution\n",
+ " is to minimize $||{M} {\\vec{p}} - {\\tilde{p}}||^2$ subject to the constraint that\n",
+ " each entry in ${\\vec{p}}$ be non-negative. However, the runtime of such a\n",
+ " method is orders of magnitude longer than directly solving equation (2).\n",
+ "- This mitigation procedure works on the level of a probability distribution\n",
+ " over bitstrings. In particular, it cannot correct an error in an individual\n",
+ " observed bitstring."
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "c98d8409-66ad-452f-b87a-d8df886d77d4",
+ "metadata": {},
+ "source": [
+ "### Qiskit M3 addon: Scaling to longer bitstrings\n",
+ "\n",
+ "The M3 method is an alternative to direct solution of equation (2)\n",
+ "for mitigating readout error in the experimental distribution over bitstrings.\n",
+ "While the direct method is limited to bitstrings no longer than about $10$ bits,\n",
+ "M3 can handle much longer bitstrings.\n",
+ "The two important features of M3 are\n",
+ "- Correlations in readout error of order three and higher among collections of bits\n",
+ " are assumed to be negligible and are ignored. In principle, at the cost of more shots,\n",
+ " one could estimate higher correlations as well.\n",
+ "- Rather than constructing ${M}$ explicitly, we use a much smaller effective matrix that records\n",
+ " probabilities only for bitstrings collected when constructing ${\\tilde{p}}$.\n",
+ "\n",
+ "At a high level, the procedure works as follows.\n",
+ "\n",
+ "First, we construct building blocks from which we can construct a simplified, effective, description of ${M}$.\n",
+ "Then, we repeatedly run the circuit of interest and collect bitstrings which we use to construct\n",
+ "both ${\\tilde{p}}$ and, with the help of the building blocks, an effective ${M}$.\n",
+ "\n",
+ "More precisely,\n",
+ "- Single-qubit assignment matrices are estimated for each qubit. To do this, we repeatedly\n",
+ " prepare the qubit register in the all-zero state $|0 ... 0 \\rangle$ and then in the all-one\n",
+ " state $|1 ... 1 \\rangle$, and record the probability for each qubit that it is read\n",
+ " incorrectly.\n",
+ "- Correlations of order three and higher are assumed to be negligible and are ignored.\n",
+ "\n",
+ " Instead we construct a number $n$ of $2 \\times 2$ single qubit\n",
+ " assignment matrices, and a number $n(n-1)/2$ of $4 \\times 4$ two-qubit assignment\n",
+ " matrices. These one- and two-qubit assignment matrices are stored for later\n",
+ " use.\n",
+ "- After repeatedly sampling a circuit to construct ${\\tilde{p}}$,\n",
+ " we construct an effective approximation to ${M}$ using only\n",
+ " bitstrings that are sampled when constructing ${\\tilde{p}}$. This effective matrix\n",
+ " is built using the one- and two-qubit matrices described in the previous item.\n",
+ " The linear dimension of this matrix is at most of the order of the number\n",
+ " of shots used in constructing ${\\tilde{p}}$, which is much smaller than\n",
+ " the dimension $2^n$ of the full assignment matrix ${M}$ .\n",
+ "\n",
+ "For technical details on M3, you can refer to [*Scalable Mitigation of Measurement Errors on Quantum Computers*](https://journals.aps.org/prxquantum/abstract/10.1103/PRXQuantum.2.040326)."
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "14f04ac4-5142-4436-93e1-f712d47dde1a",
+ "metadata": {},
+ "source": [
+ "### Application of M3 to a quantum algorithm\n",
+ "We'll apply M3's readout mitigation to the hidden shift problem. The hidden shift problem, and closely related problems such as the [hidden subgroup problem](https://en.wikipedia.org/wiki/Hidden_subgroup_problem), were originally conceived in a fault-tolerant setting (more precisely, before fault-tolerant QPUs were proved to be possible!). But they are studied with available processors as well. An example of algorithmic exponential speedup obtained for a variant of the hidden shift problem obtained on 127-qubit IBM QPU's can be found in [this paper](https://journals.aps.org/prx/accepted/a9074K06A8e1590147da9c69f8c4b64c28247be5a) ([arXiv version](https://arxiv.org/abs/2401.07934)).\n",
+ "\n",
+ "In the following, all arithmetic is Boolean.\n",
+ "That is, for $a, b \\in \\mathbb{Z}_2 = \\{0, 1\\}$, addition, $a + b$ is the logical XOR function.\n",
+ "Furthermore, multiplication $a \\times b$ (or $a b$) is the logical AND function. For $x, y \\in \\{0, 1\\}^n$,\n",
+ "$x + y$ is defined by bitwise application of XOR.\n",
+ "The dot product $\\cdot: {\\mathbb{Z}_2^n} \\rightarrow \\mathbb{Z}_2$ is defined\n",
+ "by $x \\cdot y = \\sum_i x_i y_i$."
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "81cb34b4-566b-4ee1-8309-7464a5fd3bb7",
+ "metadata": {},
+ "source": [
+ "#### Hadamard operator and Fourier transform\n",
+ "\n",
+ "In implementing quantum algorithms, it is very common to use the Hadamard operator as a Fourier transform.\n",
+ "The computational basis states are sometimes called _classical states_. They stand in\n",
+ "a one-to-one relation to the classical bitstrings.\n",
+ "The $n$-qubit Hadamard operator on classical states can be viewed as a Fourier transform on the Boolean hypercube:\n",
+ "$$\n",
+ "H^{\\otimes n} = \\frac{1}{\\sqrt{2^n}} \\sum_{x,y \\in {\\mathbb{Z}_2^n}} (-1)^{x \\cdot y} {|{y}\\rangle}{\\langle{x}|}.\n",
+ "$$\n",
+ "Consider a state ${|{s}\\rangle}$ corresponding to fixed bitstring $s$.\n",
+ "Applying $H^{\\otimes n}$, and using ${\\langle {x}|{s}\\rangle} = \\delta_{x,s}$,\n",
+ "we see that the Fourier transform of ${|{s}\\rangle}$ can be written\n",
+ "$$\n",
+ " H^{\\otimes n} {|{s}\\rangle} = \\frac{1}{\\sqrt{2^n}} \\sum_{y \\in {\\mathbb{Z}_2^n}} (-1)^{s \\cdot y} {|{y}\\rangle}.\n",
+ "$$\n",
+ "\n",
+ "The Hadamard is its own inverse, that is,\n",
+ " $H^{\\otimes n} H^{\\otimes n} = (H H)^{\\otimes n} = I^{\\otimes n}$.\n",
+ "Thus, the inverse Fourier transform is the same operator, $H^{\\otimes n}$.\n",
+ "Explicitly, we have,\n",
+ "$$\n",
+ " {|{s}\\rangle} = H^{\\otimes n} H^{\\otimes n} {|{s}\\rangle} = H^{\\otimes n} \\frac{1}{\\sqrt{2^n}} \\sum_{y \\in {\\mathbb{Z}_2^n}} (-1)^{s \\cdot y} {|{y}\\rangle}.\n",
+ "$$"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "d10fce93-3e08-4a27-b641-eda5e7b05ce0",
+ "metadata": {},
+ "source": [
+ "#### The hidden shift problem\n",
+ "\n",
+ "We consider a simple example of a _hidden shift problem_.\n",
+ "The problem is to identify a constant shift in the input to a function.\n",
+ "The function we consider is the dot product. It is the simplest member\n",
+ "of a large class of functions that admit a quantum speedup for the hidden shift\n",
+ "problem via techniques similar to those presented below.\n",
+ "\n",
+ "Let $x,y \\in {\\mathbb{Z}_2^m}$ be bitstrings of length $m$.\n",
+ "We define ${f}: {\\mathbb{Z}_2^m} \\times {\\mathbb{Z}_2^m} \\rightarrow \\{-1,1\\}$ by\n",
+ "$$\n",
+ " {f}(x, y) = (-1)^{x \\cdot y}.\n",
+ "$$\n",
+ " Let $a,b \\in {\\mathbb{Z}_2^m}$ be fixed bitstrings of length $m$.\n",
+ " We furthermore define $g: {\\mathbb{Z}_2^m} \\times {\\mathbb{Z}_2^m} \\rightarrow \\{-1,1\\}$ by\n",
+ "$$\n",
+ " g(x, y) = {f}(x+a, y+b) = (-1)^{(x+a) \\cdot (y+b)},\n",
+ " $$\n",
+ " where $a$ and $b$ are (hidden) parameters.\n",
+ " We are given two black boxes, one implementing $f$, and the other $g$.\n",
+ " We suppose that we know that they compute the functions defined above, except that we know\n",
+ " neither $a$ nor $b$. The game is to determine the hidden bitstrings (shifts)\n",
+ " $a$ and $b$ by making queries to $f$ and $g$. It's clear that if we play the game classically,\n",
+ " we need $O(2m)$ queries to determine $a$ and $b$. For example, we can query $g$\n",
+ " with all pairs of strings that combined have exactly one element set to $1$.\n",
+ " On each query, we learn one element of either $a$ or $b$.\n",
+ " However, we will see that, if the black boxes are implemented as quantum circuits, we can\n",
+ " determine $a$ and $b$ with a single query to each of $f$ and $g$.\n",
+ "\n",
+ " In the context of algorithmic complexity, a black box is called an _oracle_.\n",
+ " In addition to being opaque, an oracle has the property that it consumes the input and\n",
+ " produces the output instantly, adding nothing to the complexity budget of the algorithm\n",
+ " in which it is embedded. In fact, in the case at hand, the oracles implementing $f$ and\n",
+ " $g$ will be seen to be efficient."
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "d5e4c1e2-d6f5-4093-afdb-f97e9f062179",
+ "metadata": {
+ "jp-MarkdownHeadingCollapsed": true
+ },
+ "source": [
+ "#### Quantum circuits for $f$ and $g$\n",
+ "\n",
+ "We will need the following ingredients in order to implement $f$ and $g$ as quantum circuits.\n",
+ "\n",
+ "For single qubit classical states ${|{x_1}\\rangle}, {|{y_1}\\rangle}$, with $x_1,y_1 \\in \\mathbb{Z}_2$,\n",
+ "the controlled $Z$ gate ${CZ}$ may be written\n",
+ "$$\n",
+ "{CZ} {|{x_1}\\rangle}{|{y_1}\\rangle}{x_1} = (-1)^{x_1 y_1} {|{x_1}\\rangle}{x_1}{|{y_1}\\rangle}.\n",
+ "$$\n",
+ "We will operate with $m$ CZ gates, one on $(x_1, y_1)$, and one on $(x_2, y_2)$, and so on, through $(x_m, y_m)$.\n",
+ "We call this operator ${CZ}_{x,y}$.\n",
+ "\n",
+ "$U_f = {CZ}_{x,y}$ is a quantum version of ${f} = {f}(x,y)$:\n",
+ "$$\n",
+ "%\\CZ_{x,y} {|#1\\rangle}{z} =\n",
+ "U_f {|{x}\\rangle}{|{y}\\rangle} = {CZ}_{x,y} {|{x}\\rangle}{|{y}\\rangle} = (-1)^{x \\cdot y} {|{x}\\rangle}{|{y}\\rangle}.\n",
+ "$$\n",
+ "\n",
+ "We also need to implement a bitstring shift.\n",
+ "We denote the operator on the $x$ register $X^{a_1}\\cdots X^{a_m}$ by $X_a$\n",
+ "and likewise on the $y$ register $X_b = X^{b_1}\\cdots X^{b_m}$.\n",
+ "These operators apply $X$ wherever a single bit is $1$, and the identity $I$ wherever it is $0$.\n",
+ "Then we have\n",
+ "$$\n",
+ " X_a X_b {|{x}\\rangle}{|{y}\\rangle} = {|{x+a}\\rangle}{|{y+b}\\rangle}.\n",
+ "$$\n",
+ "\n",
+ "The second black box $g$ is implemented by the unitary $U_g$, given by\n",
+ "$$\n",
+ "%U_g {|{x}\\rangle}{|{y}\\rangle} = X_aX_b \\CZ_{x,y} X_aX_b {|{x}\\rangle}{|{y}\\rangle}.\n",
+ "U_g = X_aX_b {CZ}_{x,y} X_aX_b.\n",
+ "$$\n",
+ "To see this, we apply the operators from right to left to the state ${|{x}\\rangle}{|{y}\\rangle}$.\n",
+ "First\n",
+ "\n",
+ "$$\n",
+ " X_a X_b {|{x}\\rangle}{|{y}\\rangle} = {|{x+a}\\rangle}{|{y+b}\\rangle}.\n",
+ "$$\n",
+ "\n",
+ "Then,\n",
+ "$$\n",
+ " {CZ}_{x,y} {|{x+a}\\rangle}{|{y+b}\\rangle} = (-1)^{(x+a)\\cdot (y+b)} {|{x+a}\\rangle}{|{y+b}\\rangle}\n",
+ "$$\n",
+ "\n",
+ "Finally,\n",
+ "\n",
+ "$$\n",
+ " X^a X^b (-1)^{(x+a)\\cdot (y+b)} {|{x+a}\\rangle}{|{y+b}\\rangle} = (-1)^{(x+a)\\cdot (y+b)} {|{x}\\rangle}{|{y}\\rangle},\n",
+ "$$\n",
+ "\n",
+ "which is indeed the quantum version of $f(x+a, y+b)$."
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "de0bb5f9-517e-4bea-91ab-f65e012e8ae9",
+ "metadata": {
+ "editable": true,
+ "jp-MarkdownHeadingCollapsed": true,
+ "slideshow": {
+ "slide_type": ""
+ },
+ "tags": []
+ },
+ "source": [
+ "#### The hidden shift algorithm\n",
+ "\n",
+ "Now we put the pieces constructed above together to solve the hidden shift problem.\n",
+ "We begin by applying Hadamards to the registers initialized to the all-zero state.\n",
+ "$$\n",
+ "H^{\\otimes 2m} = H^{\\otimes m} \\otimes H^{\\otimes m} {{|{0}\\rangle}^{\\otimes m}}{{|{0}\\rangle}^{\\otimes m}} = \\frac{1}{\\sqrt{2^{2m}}} \\sum_{x, y \\in {\\mathbb{Z}_2^m}} (-1)^{x \\cdot y} {|{x}\\rangle}{|{y}\\rangle}.\n",
+ "$$\n",
+ "\n",
+ "Next, we query the oracle $g$ to arrive at\n",
+ "$$\n",
+ "U_g H^{\\otimes 2m} {{|{0}\\rangle}^{\\otimes m}}{{|{0}\\rangle}^{\\otimes m}}\n",
+ "= \\frac{1}{\\sqrt{2^{2m}}} \\sum_{x, y \\in {\\mathbb{Z}_2^m}} (-1)^{(x+a) \\cdot (y+b)} {|{x}\\rangle}{|{y}\\rangle}\n",
+ "$$\n",
+ "$$\n",
+ "\\approx \\frac{1}{\\sqrt{2^{2m}}} \\sum_{x, y \\in {\\mathbb{Z}_2^m}} (-1)^{x \\cdot y + x \\cdot b + y \\cdot a} {|{x}\\rangle}{|{y}\\rangle}.\n",
+ "$$\n",
+ "In the last line, we omitted the constant global phase factor $(-1)^{a \\cdot b}$,\n",
+ "and denote equality up to a phase by $\\approx$.\n",
+ "Next, applying the oracle $f$ introduces another factor of $(-1)^{x \\cdot y}$, canceling the one already\n",
+ "present. We then have,\n",
+ "$$\n",
+ "U_f U_g H^{\\otimes 2m} {{|{0}\\rangle}^{\\otimes m}}{{|{0}\\rangle}^{\\otimes m}}\n",
+ "\\approx \\frac{1}{\\sqrt{2^{2m}}} \\sum_{x, y \\in {\\mathbb{Z}_2^m}} (-1)^{x \\cdot b + y \\cdot a} {|{x}\\rangle}{|{y}\\rangle}.\n",
+ "$$\n",
+ "The final step is to apply the inverse Fourier transform, $H^{\\otimes 2m} = H^{\\otimes m} \\otimes H^{\\otimes m}$,\n",
+ "resulting in\n",
+ "$$\n",
+ "H^{\\otimes 2m} U_f U_g H^{\\otimes 2m} {{|{0}\\rangle}^{\\otimes m}}{{|{0}\\rangle}^{\\otimes m}}\n",
+ "\\approx {|{b}\\rangle}{|{a}\\rangle}.\n",
+ "$$\n",
+ "The circuit is finished. In the absence of noise, sampling the quantum registers will\n",
+ "return the bitstrings $b, a$ with probability one.\n",
+ "\n",
+ "The Boolean inner product is an example of the so-called bent functions.\n",
+ "We will not define bent functions here\n",
+ "but merely note that they\n",
+ "\"are maximally resistant against attacks that seek to exploit a dependence of\n",
+ "the outputs on some linear subspace of the inputs.\"\n",
+ "This quote is from the article [_Quantum algorithms for highly non-linear Boolean functions_](https://arxiv.org/abs/0811.3208), which\n",
+ "gives efficient hidden-shift algorithms for several classes of bent functions.\n",
+ "The algorithm in this tutorial appears in Section 3.1 of the article.\n",
+ "\n",
+ "In the more general case, the circuit for finding a hidden shift $s \\in \\mathbb{Z}^n$ is\n",
+ "$$\n",
+ " H^{\\otimes n} U_{\\tilde{f}} H^{\\otimes n} U_g H^{\\otimes n} {|{0}\\rangle}^{\\otimes n} = {|{s}\\rangle}.\n",
+ "$$\n",
+ " In the general case, $f$ and $g$ are functions of a single variable.\n",
+ " Our example of the inner product has this form if we let $f(x, y) \\to f(z)$,\n",
+ " with $z$ equal to the concatenation of $x$ and $y$, and $s$ equal to the concatenation\n",
+ " of $a$ and $b$.\n",
+ " The general case requires exactly two oracles: One oracle for $g$ and one for $\\tilde{f}$,\n",
+ " where the latter is a function known as the _dual_ of the bent function $f$.\n",
+ " The inner product function has the self-dual property $\\tilde{f}=f$.\n",
+ "\n",
+ " In our circuit for the hidden shift on the inner product we omitted the middle layer\n",
+ " of Hadamards that appears in the circuit for the general case. While in the general case\n",
+ " this layer is necessary, we saved a bit of depth by omitting it, at the expense of a bit\n",
+ " of post-processing because the output is ${|{b}\\rangle}{|{a}\\rangle}$ instead of the desired ${|{a}\\rangle}{|{b}\\rangle}$."
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "ccc0aec5-63d3-4a43-8fad-fea56ea24ec7",
+ "metadata": {},
+ "source": [
+ "## Requirements\n",
+ "\n",
+ "Before starting this tutorial, ensure that you have the following installed:\n",
+ "\n",
+ "- Qiskit SDK 2.1 or later with visualization support (`pip install 'qiskit[visualization]'`)\n",
+ "- Qiskit Runtime 0.41 or later (`pip install qiskit-ibm-runtime`)\n",
+ "- M3 Qiskit addon 3.0 (`pip install mthree`)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "c090b8dd-754f-4390-93e5-663f818e61d0",
+ "metadata": {
+ "jp-MarkdownHeadingCollapsed": true
+ },
+ "source": [
+ "## Setup"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "233dfc67-d280-4827-910d-8e7170ee95ad",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "from collections.abc import Iterator, Sequence\n",
+ "from random import Random\n",
+ "from qiskit.circuit import (\n",
+ " CircuitInstruction,\n",
+ " QuantumCircuit,\n",
+ " QuantumRegister,\n",
+ " Qubit,\n",
+ ")\n",
+ "from qiskit.circuit.library import CZGate, HGate, XGate\n",
+ "from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager\n",
+ "from qiskit_ibm_runtime import QiskitRuntimeService\n",
+ "import timeit\n",
+ "import matplotlib.pyplot as plt\n",
+ "from qiskit_ibm_runtime import SamplerV2 as Sampler\n",
+ "import mthree"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "5cf89647-261d-4af1-af54-75548345df6e",
+ "metadata": {},
+ "source": [
+ "## Step 1: Map classical inputs to a quantum problem"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "940cc113-5cf3-4688-9747-93268be10088",
+ "metadata": {},
+ "source": [
+ "First, we write the functions to implement the hidden shift problem as a `QuantumCircuit`."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "675fe72b-6a86-42ea-a65f-b29a60725ecb",
+ "metadata": {
+ "editable": true,
+ "scrolled": true,
+ "slideshow": {
+ "slide_type": ""
+ },
+ "tags": []
+ },
+ "outputs": [],
+ "source": [
+ "def apply_hadamards(qubits: Sequence[Qubit]) -> Iterator[CircuitInstruction]:\n",
+ " \"\"\"Apply a Hadamard gate to every qubit.\"\"\"\n",
+ " for q in qubits:\n",
+ " yield CircuitInstruction(HGate(), [q], [])\n",
+ "\n",
+ "\n",
+ "def apply_shift(\n",
+ " qubits: Sequence[Qubit], shift: int\n",
+ ") -> Iterator[CircuitInstruction]:\n",
+ " \"\"\"Apply X gates where the bits of the shift are equal to 1.\"\"\"\n",
+ " for i, q in zip(range(shift.bit_length()), qubits):\n",
+ " if shift >> i & 1:\n",
+ " yield CircuitInstruction(XGate(), [q], [])\n",
+ "\n",
+ "\n",
+ "def oracle_f(qubits: Sequence[Qubit]) -> Iterator[CircuitInstruction]:\n",
+ " \"\"\"Apply the f oracle.\"\"\"\n",
+ " for i in range(0, len(qubits) - 1, 2):\n",
+ " yield CircuitInstruction(CZGate(), [qubits[i], qubits[i + 1]])\n",
+ "\n",
+ "\n",
+ "def oracle_g(\n",
+ " qubits: Sequence[Qubit], shift: int\n",
+ ") -> Iterator[CircuitInstruction]:\n",
+ " \"\"\"Apply the g oracle.\"\"\"\n",
+ " yield from apply_shift(qubits, shift)\n",
+ " yield from oracle_f(qubits)\n",
+ " yield from apply_shift(qubits, shift)\n",
+ "\n",
+ "\n",
+ "def determine_hidden_shift(\n",
+ " qubits: Sequence[Qubit], shift: int\n",
+ ") -> Iterator[CircuitInstruction]:\n",
+ " \"\"\"Determine the hidden shift.\"\"\"\n",
+ " yield from apply_hadamards(qubits)\n",
+ " yield from oracle_g(qubits, shift)\n",
+ " # We omit this layer in exchange for post processing\n",
+ " # yield from apply_hadamards(qubits)\n",
+ " yield from oracle_f(qubits)\n",
+ " yield from apply_hadamards(qubits)\n",
+ "\n",
+ "\n",
+ "def run_hidden_shift_circuit(n_qubits, rng):\n",
+ " hidden_shift = rng.getrandbits(n_qubits)\n",
+ "\n",
+ " qubits = QuantumRegister(n_qubits, name=\"q\")\n",
+ " circuit = QuantumCircuit.from_instructions(\n",
+ " determine_hidden_shift(qubits, hidden_shift), qubits=qubits\n",
+ " )\n",
+ " circuit.measure_all()\n",
+ " # Format the hidden shift as a string.\n",
+ " hidden_shift_string = format(hidden_shift, f\"0{n_qubits}b\")\n",
+ " return (circuit, hidden_shift, hidden_shift_string)\n",
+ "\n",
+ "\n",
+ "def display_circuit(circuit):\n",
+ " return circuit.remove_final_measurements(inplace=False).draw(\n",
+ " \"mpl\", idle_wires=False, scale=0.5, fold=-1\n",
+ " )"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "6960aa96-f612-40b3-9f26-3d3afb07262b",
+ "metadata": {},
+ "source": [
+ "We'll start with a fairly small example:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 2,
+ "id": "8297843e-00c3-4bb5-9d33-a7e558d1698c",
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Hidden shift string 011010\n"
+ ]
+ },
+ {
+ "data": {
+ "text/plain": [
+ ""
+ ]
+ },
+ "execution_count": 2,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "n_qubits = 6\n",
+ "random_seed = 12345\n",
+ "rng = Random(random_seed)\n",
+ "circuit, hidden_shift, hidden_shift_string = run_hidden_shift_circuit(\n",
+ " n_qubits, rng\n",
+ ")\n",
+ "\n",
+ "print(f\"Hidden shift string {hidden_shift_string}\")\n",
+ "\n",
+ "display_circuit(circuit)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "c136f784-e105-4244-a483-b8221577afc3",
+ "metadata": {
+ "editable": true,
+ "slideshow": {
+ "slide_type": ""
+ },
+ "tags": []
+ },
+ "source": [
+ "## Step 2: Optimize circuits for quantum hardware execution"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 3,
+ "id": "ee4d384a-15e2-4a3f-a52a-def3d184c4fa",
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "['shift 011010', 'n_qubits 6', 'seed = 12345']"
+ ]
+ },
+ "execution_count": 3,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "job_tags = [\n",
+ " f\"shift {hidden_shift_string}\",\n",
+ " f\"n_qubits {n_qubits}\",\n",
+ " f\"seed = {random_seed}\",\n",
+ "]\n",
+ "job_tags"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "f2b77d93-c34a-43a4-b436-e7a25024a94a",
+ "metadata": {
+ "editable": true,
+ "scrolled": true,
+ "slideshow": {
+ "slide_type": ""
+ },
+ "tags": []
+ },
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Using backend ibm_kingston\n"
+ ]
+ },
+ {
+ "data": {
+ "text/plain": [
+ ""
+ ]
+ },
+ "execution_count": 4,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "# Uncomment this to run the circuits on a quantum computer on IBMCloud.\n",
+ "service = QiskitRuntimeService()\n",
+ "backend = service.backend(\"ibm_kingston\")\n",
+ "\n",
+ "# from qiskit_ibm_runtime.fake_provider import FakeMelbourneV2\n",
+ "# backend = FakeMelbourneV2()\n",
+ "# backend.refresh(service)\n",
+ "\n",
+ "print(f\"Using backend {backend.name}\")\n",
+ "\n",
+ "\n",
+ "def get_isa_circuit(circuit, backend):\n",
+ " pass_manager = generate_preset_pass_manager(\n",
+ " optimization_level=3, backend=backend, seed_transpiler=1234\n",
+ " )\n",
+ " isa_circuit = pass_manager.run(circuit)\n",
+ " return isa_circuit\n",
+ "\n",
+ "\n",
+ "isa_circuit = get_isa_circuit(circuit, backend)\n",
+ "display_circuit(isa_circuit)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "f2f426da-3601-4442-9593-16bb19abf91e",
+ "metadata": {},
+ "source": [
+ "## Step 3: Execute circuits using Qiskit Primitives"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "6ebaea63-e352-4d7c-a7ee-57a2f0a66306",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# submit job for solving the hidden shift problem using the Sampler primitive\n",
+ "NUM_SHOTS = 50_000\n",
+ "\n",
+ "\n",
+ "def run_sampler(backend, isa_circuit, num_shots):\n",
+ " sampler = Sampler(mode=backend)\n",
+ " sampler.options.environment.job_tags\n",
+ " pubs = [(isa_circuit, None, NUM_SHOTS)]\n",
+ " job = sampler.run(pubs)\n",
+ " return job\n",
+ "\n",
+ "\n",
+ "def setup_mthree_mitigation(isa_circuit, backend):\n",
+ " # retrieve the final qubit mapping so mthree knows which qubits to calibrate\n",
+ " qubit_mapping = mthree.utils.final_measurement_mapping(isa_circuit)\n",
+ "\n",
+ " # submit jobs for readout error calibration\n",
+ " mit = mthree.M3Mitigation(backend)\n",
+ " mit.cals_from_system(qubit_mapping, rep_delay=None)\n",
+ "\n",
+ " return mit, qubit_mapping"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 6,
+ "id": "9b76b943-fd71-4287-804b-980cacc078f8",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "job = run_sampler(backend, isa_circuit, NUM_SHOTS)\n",
+ "mit, qubit_mapping = setup_mthree_mitigation(isa_circuit, backend)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "674ad6e1-c496-44ca-9598-470edf4c15dc",
+ "metadata": {},
+ "source": [
+ "## Step 4: Post-process and return results in classical format\n",
+ "\n",
+ "In the theoretical discussion above, we determined that for input $ab$, we expect output $ba$.\n",
+ "An additional complication is that, in order to have a simpler (pre-transpiled) circuit, we inserted the required CZ gates between\n",
+ "neighboring pairs of qubits. This amounts to interleaving the bitstrings $a$ and $b$ as $a1 b1 a2 b2 \\ldots$.\n",
+ "The output string $ba$ will be interleaved in a similar way: $b1 a1 b2 a2 \\ldots$. The function `unscramble` below\n",
+ "transforms the output string from $b1 a1 b2 a2 \\ldots$ to $a1 b1 a2 b2 \\ldots$ so that the input and output strings may be compared directly."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 7,
+ "id": "aada9b41-0600-4814-a5a3-917d7fb171f9",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "# retrieve bitstring counts\n",
+ "def get_bitstring_counts(job):\n",
+ " result = job.result()\n",
+ " pub_result = result[0]\n",
+ " counts = pub_result.data.meas.get_counts()\n",
+ " return counts, pub_result"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 8,
+ "id": "3930e290-89ff-404d-946b-cdf90119cb48",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "counts, pub_result = get_bitstring_counts(job)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "138ef822-8e66-4f01-94f6-5a634465d2ac",
+ "metadata": {},
+ "source": [
+ "The Hamming distance between two bit strings is the number of indices at which the bits differ."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 9,
+ "id": "3488f07f-7705-4ef2-b033-ebc2cbe19d45",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "def hamming_distance(s1, s2):\n",
+ " weight = 0\n",
+ " for c1, c2 in zip(s1, s2):\n",
+ " (c1, c2) = (int(c1), int(c2))\n",
+ " if (c1 == 1 and c2 == 1) or (c1 == 0 and c2 == 0):\n",
+ " weight += 1\n",
+ "\n",
+ " return weight"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 10,
+ "id": "9e2f66e7-cad3-4893-9fde-ddaec077d938",
+ "metadata": {
+ "scrolled": true
+ },
+ "outputs": [],
+ "source": [
+ "# Replace string of form a1b1a2b2... with b1a1b2a1...\n",
+ "# That is, reverse order of successive pairs of bits.\n",
+ "def unscramble(bitstring):\n",
+ " ps = [bitstring[i : i + 2][::-1] for i in range(0, len(bitstring), 2)]\n",
+ " return \"\".join(ps)\n",
+ "\n",
+ "\n",
+ "def find_hidden_shift_bitstring(counts, hidden_shift_string):\n",
+ " # convert counts to probabilities\n",
+ " probs = {\n",
+ " unscramble(bitstring): count / NUM_SHOTS\n",
+ " for bitstring, count in counts.items()\n",
+ " }\n",
+ "\n",
+ " # Retrieve the most probable bitstring.\n",
+ " most_probable = max(probs, key=lambda x: probs[x])\n",
+ "\n",
+ " print(f\"Expected hidden shift string: {hidden_shift_string}\")\n",
+ " if most_probable == hidden_shift_string:\n",
+ " print(\"Most probable bitstring matches hidden shift π.\")\n",
+ " else:\n",
+ " print(\"Most probable bitstring didn't match hidden shift βΉοΈ.\")\n",
+ " print(\"Top 10 bitstrings and their probabilities:\")\n",
+ " display(\n",
+ " {\n",
+ " k: (v, hamming_distance(hidden_shift_string, k))\n",
+ " for k, v in sorted(\n",
+ " probs.items(), key=lambda x: x[1], reverse=True\n",
+ " )[:10]\n",
+ " }\n",
+ " )\n",
+ "\n",
+ " return probs, most_probable"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 11,
+ "id": "9d3dd975-b628-48fc-8565-86b0ed88c973",
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Expected hidden shift string: 011010\n",
+ "Most probable bitstring matches hidden shift π.\n",
+ "Top 10 bitstrings and their probabilities:\n"
+ ]
+ },
+ {
+ "data": {
+ "text/plain": [
+ "{'011010': (0.9743, 6),\n",
+ " '001010': (0.00812, 5),\n",
+ " '010010': (0.0063, 5),\n",
+ " '011000': (0.00554, 5),\n",
+ " '011011': (0.00492, 5),\n",
+ " '011110': (0.00044, 5),\n",
+ " '001000': (0.00012, 4),\n",
+ " '010000': (8e-05, 4),\n",
+ " '001011': (6e-05, 4),\n",
+ " '000010': (6e-05, 4)}"
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ }
+ ],
+ "source": [
+ "probs, most_probable = find_hidden_shift_bitstring(\n",
+ " counts, hidden_shift_string\n",
+ ")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "67a88958-c900-4bd1-8329-d530ba3427ea",
+ "metadata": {},
+ "source": [
+ "Let's record the probability of the most probable bitstring before applying readout error mitigation with M3."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 12,
+ "id": "a95daf38-009e-44ea-a9c5-8041bf7b22e2",
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "0.9743"
+ ]
+ },
+ "execution_count": 12,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "max_probability_before_M3 = probs[most_probable]\n",
+ "max_probability_before_M3"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "9531e7bd-e134-4110-ba00-d5a08fe612d9",
+ "metadata": {},
+ "source": [
+ "Now we apply the readout correction learned by M3 to the counts.\n",
+ "The function `apply_corrections` returns a quasi-probability distribution. This is a list of `float`s that sum to one. But some values might be negative."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 13,
+ "id": "3a24d85c-a980-48cf-96ab-5237d370d7af",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "def perform_mitigation(mit, counts, qubit_mapping):\n",
+ " # mitigate readout error\n",
+ " quasis = mit.apply_correction(counts, qubit_mapping)\n",
+ "\n",
+ " # print results\n",
+ " most_probable_after_m3 = unscramble(max(quasis, key=lambda x: quasis[x]))\n",
+ "\n",
+ " is_hidden_shift_identified = most_probable_after_m3 == hidden_shift_string\n",
+ " if is_hidden_shift_identified:\n",
+ " print(\"Most probable bitstring matches hidden shift π.\")\n",
+ " else:\n",
+ " print(\"Most probable bitstring didn't match hidden shift βΉοΈ.\")\n",
+ " print(\"Top 10 bitstrings and their quasi-probabilities:\")\n",
+ " topten = {\n",
+ " unscramble(k): f\"{v:.2e}\"\n",
+ " for k, v in sorted(quasis.items(), key=lambda x: x[1], reverse=True)[\n",
+ " :10\n",
+ " ]\n",
+ " }\n",
+ " max_probability_after_M3 = float(topten[most_probable_after_m3])\n",
+ " display(topten)\n",
+ "\n",
+ " return max_probability_after_M3, is_hidden_shift_identified"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 14,
+ "id": "d304c727-f163-4842-a96a-dad8a067c8d1",
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Expected hidden shift string: 011010\n",
+ "Most probable bitstring matches hidden shift π.\n",
+ "Top 10 bitstrings and their quasi-probabilities:\n"
+ ]
+ },
+ {
+ "data": {
+ "text/plain": [
+ "{'011010': '1.01e+00',\n",
+ " '001010': '8.75e-04',\n",
+ " '001000': '7.38e-05',\n",
+ " '010000': '4.51e-05',\n",
+ " '111000': '2.18e-05',\n",
+ " '001011': '1.74e-05',\n",
+ " '000010': '6.42e-06',\n",
+ " '011001': '-7.18e-06',\n",
+ " '011000': '-4.53e-04',\n",
+ " '010010': '-1.28e-03'}"
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ }
+ ],
+ "source": [
+ "print(f\"Expected hidden shift string: {hidden_shift_string}\")\n",
+ "max_probability_after_M3, is_hidden_shift_identified = perform_mitigation(\n",
+ " mit, counts, qubit_mapping\n",
+ ")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "b62fac82-e883-4030-9793-86af8299abaa",
+ "metadata": {},
+ "source": [
+ "#### Comparing identifying the hidden shift string before and after applying M3 correction"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 15,
+ "id": "a3864692-582e-4fe0-8059-f0be59cfd229",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "def compare_before_and_after_M3(\n",
+ " max_probability_before_M3,\n",
+ " max_probability_after_M3,\n",
+ " is_hidden_shift_identified,\n",
+ "):\n",
+ " is_probability_improved = (\n",
+ " max_probability_after_M3 > max_probability_before_M3\n",
+ " )\n",
+ " print(f\"Most probable probability before M3: {max_probability_before_M3}\")\n",
+ " print(f\"Most probable probability after M3: {max_probability_after_M3}\")\n",
+ " if is_hidden_shift_identified and is_probability_improved:\n",
+ " print(\"Readout error mitigation effective! π\")\n",
+ " else:\n",
+ " print(\"Readout error mitigation not effective. βΉοΈ\")"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 16,
+ "id": "98f92ec2-50b5-4456-b40a-3f619cc42ff1",
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Most probable probability before M3: 0.9743\n",
+ "Most probable probability after M3: 1.01\n",
+ "Readout error mitigation effective! π\n"
+ ]
+ }
+ ],
+ "source": [
+ "compare_before_and_after_M3(\n",
+ " max_probability_before_M3,\n",
+ " max_probability_after_M3,\n",
+ " is_hidden_shift_identified,\n",
+ ")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "90c3bec4-9a3d-43b7-8160-bcb615230d01",
+ "metadata": {},
+ "source": [
+ "### Plot how CPU time required by M3 scales with shots"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "33addc38-f738-48ed-a29d-9790f446c036",
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Applying M3 correction to 5000 shots...\n",
+ "\tDone in 0.003321983851492405 seconds.\n",
+ "Applying M3 correction to 7500 shots...\n",
+ "\tDone in 0.004425413906574249 seconds.\n",
+ "Applying M3 correction to 10000 shots...\n",
+ "\tDone in 0.006366567220538855 seconds.\n",
+ "Applying M3 correction to 12500 shots...\n",
+ "\tDone in 0.0071477219462394714 seconds.\n",
+ "Applying M3 correction to 15000 shots...\n",
+ "\tDone in 0.00860048783943057 seconds.\n",
+ "Applying M3 correction to 17500 shots...\n",
+ "\tDone in 0.010026784148067236 seconds.\n",
+ "Applying M3 correction to 20000 shots...\n",
+ "\tDone in 0.011459112167358398 seconds.\n",
+ "Applying M3 correction to 22500 shots...\n",
+ "\tDone in 0.012727141845971346 seconds.\n",
+ "Applying M3 correction to 25000 shots...\n",
+ "\tDone in 0.01406092382967472 seconds.\n",
+ "Applying M3 correction to 27500 shots...\n",
+ "\tDone in 0.01546052098274231 seconds.\n",
+ "Applying M3 correction to 30000 shots...\n",
+ "\tDone in 0.016769016161561012 seconds.\n",
+ "Applying M3 correction to 32500 shots...\n",
+ "\tDone in 0.019537431187927723 seconds.\n",
+ "Applying M3 correction to 35000 shots...\n",
+ "\tDone in 0.019739801064133644 seconds.\n",
+ "Applying M3 correction to 37500 shots...\n",
+ "\tDone in 0.021093040239065886 seconds.\n",
+ "Applying M3 correction to 40000 shots...\n",
+ "\tDone in 0.022840639110654593 seconds.\n",
+ "Applying M3 correction to 42500 shots...\n",
+ "\tDone in 0.023974396288394928 seconds.\n",
+ "Applying M3 correction to 45000 shots...\n",
+ "\tDone in 0.026412792038172483 seconds.\n",
+ "Applying M3 correction to 47500 shots...\n",
+ "\tDone in 0.026364430785179138 seconds.\n",
+ "Applying M3 correction to 50000 shots...\n",
+ "\tDone in 0.02820305060595274 seconds.\n"
+ ]
+ },
+ {
+ "data": {
+ "text/plain": [
+ "Text(0.5, 1.0, 'Time to apply M3 correction')"
+ ]
+ },
+ "execution_count": 17,
+ "metadata": {},
+ "output_type": "execute_result"
+ },
+ {
+ "data": {
+ "text/plain": [
+ ""
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ }
+ ],
+ "source": [
+ "# Collect samples for numbers of shots varying from 5000 to 25000.\n",
+ "shots_range = range(5000, NUM_SHOTS + 1, 2500)\n",
+ "times = []\n",
+ "for shots in shots_range:\n",
+ " print(f\"Applying M3 correction to {shots} shots...\")\n",
+ " t0 = timeit.default_timer()\n",
+ " _ = mit.apply_correction(\n",
+ " pub_result.data.meas.slice_shots(range(shots)).get_counts(),\n",
+ " qubit_mapping,\n",
+ " )\n",
+ " t1 = timeit.default_timer()\n",
+ " print(f\"\\tDone in {t1 - t0} seconds.\")\n",
+ " times.append(t1 - t0)\n",
+ "\n",
+ "fig, ax = plt.subplots()\n",
+ "ax.plot(shots_range, times, \"o--\")\n",
+ "ax.set_xlabel(\"Shots\")\n",
+ "ax.set_ylabel(\"Time (s)\")\n",
+ "ax.set_title(\"Time to apply M3 correction\")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "bdfbeddd-e90e-4a5c-ad19-810603c4f695",
+ "metadata": {},
+ "source": [
+ "#### Interpreting the plot\n",
+ "\n",
+ "The plot above shows that the time required to apply M3 correction scales linearly in the number of shots."
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "4988be5c-46ad-473d-b1b1-c8280d8f6324",
+ "metadata": {},
+ "source": [
+ "# Scaling up!"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 18,
+ "id": "1bf46124-9261-472e-864f-ee2ed0209079",
+ "metadata": {
+ "jp-MarkdownHeadingCollapsed": true
+ },
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Hidden shift string 00000010100110101011101110010001010000110011101001101010101001111001100110000111\n"
+ ]
+ }
+ ],
+ "source": [
+ "n_qubits = 80\n",
+ "rng = Random(12345)\n",
+ "circuit, hidden_shift, hidden_shift_string = run_hidden_shift_circuit(\n",
+ " n_qubits, rng\n",
+ ")\n",
+ "\n",
+ "print(f\"Hidden shift string {hidden_shift_string}\")"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 19,
+ "id": "53cd9086-27f8-481b-ba96-9bac9b146188",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "isa_circuit = get_isa_circuit(circuit, backend)"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 20,
+ "id": "b176b2b3-de30-4a9b-905d-a06c7365376a",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "job = run_sampler(backend, isa_circuit, NUM_SHOTS)\n",
+ "mit, qubit_mapping = setup_mthree_mitigation(isa_circuit, backend)"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 21,
+ "id": "c4dda06e-e2c6-46b9-9e74-d90bb318419b",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "counts, pub_result = get_bitstring_counts(job)"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 22,
+ "id": "3421d2f3-083f-44a9-bee4-5f05e43ebb4f",
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Expected hidden shift string: 00000010100110101011101110010001010000110011101001101010101001111001100110000111\n",
+ "Most probable bitstring matches hidden shift π.\n",
+ "Top 10 bitstrings and their probabilities:\n"
+ ]
+ },
+ {
+ "data": {
+ "text/plain": [
+ "{'00000010100110101011101110010001010000110011101001101010101001111001100110000111': (0.50402,\n",
+ " 80),\n",
+ " '00000010100110101011101110010001010000110011100001101010101001111001100110000111': (0.0396,\n",
+ " 79),\n",
+ " '00000010100110101011101110010001010000110011101001101010101001111001100100000111': (0.0323,\n",
+ " 79),\n",
+ " '00000010100110101011101110010001010000110011101001101010101001101001100110000111': (0.01936,\n",
+ " 79),\n",
+ " '00000010100110101011101110010011010000110011101001101010101001111001100110000111': (0.01432,\n",
+ " 79),\n",
+ " '00000010100110101011101110010001010000110011101001101010101001011001100110000111': (0.0101,\n",
+ " 79),\n",
+ " '00000010100110101011101110010001010000110011101001101010101001110001100110000111': (0.00924,\n",
+ " 79),\n",
+ " '00000010100110101011101110010001010000010011101001101010101001111001100110000111': (0.00908,\n",
+ " 79),\n",
+ " '00000010100110101011100110010001010000110011101001101010101001111001100110000111': (0.00888,\n",
+ " 79),\n",
+ " '00000010100110101011101110010001010000110011101001100010101001111001100110000111': (0.0082,\n",
+ " 79)}"
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ }
+ ],
+ "source": [
+ "probs, most_probable = find_hidden_shift_bitstring(\n",
+ " counts, hidden_shift_string\n",
+ ")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "2ec8335f-d3ef-4951-8bbd-b4d01f899be1",
+ "metadata": {},
+ "source": [
+ "We see that the correct hidden shift string is found. Furthermore, the nine next-most-probably bitstrings are wrong in only one position."
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "98d1fb7b-e696-4cef-a8b8-18e44b96eece",
+ "metadata": {},
+ "source": [
+ "Record the most probable probability"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 23,
+ "id": "c5a9dbdc-ebb7-470a-ac31-10b41217eeed",
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/plain": [
+ "0.50402"
+ ]
+ },
+ "execution_count": 23,
+ "metadata": {},
+ "output_type": "execute_result"
+ }
+ ],
+ "source": [
+ "max_probability_before_M3 = probs[most_probable]\n",
+ "max_probability_before_M3"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 24,
+ "id": "b4514b88-3f2d-4e96-88ce-6d291576568d",
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Expected hidden shift string: 00000010100110101011101110010001010000110011101001101010101001111001100110000111\n",
+ "Most probable bitstring matches hidden shift π.\n",
+ "Top 10 bitstrings and their quasi-probabilities:\n"
+ ]
+ },
+ {
+ "data": {
+ "text/plain": [
+ "{'00000010100110101011101110010001010000110011101001101010101001111001100110000111': '9.85e-01',\n",
+ " '00000010100110101011101110010001010000110011100001101010101001111001100110000111': '6.84e-03',\n",
+ " '00000010100110101011100110010001010000110011101001101010101001111001100110000111': '3.87e-03',\n",
+ " '00000010100110101011101110010011010000110011101001101010101001111001100110000111': '3.42e-03',\n",
+ " '00000010100110101011101110010001010000110011101001101010101001111001100100000111': '3.30e-03',\n",
+ " '00000010100110101011101110010001010000110011101001101010101001110001100110000111': '3.28e-03',\n",
+ " '00000010100010101011101110010001010000110011101001101010101001111001100110000111': '2.62e-03',\n",
+ " '00000010100110101011101110010001010000110011101001101010101001101001100110000111': '2.43e-03',\n",
+ " '00000010100110101011101110010000010000110011101001101010101001111001100110000111': '1.73e-03',\n",
+ " '00000010100110101011101110010001010000110011101001101010101001111001000110000111': '1.63e-03'}"
+ ]
+ },
+ "metadata": {},
+ "output_type": "display_data"
+ }
+ ],
+ "source": [
+ "print(f\"Expected hidden shift string: {hidden_shift_string}\")\n",
+ "max_probability_after_M3, is_hidden_shift_identified = perform_mitigation(\n",
+ " mit, counts, qubit_mapping\n",
+ ")"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 24,
+ "id": "ccd14379-6c34-4be8-93ab-728356bf81e0",
+ "metadata": {},
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Most probable probability before M3: 0.54348\n",
+ "Most probable probability after M3: 0.99\n",
+ "Readout error mitigation effective! π\n"
+ ]
+ }
+ ],
+ "source": [
+ "compare_before_and_after_M3(\n",
+ " max_probability_before_M3,\n",
+ " max_probability_after_M3,\n",
+ " is_hidden_shift_identified,\n",
+ ")"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "96b6648b-5c26-4a0c-92f2-2d1810fa14b8",
+ "metadata": {},
+ "source": [
+ "The results show that readout error was the dominant source of error and the M3 mitigation was effective."
+ ]
+ }
+ ],
+ "metadata": {
+ "description": "Use the M3 readout mitigation addon with the Sampler primitive",
+ "kernelspec": {
+ "display_name": "Python 3",
+ "language": "python",
+ "name": "python3"
+ },
+ "language_info": {
+ "codemirror_mode": {
+ "name": "ipython",
+ "version": 3
+ },
+ "file_extension": ".py",
+ "mimetype": "text/x-python",
+ "name": "python",
+ "nbconvert_exporter": "python",
+ "pygments_lexer": "ipython3",
+ "version": "3"
+ },
+ "title": "Readout error mitigation for the Sampler primitive using M3",
+ "toc": {
+ "base_numbering": 0
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 4
+}
diff --git a/public/docs/images/tutorials/readout-error-mitigation-sampler/extracted-outputs/33addc38-f738-48ed-a29d-9790f446c036-2.avif b/public/docs/images/tutorials/readout-error-mitigation-sampler/extracted-outputs/33addc38-f738-48ed-a29d-9790f446c036-2.avif
new file mode 100644
index 00000000000..8b3aa0fc4ae
Binary files /dev/null and b/public/docs/images/tutorials/readout-error-mitigation-sampler/extracted-outputs/33addc38-f738-48ed-a29d-9790f446c036-2.avif differ
diff --git a/public/docs/images/tutorials/readout-error-mitigation-sampler/extracted-outputs/8297843e-00c3-4bb5-9d33-a7e558d1698c-1.avif b/public/docs/images/tutorials/readout-error-mitigation-sampler/extracted-outputs/8297843e-00c3-4bb5-9d33-a7e558d1698c-1.avif
new file mode 100644
index 00000000000..a6c77b0c357
Binary files /dev/null and b/public/docs/images/tutorials/readout-error-mitigation-sampler/extracted-outputs/8297843e-00c3-4bb5-9d33-a7e558d1698c-1.avif differ
diff --git a/public/docs/images/tutorials/readout-error-mitigation-sampler/extracted-outputs/f2b77d93-c34a-43a4-b436-e7a25024a94a-1.avif b/public/docs/images/tutorials/readout-error-mitigation-sampler/extracted-outputs/f2b77d93-c34a-43a4-b436-e7a25024a94a-1.avif
new file mode 100644
index 00000000000..a8b1ed57412
Binary files /dev/null and b/public/docs/images/tutorials/readout-error-mitigation-sampler/extracted-outputs/f2b77d93-c34a-43a4-b436-e7a25024a94a-1.avif differ
diff --git a/qiskit_bot.yaml b/qiskit_bot.yaml
index 7218a9cb80e..1c035202747 100644
--- a/qiskit_bot.yaml
+++ b/qiskit_bot.yaml
@@ -683,6 +683,8 @@ notifications:
- "@annaliese-estes"
"docs/tutorials/quantum-phase-estimation-qctrl":
- "@dsierrasosa"
+ "docs/tutorials/readout-error-mitigation-sampler":
+ - "@jlapeyre"
"learning/courses/quantum-computing-in-practice/introduction":
- "@livlanes"
"learning/courses/quantum-computing-in-practice/running-quantum-circuits":
diff --git a/scripts/config/notebook-testing.toml b/scripts/config/notebook-testing.toml
index 12e4465dd70..b09a59942e8 100644
--- a/scripts/config/notebook-testing.toml
+++ b/scripts/config/notebook-testing.toml
@@ -162,6 +162,7 @@ notebooks = [
"docs/guides/qiskit-transpiler-service.ipynb",
# We never run tutorials notebooks
+ "docs/tutorials/readout-error-mitigation-sampler.ipynb",
"docs/tutorials/quantum-phase-estimation-qctrl.ipynb",
"docs/tutorials/transverse-field-ising-model.ipynb",
"docs/tutorials/fractional-gates.ipynb",