# Near-Term Quantum Computing Techniques: Variational Quantum Algorithms,
Error Mitigation, Circuit Compilation, Benchmarking and Classical Simulation

He-Liang Huang#<EMAIL_ADDRESS>Henan Key Laboratory of Quantum
Information and Cryptography, Zhengzhou, Henan 450000, China  Xiao-Yue Xu#
Henan Key Laboratory of Quantum Information and Cryptography, Zhengzhou, Henan
450000, China  Chu Guo# Henan Key Laboratory of Quantum Information and
Cryptography, Zhengzhou, Henan 450000, China  Guojing Tian# State Key Lab of
Processors, Institute of Computing Technology, Chinese Academy of Sciences,
100190 Beijing, China  School of Computer Science and Technology, University
of Chinese Academy of Sciences, Beijing 100049, China  Shi-Jie Wei  Beijing
Academy of Quantum Information Sciences, Beijing 100193, China  Xiaoming Sun
<EMAIL_ADDRESS>State Key Lab of Processors, Institute of Computing
Technology, Chinese Academy of Sciences, 100190 Beijing, China  School of
Computer Science and Technology, University of Chinese Academy of Sciences,
Beijing 100049, China  Wan-Su Bao<EMAIL_ADDRESS>Henan Key Laboratory of
Quantum Information and Cryptography, Zhengzhou, Henan 450000, China  Gui-Lu
Long<EMAIL_ADDRESS>Department of Physics, Tsinghua University,
Beijing 100084, China  Beijing Academy of Quantum Information Sciences,
Beijing 100193, China

###### Abstract

Quantum computing is a game-changing technology for global academia, research
centers and industries including computational science, mathematics, finance,
pharmaceutical, materials science, chemistry and cryptography. Although it has
seen a major boost in the last decade, we are still a long way from reaching
the maturity of a full-fledged quantum computer. That said, we will be in the
Noisy-Intermediate Scale Quantum (NISQ) era for a long time, working on dozens
or even thousands of qubits quantum computing systems. An outstanding
challenge, then, is to come up with an application that can reliably carry out
a nontrivial task of interest on the near-term quantum devices with non-
negligible quantum noise. To address this challenge, several near-term quantum
computing techniques, including variational quantum algorithms, error
mitigation, quantum circuit compilation and benchmarking protocols, have been
proposed to characterize and mitigate errors, and to implement algorithms with
a certain resistance to noise, so as to enhance the capabilities of near-term
quantum devices and explore the boundaries of their ability to realize useful
applications. Besides, the development of near-term quantum devices is
inseparable from the efficient classical simulation, which plays a vital role
in quantum algorithm design and verification, error-tolerant verification and
other applications. This review will provide a thorough introduction of these
near-term quantum computing techniques, report on their progress, and finally
discuss the future prospect of these techniques, which we hope will motivate
researchers to undertake additional studies in this field.

# #footnotetext:  These four authors contributed equally.

##  I Introduction

Quantum computing is a rapidly emerging new-generation computing paradigm that
harnesses the laws of quantum mechanics, offering the potential of exponential
speedup over classical computation for certain problems [1, 2, 3, 4]. Over the
last decade, quantum computing technologies have advanced by leaps and bounds
into the NISQ era [5, 6], an important sign of which is the significant
achievement of quantum computational advantage on quantum sampling problems
has been realized in practice [7, 8, 9, 10, 11, 12]. The near-term quantum
processors produced during this period contain only a few dozen or even a few
thousand qubits, falling well short of the specifications for fault-tolerant
quantum computing [13, 14, 15, 16, 17]. They do, however, serve as fantastic
testing grounds for investigating a variety of applications, including machine
learning [18], secure cloud quantum computing [19, 20, 21, 22, 23],
computational science, and complicated quantum chemistry [24] and many-body
quantum systems [25] that are not feasible to be simulated with the state-of-
the-art supercomputers.

To fully exploit the potential of near-term quantum devices,
algorithms/protocols have to be tailored to the constraint of present quantum
hardware, especially a modest of qubits with non-negligible errors and limited
qubit connectivity. To mitigate quantum errors and achieve valid computational
results, several prominent classes of algorithms and tools have been developed
specifically for near-term quantum computers, including:

  * •

Variational Quantum Algorithms (VQAs) [26, 27]. VQAs have emerged as one of
the leading candidates to achieve application-oriented quantum computational
advantage on NISQ devices, owing to their hybrid quantum-classical approach
which has potential noise resilience.

  * •

Quantum Error Mitigation (QEM) [28]. QEM refers to a series of techniques that
allow us to reduce the computational errors and then evaluate accurate results
from noisy quantum circuits, although we still lack practical quantum error-
correction technique.

  * •

Quantum Circuit Compilation (QCC) [29, 30, 31]. QCC is a key technique to
transform the nonconforming quantum circuit to an executable circuit on the
target quantum platform according to its constraints, including the native
gateset, connectivity and so on.

  * •

Benchmarking Protocols [32, 33]. Quantum benchmarking helps to evaluate the
basic performance of a quantum computer and even the capacity to solve real-
world problems.

  * •

Classical Simulation [34, 35, 36, 37, 38]. Classical simulation of quantum
circuits is one of the core tools for designing quantum algorithms and
validating quantum devices.

These techniques work together to accomplish computational tasks rather than
independently of each other (see a rough relationship provided in Fig. 1). For
instance, variational quantum algorithms, quantum error mitigation, quantum
circuit compilation, and quantum benchmarking may all require the help of
classical simulation for verification or algorithm design. These near-term
quantum computing techniques offer the opportunity to attain practical quantum
advantages that can be applied to significant applications such as chemistry
[39, 40] and machine learning [41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52,
53, 54], and they also provide a continuous path for the advancement of
quantum computing from today’s near-term quantum hardware to tomorrow’s fault-
tolerant quantum computers.

This review aims to present an overview of these crucial techniques, part of
which may serve as a supplement and update to related published review works
[26, 27, 28, 32], and some algorithms and protocols are step-by-step taught to
help the reader rapidly understand the specifics. We also discuss the
challenges and opportunities, and offer some insights on potential future
advances.

Figure 1: Relationship between near-term quantum techniques. These techniques
collaborate to execute various computing tasks.

##  II Variational Quantum Algorithms

Although some quantum algorithms, such as Shor’s algorithm [3] and Grover’s
search algorithm [4], promise to be much faster than classical algorithms,
these algorithms are unattainable in the NISQ era due to the lack of error
correction capability. Variational quantum algorithms (VQAs) have been shown
to be naturally resistant [55] to, and even benefit from, noise, making them
particularly suitable for near-term quantum devices, and are therefore
considered the most promising route to quantum advantage on practical problems
in the NISQ era [26]. In this section, we will first introduce the basic
concepts, architectures and applications of VQAs, and then discuss the
opportunities and challenges for future industrial applications or scientific
research.

###  II.1 Basic concepts

Similar to the classical neural networks, the VQAs is implemented by training
a parametrized quantum circuit (PQC) to minimize a problem-specific cost
function. Thus, the basic components of VQAs include cost function, PQC and
optimization algorithms, as shown in Fig. 2.

Figure 2: A high-level overview of variational quantum algorithms (VQAs). VQAs
are a hybrid quantum-classical optimization algorithm. During the training
procedure, the quantum processing unit (QPU) estimates the cost function via
the parameterized quantum circuit (PQC). Then a classical processing unit
(CPU) is employed to implement the parameter optimization. VQAs can be
analogous to the quantum counterpart for deep learning, and thus have the
potential to be applied to a range of applications, such as find the ground
state, combinatorial optimization, machine learning, etc.

####  II.1.1 Cost function

The cost function (sometimes referred to as the loss or objective function) is
an integral part of VQAs for encoding problems, and serves essentially the
same role as the cost function in classical machine learning. During the
training process, the cost function can be considered as a measure of the
performance of a VQA with respect to the given training samples and the
expected output, which helps to find the global minima. Without loss of
generality, the cost function can be expressed as

$C(\bm{\theta})=f(O(U(\bm{\theta}),\rho_{\text{in}}))$ (1)

where $O$ are a set of observables on the output state obtained by the input
state $\rho_{\text{in}}$ under the action of the PQC $U(\bm{\theta})$ with
tunable parameters $\bm{\theta}$, and $f$ is some function, which is designed
according to the specific problem. For example, the cost function of
variational quantum eigensolver is often set as

$\langle H\rangle_{O(U(\bm{\theta}),\rho_{\text{in}})},$ (2)

which is expectation value of qubit Hamiltonian $H$. In addition, for quantum
machine learning tasks, one can usually follow the cost functions commonly
used in classical community, such as using mean squared error (MSE) or cross
entropy as the cost function for the classification task.

The choice of cost function could affect the trainability of VQAs. Cerezo et
al. proved that barren plateaus are cost-function-dependent in shallow PQCs
[56]. Specifically, the cost function with global observables leads to an
exponentially vanishing gradient (i.e., barren plateau), while the gradient
vanishes polynomially when the cost function is defined in terms of local
observables. Recently, Liu et al. found a similar phenomenon in tensor-network
based machine learning [57].

####  II.1.2 Parameterized quantum circuits (PQCs)

PQCs are the main part that differentiates VQAs from classical neural
networks. A PQC consists of a number of fixed quantum gates and trainable
quantum gates, and may even includes some measurement and feedback operations.
The trainable quantum gate is usually a single-qubit rotation,
$R_{x}(\theta)$, $R_{y}(\theta)$, and $R_{z}(\theta)$ rotations, and the
trainable parameter is its rotation angle $\theta$. These trainable parameters
are trained to minimize the cost function.

The structure of the variational ansatz plays an essential role in the
performance of VQA, including the convergence speed and the closeness of the
output to the optimal solution of the problem. The design of an effective
structure is a challenge for VQAs, which requires optimization in many
aspects, such as strong expressibility, shallow circuit depth, and small
number of parameters. In the following, we will present some common ansatzes.

Problem-Inspired Ansatzes. The problem-inspired ansatz is constructed by using
knowledge about a specific problem. A typical example is the quantum
approximate optimization algorithm (QAOA) [58], whose ansatz is set as

$U(\bm{\gamma},\bm{\beta})=\prod\nolimits_{l=1}^{p}{{e^{-i{\beta_{l}}{H_{M}}}}{e^{-i{\gamma_{l}}{H_{C}}}}},$
(3)

where $H_{C}$ is the problem Hamiltonian,
$H_{M}=\sum\limits_{i=1}^{n}{{\sigma_{x}^{i}}}$ is the mixing Hamiltonian, and
$\sigma_{x}^{i}$ denotes the Pauli $X$ operator acting on qubit $i$. This
ansatz is inspired by quantum annealing to map the initial state to the ground
state of the problem Hamiltonian $H_{C}$, by optimizing the parameters
$\bm{\gamma}=({\gamma_{1}},{\gamma_{2}},...,{\gamma_{p}})$ and
$\bm{\beta}=({\beta_{1}},{\beta_{2}},...,{\beta_{p}})$.

The problem-inspired ansatz is also widely used in quantum chemistry problems.
The unitary coupled cluster (UCC) [59] ansatz is a commonly used structure for
variational quantum eigensolver (VQE). And historically, the first VQE
experiment by Peruzzo et al. [60] utilized the UCC with singles and doubles
(UCCSD) ansatz. Subsequently, some improved ansatzes enable shallower depths
and higher accuracy. Lee et al. proposed the unitary pair coupled cluster with
generalized singles and doubles (k-UpCCGSD) method to reduce circuit depth
[61]. The orbital optimized UCC (OO-UCC) ansatz [62] is another variant of UCC
to reduce the number of parameters and circuit depth, while maintains a
similar level of accuracy to that of UCCSD. Besides reducing circuit depth,
Metcalf et al. employed the double unitary coupled-cluster (DUCC) method to
effectively reduce the required number of qubits [63]. Here we only make a
list of these ansatzes, more detailed introduction can be found in the Ref.
[64].

In general, the problem-inspired ansatz is designed according to the
characteristics of the problem, so that the output has a high accuracy.
However, limited by the connectivity of quantum processors, implementing
problem-inspired ansatz on real quantum devices may necessitate a relatively
deep PQC, which might be a challenge.

Hardware-Efficient Ansatzes. Unlike the problem-inspired ansatz, the hardware-
efficient ansatz [39] is designed around properties of quantum hardware, such
as the qubit connectivity and restricted gate sets, to ensure efficient
implementation on near-term quantum devices.

The trainable unitary operator $U(\bm{\theta})$, represented by the PQC of
this anstz, is composed of $L$ layers and each layer merges trainable single-
qubit gates with a entangling block. Mathematically, we have
$U(\bm{\theta}):=\prod_{l=1}^{L}(U_{E}U_{l}(\bm{\theta}))$, where
$U_{l}(\bm{\theta})$ is the $l$-th trainable layer and $U_{E}$ is the
entanglement layer. In particular, we have
$U_{l}(\bm{\theta})=\bigotimes_{i=1}^{N}(U_{S}(\bm{\theta}^{(i,l)}))$, where
$\bm{\theta}^{(i,l)}$ represents the $(i,l)$-th entry of
$\bm{\theta}\in\mathcal{R}^{N\times L}$, $U_{S}$ is the trainable unitary with
$U_{S}\in SU(2)$, e.g., the rotation single qubit gates $R_{x}$, $R_{y}$, and
$R_{z}$. The entangle layer $U_{E}$ is constructed from the entanglement
operations that the hardware can provide, such as CNOT, CZ gates, iSWAP gate,
and even many-body quantum evolution. These entanglement operations should
follow the connectivity of quantum hardware, such as one-dimensional chains,
lattices, etc.

Variable-Structure Ansatzes. Unlike the problem-inspired ansatz and hardware-
efficient ansatz, which take a fixed structure, the circuit structure in the
variable-structure ansatz is no longer fixed and trained as a parameter with
the aim of further reducing the circuit depth and the number of gates. A
typical example is ADAPT-VQE [65] and qubit-ADAPT-VQE [66], which are proposed
to optimize the circuit structure of the ansatz to reduce the circuit depth in
quantum chemistry applications. Zhu et al. also developed an adaptive QAOA for
solving combinatorial problems [67]. Interested readers can refer to the Refs.
[68, 69, 70, 71] for other variants or optimizations, and we will not
introduce them one by one.

####  II.1.3 Parameter optimization

The procedure of parameter optimization is implemented after the cost function
and ansatz are determined, to minimize the cost function. An appropriate
optimizer should be determined by considering the requirements of the
application. There are mainly two classes of optimization methods, gradient-
based and gradient-free methods, which we will introduce below.

Gradient Descent Methods. Gradient descent is a general-purpose optimization
algorithm for finding a local minimum of the given cost function. For this
purpose, the gradient of a cost function $C(\bm{\theta})$ with respect to its
parameters $\bm{\theta}$, i.e., ${\nabla_{\bm{\theta}}}C(\bm{\theta})$, is
calculated to obtain the direction of steepest descent. Then a small updates
to the parameters is implemented according to

$\bm{\theta}\to\bm{\theta}-\eta\nabla_{\bm{\theta}}C(\bm{\theta}),$ (4)

where $\eta$ is the step size, to iteratively approach the minimum value of
the cost function.

The most commonly used gradient descent method for VQAs is parameter-shift
rule [72, 73, 74, 75]. As shown in Eq. 1, the cost functions are usually
phrased in terms of expectation values of some observables $O$, evaluated on a
output quantum state $\rho_{\text{out}}$ of the PQC $U(\bm{\theta})$. In the
case that $\rho_{\text{out}}$ depends on a parameter $\theta_{i}$ that
parametrizes a Pauli-rotation gate $e^{-i\theta_{i}P/2}$, where $P$ is a Pauli
operator, we can compute the derivative with respect to $\theta_{i}$ as

${\nabla_{{\theta_{i}}}}\langle
O\rangle({\theta_{i}})=\frac{1}{2}\left({\langle
O\rangle({\theta_{i}}+\frac{\pi}{2})-\langle
O\rangle({\theta_{i}}-\frac{\pi}{2})}\right).$ (5)

The original parameter-shift rule is design for single-parameter gates.
Wierichs et al. then extended this method to be applicable to multi-parameter
quantum gates [76].

The gradient descent methods can also be transformed into the quantum version
[77, 78, 79, 80, 81, 82]. The quantum versions of the gradient descent was
firstly proposed by Rebentrost et al. for high-dimensional optimization
problems [83]. The gradient algorithm in Ref. [83] involves phase estimation,
which requires substantial circuit depth and is difficult to implement with
current quantum hardwares. Instead of using phase estimation, Li et al.
proposed an experimental friendly quantum gradient algorithm [77] using the
linear combination of unitaries (LCU).

In addition to above methods, some other gradient descent methods have been
proposed, such as quantum analytic descent [84] and stochastic gradient
descent [85, 86, 75], etc.

Gradient-free Methods. Similar to the classical machine learning, gradient-
free optimization methods can also be used to optimize PQCs. The Nelder-Mead
method is one of commonly used gradient-free optimization for finding a local
minimum of a function of several variables, and thus is naturally suitable for
VQAs as well. Besides, evolutionary algorithms and reinforcement learning can
be also used to train VQAs [87, 88, 89, 90, 91]. These optimization methods
are standard in the classical machine learning community and we will not go
into detail in this review.

###  II.2 Architectures and applications

VQAs can be applied to a wide range of fields, including finding ground states
of molecules, combinatorial optimization, machine learning, etc. Here, we will
introduce some typical applications and their experimental progress.

####  II.2.1 Variational quantum eigensolver (VQE)

The VQE [60, 55] is a flagship algorithm for NISQ applications. This type of
algorithms allows one to find the ground state of a given Hamiltonian $H$,
which may be used for simulating molecules and chemical reactions. According
to the Rayliegh-Ritz variational principle [92, 93, 94], the ground state
energy $E_{0}$ associated with this Hamiltonian $H$ is bounded by

${E_{0}}\leq\frac{{\langle\psi|H|\psi\rangle}}{{\langle\psi|\psi\rangle}},$
(6)

where $|\psi\rangle$ is a trial quantum state. The original VQE is to find a
quantum state approximating the ground state by training a PQC, such that the
expected value of the Hamiltonian is constantly approaching the minimum. In
addition to its original purpose, VQE has been widely extended to other
objectives, such as finding excited states or the spectrum of Hamiltonian [95,
96, 97, 98, 99, 100, 101], nonequilibrium steady state [102, 103, 104], and
calculating energy derivatives [105, 106, 107]. Recently, Wei et al. propose a
full quantum eigensolver (FQE) [108] which treats the gradient descent part in
a full quantum mechanical manner.

VQE has been demonstrated on various quantum architectures such as
superconducting qubits [109, 110, 39, 111, 112, 113, 40, 114, 115], photonic
system [60, 116], and trapped ions [117, 118, 119]. These experiments are
still boil down to proof-of-principle demonstrations on small-molecule
systems, and thus further experimental efforts are required to scale up this
approach to larger molecular systems of chemical interest. Meanwhile, the
algorithmic approaches also needs to be further developed to relax the
hardware limitation, such as reducing the number of qubits [120, 121, 95, 122,
62, 123, 66, 124, 125] and circuit depth [66, 126, 127, 128, 129], error-
mitigation techniques [130, 131, 132], accelerated VQE [133], measurement-
based VQE [134] and so on.

####  II.2.2 Combinatorial optimization

Combinatorial optimization is a promising path to demonstrating practical
quantum advantage on near-term quantum devices, most notably the QAOA [58,
135] for solving the combinatorial optimization problem MaxCut. Max-Cut is the
NP-complete problem, which is defined to partition the nodes of a graph into
two distinct sets $A$ and $B$ that maximizes the number of edges connecting
nodes in opposite sets. It has been shown that to achieve an approximation
ratio of $r\geq 16/17\approx 0.9412$ for Max-Cut on all graphs is NP-Hard
[136].

Mathematically, consider a graph with $m$ edges and $n$ vertices, we can use
bitstring $z={z_{1}}\ldots{z_{n}}$ to represent the assignment of vertices to
the two sets, where ${z_{i}}=0$ if $i$-th vertex is belongs to set $A$ and
${z_{i}}=1$ if it belongs to the other set $B$. The objective is to maximize
the number of edges cut, denoted as $C(z)$. When implementing the QAOA to find
approximate solutions to the MaxCut problem, we denote the bitstring using
computational basis states $|z\rangle$, and define the objective function to
maximize as

$C(z)=\sum\limits_{\alpha=1}^{m}{{C_{\alpha}}(z)}=\sum\limits_{\text{edge}(j,k)}^{m}{\frac{1}{2}(1-\sigma_{z}^{j}\sigma_{z}^{k})}|z\rangle,$
(7)

where $(j,k)$ vertices connect the $\alpha$ edge. ${C_{\alpha}}$ has the value
1 only if the $j$-th and $k$-th qubits have different measurement results on
the $Z$ basis, representing separate partitions. The $p$-layer QAOA
parametered circuit is usually governed by the problem Hamiltonian
$H_{C}=\sum\limits_{\text{edge}(j,k)}^{m}{\frac{1}{2}(1-\sigma_{z}^{j}\sigma_{z}^{k})}$
and a mixing Hamiltonian $H_{M}=\sum\limits_{i=1}^{n}{\sigma_{x}^{i}}$
alternatively in each layer, as

$U({\bm{\gamma}},{\bm{\beta}})=\prod\nolimits_{l=1}^{p}{{e^{-i{\beta_{l}}{H_{M}}}}{e^{-i{\gamma_{l}}{H_{C}}}}},$
(8)

where $\bm{\gamma}=(\gamma_{1},\ldots,\gamma_{p})$ and
$\bm{\beta}=(\beta_{1},\ldots,\beta_{p})$ are variational parameters. The
input state of QAOA circuit is set as a $n$-qubit uniform superposition state
$|+\rangle^{n}$, and the final state
$\varphi({\bm{\gamma}},{\bm{\beta}})=U({\bm{\gamma}},{\bm{\beta}})|+\rangle^{n}$
is measured to output the bitstring $z$. We train the QAOA parametered circuit
to evolve the initial state into the bitstring $z$ with the maximum partition.

There has been a lot of theoretical works discussing the performance of QAOA
and evaluating the resources needed to achieve quantum computational supremacy
[50, 137, 138, 139, 140, 141, 142, 143]. In Ref. [138], Streif et al.
presented a comparison between the QAOA with competing methods, quantum
annealing and simulated annealing. Guerreschi et al. claimed that at least
several hundreds of qubits are required for QAOA to achieve quantum speed-up
[139]. Dalzell et al. concluded that 420 qubits QAOA would be sufficient for
quantum computational supremacy [141]. Furthermore, it was proven that
classical Goemans-Williamson algorithm outperforms the QAOA for certain
instances of MaxCut at any constant level [144]. However, Bravyi et al. [142]
and Egger et al. [143] suggested that higher-level Recursive-QAOA is
competitive (and often better than) to best known generic classical algorithm
based on rounding an SDP relaxation. Overall, whether QAOA has advantages over
classical algorithms, and in which issues, is still a controversial issue that
needs more research. Besides, to improve the performance of the QAOA, great
efforts have been made, such as finding a good classical optimizer [140, 145,
146, 147, 148], ansatz design [67, 149, 150, 151], reducing circuit depth
[152, 153] and number of qubits [154], and so on.

In terms of experimental progress, QAOA has been experimentally implemented in
superconducting [155, 49], trapped-ion [156], and Rydberg atom [48] platform.
In particular, in 2021, Harrigan et al. realized the QAOA with up to 23 qubits
on the Google’s Sycamore superconducting quantum processor [49]. Ebadi et al.
experimentally investigated quantum algorithms for solving the maximum
independent set problem using Rydberg atom arrays with up to 289 qubits in two
spatial dimensions [48].

####  II.2.3 Machine learning

Nowadays, artificial intelligence and machine learning have become an integral
part of modern life. As quantum computing promises to enhance our capacity to
do some crucial computational tasks, it is natural to look for linkages
between these two cutting-edge techniques in the goal of gaining a variety of
advantages. In fact, VQAs can be seen as quantum analogs of classical neural
networks, and thus can easily be used for various machine learning tasks, such
as classification and generative models.

QCNN. Convolutional neural networks (CNNs) stand out from other neural
networks for their outstanding capabilities in computer vision [157, 158, 159,
160], and natural language processing [161, 162, 163, 164, 165, 166]. A CNN
usually consists of three main types of layers: 1) convolutional layer, 2)
pooling layer, and 3) fully-connected layer, where convolutional layers use
filters that perform convolution operations on the input, and pooling layers
are a type of downsampling operation.

Inspired by CNNs, Cong et al. proposed a quantum convolutional neural network
(QCNN), and showed its ability to solve quantum many-body problems [41]. In
the QCNN architecture, the convolutional and pooling layers, as well as the
fully-connected layer, are designed as quantum circuits, although their actual
functions do not correspond exactly to the classical ones. The analysis by
Pesah et al. implied that QCNNs do not exhibit barren plateaus [167], which
are important for trainability as the system scales. The QCNN has been
realized on a 7-qubit superconducting quantum processor to identify symmetry-
protected topological phases of a spin model [168]. QCNN is more suitable for
quantum problems, in contrast, the hybrid quantum-classical convolutional
neural network (QCCNN) [43] is specially proposed for real-world problems. The
central idea of QCCNN is to implement the feature map in the convolutional
layer with a PQC, while other layers remain the classical operation. Such a
design not only introduces a quantum-enhanced feature map, but also preserves
important features of classical CNN, such as nonlinearity and scalability.
Numerical experiments on the Tetris dataset show that QCCNN can accomplish
classification tasks with learning accuracy surpassing that of classical CNN
with the same structure. Besides, Wei et al. proposed a QCNN framework based
on LCU [169, 170] that can transform CNN to QCNN directly [171].

Quantum GAN. Generative adversarial networks (GANs) are at the forefront of
the generative learning and have been widely used for image processing, video
processing, and molecule development [172]. GANs exploit a two-player minimax
game between a generator and a discriminator, where generator aims to output
the generated data to fool discriminator, and meanwhile, discriminator tries
to distinguish the true example from generator.

Recently, theoretical works show that quantum generative models may exhibit an
exponential advantage over classical counterparts [45, 173, 174], arousing
extensive research interest in the theory and experiment of quantum GANs [175,
176, 177, 178, 179, 180, 181, 47, 182, 183]. In particular, Huang et al.
developed a resource-efficient quantum GAN, and for the first time
accomplished the generative task of real-world handwritten digit image on a
superconducting quantum processor [46]. Being able to handle real-world data
generation tasks is an exciting thing for an “infant” quantum computer.
Subsequently, Rudolph et al. provided an experimental implementation of
generating more high-resolution images of handwritten digits with quantum GAN
on an ion-trap quantum computer [47]. Niu et al. showed that by learning a
shallow quantum circuit to generate a superposition of classical data, their
proposed entangling quantum GAN can be used to create an approximate quantum
random access memory (QRAM) [181].

Quantum Kernel Estimation. Quantum Kernel estimation is a common method for
supervised learning classification task on the NISQ device. In Ref. [184],
Schuld et al. interpreted the process of encoding classical information into a
quantum state as a nonlinear feature map which assigns data to the Hilbert
space of the quantum system, and then proposed to use a variational quantum
circuit to process the feature vectors. Havlíček et al. proposed and
experimentally implemented binary classifiers that use the quantum state space
as feature space on a superconducting processor [42], providing a possible
path to quantum advantage.

Quantum Auto-encoder. Quantum auto-encoder (QAE) is an efficient VQA for
quantum data compression. QAE starts from a large-scale quantum state
$\rho_{AB}$, utilizes a PQC $U({\bm{\theta}})$ to encode $\rho_{AB}$ into the
$A$ subsystem (latent space), and then recover the initial state $\rho_{AB}$
from this compressed state with high fidelity via a decoder constructed by
another PQC ${V({\bm{\theta^{\prime}}})}$. Usually $V$ is set as
${V({\bm{\theta^{\prime}}})}={U({\bm{\theta}})^{\dagger}}$. QAE has recently
attracted great attentions [185, 186, 187, 188, 189, 190], and its
experimental implementations have been carried out on linear optical systems
[191, 192] and superconducting system [193].

####  II.2.4 Cryptanalysis

Shor’s integer factorization algorithm is considered to be one of the most
influential algorithms that shape the research interest in quantum computing
today. However, implementing Shor’s algorithm requires a fault-tolerant
quantum computer, which is far from attainable. Variational Quantum Factoring
(VQF) algorithm is heuristic alternative to Shor’s algorithm suitable for NISQ
devices [194]. It works by reducing the factoring problem to the ground state
of the Ising Hamiltonian, which can then be found using VQA. Karamlou et al.
reported a experimental demonstrations of VQF using a superconducting quantum
processor [195], the integers 1099551473989, 3127, and 6557 are factored with
3, 4, and 5 qubits, respectively. Besides, the performance and resource
analysis of VQF are presented in Refs. [196, 197].

In addition to the factoring task, VQAs can also be used to attack advanced
encryption standard (AES)-like symmetric cryptography [198]. In their
simulation results, sometimes the variational quantum attack algorithm is even
faster than Grover’s algorithm. Moreover, Coyle et al. proposed cryptanalysis
algorithm based on variational quantum cloning (VarQlone), which allows an
adversary to obtain optimal approximate cloning strategies for unknown quantum
states with shallow quantum circuits [199].

####  II.2.5 Mathematical applications

Numerical solvers such as solving linear systems, non-linear differential
equations, $etc$, have a very wide range of engineering application scenarios.
Several classical-quantum hybrid algorithms derived from VQAs make it possible
for near-term quantum devices to solve large-scale numerical problems without
the use of fault-tolerant quantum computing.

Variational Quantum Linear Solver. Solving systems of linear equations is a
fundamental computational problem, and the vast majority of numerical problems
in science and engineering boil down to solving systems of linear equations.
Unlike the Harrow-Hassidim-Lloyd (HHL) algorithm [200], the variational
quantum linear solver (VQLS) proposed by Bravo-Prieto et al. can run on near-
term quantum computers [201]. By employing a hybrid quantum-classical
algorithm, the goal of VQLS is to variationally prepare $|x\rangle$ such that
$A|x\rangle\propto{\kern 1.0pt}{\kern 1.0pt}|b\rangle$. The numerical
simulations give the evidence that the run time of VQLS has linear scaling in
$\kappa$, logarithmic scaling in $1/\varepsilon$, and polylogarithmic scaling
in $N$, where $\kappa$ is the condition number of a $N\times N$ matrix $A$,
and $\varepsilon$ is the precision of the solution. In the same period, Huang
et al. [202] and Xu et al. [203] independently proposed similar quantum
algorithms.

Solving nonlinear differential equations. Differential equations are
ubiquitous in various fields of science and engineering. Lubasch et al. first
extended the concept of VQAs to introduce nonlinearity for solving nonlinear
problems [204], such as the nonlinear Schrödinger equation. Subsequently,
several approaches to solve the nonlinear differential equation using the VQAs
were proposed [205, 206, 207, 208], and the application to computational fluid
dynamics was studied [209].

####  II.2.6 Many-body physics

It is more natural to use quantum circuits to process quantum data than to
process classical data, because it avoids the input and output problems of
classical-quantum data conversion. VQAs enable us to realize the simulation of
many-body dynamics [210], as well as identifying many-body phase transition
[51, 211, 168]. In particular, Gong et al. proposed a quantum neuronal sensing
based on digital-analog variational quantum circuit, and exprimentally
realized this scheme to classify the ergodic and localized phases of matter on
a 61 qubit superconducting quantum processor [51]. This experiment
demonstrates the experimental feasibility of large-scale VQAs, and opens new
avenues for exploring quantum many-body phenomena in larger-scale systems.

###  II.3 Opportunities and Challenges

Due to their wide applications and noise resistance to NISQ devices, VQAs have
a huge potential to achieve quantum advantages for practically interesting
problems. Over the past few years, VQAs have received unprecedented attention,
and papers related to VQAs appear on the arXiv website almost every day. The
following points deserve our attention in advancing VQAs into future
industrial applications or truly meaningful scientific research.

Where is the main advantage? Although VQAs have been hotly researched for
several years, this core question has never been well answered. Perhaps using
VQAs to handle quantum problems is a natural advantage [51], since either
simulating quantum systems or measuring quantum systems using classical
methods would be exponentially resource-intensive. However, for real-world
problems, it becomes very tricky to answer why we need VQAs. Does it offer
speed or accuracy advantages over classical neural networks? Maybe we can try
to answer this question in two ways: 1) Find the provable advantages of VQAs
theoretically. This should be difficult, since today’s machine learning
algorithms are still theoretically difficult to study; 2) Experimentally
achieve milestones that beat current state-of-the-art classical machine
learning in terms of speed or accuracy on large-scale practical datasets
(rather than toy models). Either way, it will take a huge effort.

I/O problems for real-world data. To fully exploit the superposition
properties of quantum systems, we actually want the classical data to be
encoded into the amplitude of quantum system. However, this is extremely
resource-intensive, whether using QRAM [212] or pure quantum circuits [213,
214]. Thus, for real-world applications, perhaps we need to find some specfic
classical problems suitable for quantum computing. For example, the input to
this classical problem has a certain sparsity or symmetry, which is conducive
to being encoded as a quantum state. For the output, in general we need to
avoid the exponential consumption caused by a large number of measurements,
which requires us to extract and analyze only part of the information of the
output quantum state, such as the expected value of observables, principal
component analysis [215], $etc$.

Trainability and Training Efficiency. Avoiding barren plateaus is now an
important area of research, as barren plateaus in training landscapes will
appear if care is not taken enough. Some results suggest that the exponential
parameter space, the noise and decoherence, and even entanglement can induce
barren plateaus [216, 217, 218]. At present, some methods have discussed how
to mitigate barren plateaus [219, 220, 221, 222, 56], and more research is
still needed.

One of the weaknesses of quantum computers compared to classical computers is
quantum noise, which also greatly affects the training performance. Generally,
the resillience of VQAs exists for a wide class of noise, as analyzed in Refs.
[55, 223, 224]. Yet such noise resillience is highly limited as the noise
level rising. As suggested in Refs. [225, 226, 227, 217, 228], large noise may
leads to problems such as performance degradation or barren plateaus. Some
special noises, such as leakage error, also have a bad effect on the
performance of VQAs [229]. Thus, efficient error mitigation schemes for VQAs
need to be developed [230, 231, 232, 233].

Besides, the training process of VQA is very time-consuming, due to the
incompatibility with backpropagation and the cost of a large number of
measurements, posing a great challenge to the large-scale development of VQAs.
Some approaches try to use multiple QPUs for parallel training to alleviate
this deficiency, such as data parallelism [234] or parameter parallelism
[235]. However, more fundamental solutions are urgently needed.

Hardware development. For practical NISQ applications, high-quality quantum
processors with more qubits, longer coherence times and lower error rates are
prerequisites. For VQAs, iterative training requires frequent interaction
between classical and quantum computers, which imposes some new demands on the
quantum hardware. The cryogenic electronic control architecture is a hardware-
level solution for accelerating classical-quantum interactions [236, 237, 238,
239, 240], and also an advanced technique required for the future development
of quantum computing.

##  III quantum error mitigation

Figure 3: Summary of popular Quantum Error Mitigation(QEM) methods.
Probabilistic error cancellation estimates the noiseless expectation value
using a linear combination of expectation values with different noise terms,
which is based on the quasi-decomposition. The data-driven approaches
including Zero-noise extrapolation and learning-based methods collect
expectation values of circuits with different error rates, or circuits with
similar structure. Purification-based error mitigation schemes use multiple
($M$) copies of the noisy state and estimate
$\frac{\text{Tr}(\rho^{M}O)}{\text{Tr}(\rho^{M})}$ using collective
measurements. Verification methods are designed for systems with certain
symmetries and combined with classical post-processing. Quantum subspace
expansion also uses post-processing to mitigate errors for some VQE
algorithms. Other mitigation methods for certain noise types, such as
measurement error mitigation (for measurement errors) are also widely applied.

In the context of current NISQ devices, the impact of noise remains the
greatest challenge for the practical applications [6]. Although quantum error
correction (QEC) promises to enable quantum computation with arbitrary levels
of noise, it is out of reach for near-term quantum processors. Quantum error
mitigation (QEM) provides us a feasible alternative to mitigate errors of
near-term quantum processors, and is also the continuous path that will take
us from today’s quantum hardware to tomorrow’s fault-tolerant quantum
computers. Instead of active error-correction, QEM methods usually estimate
the error-free expectation value by classical post-processing of the noisy
measurement results. While it introduces the additional sampling overhead
caused by the increase in the variance of the mitigated observable, QEM
requires fewer qubits and gate resources and is therefore more suitable for
practical NISQ devices. In recent years, a wide variety of QEM methods have
been proposed [131, 28, 130, 241, 242, 243, 244], and we will discuss some
typical QEM techniques in the following subsections.

###  III.1 Probabilistic error cancellation

We first introduce probabilistic error cancellation (PEC) [245, 130]. The key
idea of PEC is to apply the quasi-probability decomposition of the inverse
noise process, leading to a linear combination of a set of noisy circuits. The
implementation of PEC requires the full knowledge of the noise model in the
target circuit and is usually combined with characterization of the noise
process. The PEC method was first introduced by Temme et al. [245] and further
detailed for a practical scheme by Endo et al. [130]. It has been shown to
mitigate Markovian noise [130] and has been further generalized to the non-
Markovian noise [246]. And it was experimentally demonstrated on
superconducting [247] and trapped ion [248] quantum computers.

####  III.1.1 Standard PEC with discrete gate-based circuits

Quasi-probability decomposition of inverse noise model. Suppose we have a
noisy quantum gate $\mathcal{E}\circ\mathcal{U}$ contaminated by a noise
channel $\mathcal{E}$, where $\mathcal{U}(\rho)=U\rho U^{\dagger}$ denotes the
noiseless gate with the ideal operator $U$ and the initial quantum state
$\rho$. With approaches to accurately characterize the noise channel
$\mathcal{E}$, we can have the mathematical form of the inverse noise channel
$\mathcal{E}^{-1}$, and then apply $\mathcal{E}^{-1}$ after the noisy gate
$\mathcal{E}\circ\mathcal{U}$ to recover the ideal gate
$\mathcal{U}=\mathcal{E}^{-1}\circ\mathcal{E}\circ\mathcal{U}$. Since the
inverse noise channel $\mathcal{E}^{-1}$ may be unphysical, we need a set of
basis operations $\\{\mathcal{B}_{k}\\}_{k}$ to decompose it as

$\mathcal{E}^{-1}=\sum_{k}q_{k}\mathcal{B}_{k},$ (9)

where $q_{k}$ is the the combination coefficients, i.e., quasi-probabilities.
So we can estimate the noiseless expectation value of the observable $O$ using
a linear combination of results obtained from circuits applied by different
basis operations

$\displaystyle\langle O\rangle_{\text{PEC}}=$
$\displaystyle\sum_{k}q_{k}\text{Tr}[\mathcal{B}_{k}\circ\mathcal{E}\circ\mathcal{U}(\rho)]$
(10) $\displaystyle=$ $\displaystyle
Q_{\mathcal{E}}\sum_{k}\text{sgn}(q_{k})\frac{|q_{k}|}{Q_{\mathcal{E}}}\text{Tr}[\mathcal{B}_{k}\circ\mathcal{E}\circ\mathcal{U}(\rho)],$

where $Q_{\mathcal{E}}=\sum_{k}|q_{k}|$. Given above quasi-probability
decomposition, we can use the Monte Carlo sampling to randomly apply the basis
operation set $\\{\mathcal{B}_{k}\\}_{k}$ with the quasi-probability
distribution $\\{\frac{|q_{k}|}{Q_{\mathcal{E}}}\\}$. Then we multiply the
measurement outcome by the corresponding sgn($q_{k}$) and weight the Monte
Carlo average by $Q_{\mathcal{E}}$. For the general circuit composed of a
$N$-gate sequences $\mathcal{U}=\prod_{n=1}^{N}\mathcal{U}_{n}$, we can just
find the quasi-probability decomposition for noise process $\mathcal{E}_{n}$
associated with each gate

$\mathcal{E}_{n}^{-1}=\sum_{k_{n}}q_{k_{n}}\mathcal{B}_{k_{n}}=Q_{n}\sum_{k_{n}}\text{sgn}(q_{k_{n}})\frac{|q_{k_{n}}|}{Q_{n}}\mathcal{B}_{k_{n}},$
(11)

in which $Q_{n}=\sum_{k}|q_{n_{k}}|$.

Then we can obtain the error mitigated value in a similar way

$\displaystyle\langle O\rangle_{\text{PEC}}=$
$\displaystyle\text{Tr}\left[\prod_{n}^{N}\left(\sum_{k_{n}}q_{k_{n}}\mathcal{B}_{k_{n}}\circ\mathcal{E}_{n}\circ\mathcal{U}_{n}\right)(\rho)\right]$
(12) $\displaystyle=$ $\displaystyle
Q_{\mathcal{E}}\sum_{\vec{k}}\text{sgn}(q_{\vec{k}})\frac{|q_{\vec{k}}|}{Q_{\mathcal{E}}}\text{Tr}\left[\prod_{n=1}^{N}\left(\mathcal{B}_{k_{n}}\circ\mathcal{E}_{n}\circ\mathcal{U}\right)(\rho)\right],$

with $\vec{k}=(k_{1},k_{2},...,k_{N})$ for
$q_{\vec{k}}=\prod_{n=1}^{N}q_{k_{n}}$ and
$Q_{\mathcal{E}}=\prod_{n=1}^{N}Q_{n}$.

Sampling overhead. Note that the total resource overhead for PEC involves both
the characterization cost of the noise channel, and the sampling overhead to
estimate the linear combinations in Eq 10. Considering that the gate noise
characterization is usually performed during the device calibration stage
[249], we here focus on the sampling overhead for the implementation of PEC.

As for the factor $Q_{\mathcal{E}}$ multiplied with the averaged expectation
value, the variance is increased by
$Q_{\mathcal{E}}^{2}=(\prod_{n=1}^{N}Q_{n})^{2}$. Therefore, it needs
$Q_{\mathcal{E}}^{2}$ times more samples to achieve the same measurement
accuracy as the unmitigated case. The square of the multiplier
$Q_{\mathcal{E}}^{2}$ can be regarded as the sampling overhead
$\Upsilon_{\text{PEC}}$ for PEC, which grows exponentially with the overall
gate error $\varepsilon N$ if we assume $Q_{n}=1+\mathcal{O}(\varepsilon)$ for
stochastic noise process [130]. This sets limitations on the efficiency of PEC
schemes to $\varepsilon N=\mathcal{O}(1)$. Moreover, the theoretical analysis
of the sample complexity for implementing PEC has been explored [250, 251] and
it has been shown the optimality of PEC among all strategies in mitigating a
certain type of noise [251].

Approaches to extract error parameters. To implement the PEC in practice, we
need to accurately characterize the noise model $\mathcal{E}$ and obtain the
analytical decomposition of the inverse noise map as Eq. 9 shown.

Gate set tomography. Endo et al. [130] has introduced gate set tomography
(GST) [252, 253] to PEC, which is free from the SPAM errors. Although GST can
reconstruct arbitrary noise process, the number of samples required for its
reconstruction grows exponentially with the size of the system $N_{q}$ and is
quite costly. If only specific types of noise are considered, there are some
more efficient characterization protocols [254, 255, 256].

Cycle error reconstruction for Pauli noise. Cycle error reconstruction (CER)
is an efficient approach to identify the Pauli noise in circuits [257, 258].
It utilizes Cycle Benchmarking [259] to characterize the Pauli noise channel
and extract the Pauli error rates after post-processing. The implementation of
PEC with CER has been shown on a 4-qubit superconducting processor [260].
Under the Pauli noise assumption, the quasi-probability decomposition can be
further simplified using sparse Pauli-Lindblad models [261]. The sparse Pauli
noise model considers Pauli operators with only one or two non-trivial terms,
but is sufficient to capture the correlated errors. Moreover, the number of
noise parameters in the sparse model scales polynomially with the system size,
so it remains efficient in large-scale quantum systems. This approach has been
experimentally demonstrated to learn the noise model on a superconducting
quantum processor of up to 20 qubits [261].

####  III.1.2 PEC with continuous time evolutions

The standard PEC scheme we discussed above is based on the discrete gate-based
quantum circuits with the gate-independent Markovian noise. To overcome this
limitation, the stochastic QEM method [262] is proposed to extend the standard
PEC to more general scenarios that may have strong gate-dependence and
complicated nonlocal effects, and general computing models such as analog
quantum simulators.

Specifically, the time evolution of a noisy quantum system can be described by
the Lindblad master equation as

$\displaystyle\frac{\partial\rho(t)}{\partial
t}=-i[H^{\text{noisy}}(t),\rho(t)]+\lambda\mathcal{L}[\rho(t)],$ (13)
$\displaystyle\mathcal{L}(\rho)=\sum_{i}(2L_{i}\rho
L_{i}^{\dagger}-L_{i}^{\dagger}L_{i}\rho-\rho L_{i}^{\dagger}L_{i})$

where $H^{\text{noisy}}=H+\delta H_{0}$ ($H$ is the ideal evolution
Hamiltonian, $H_{0}$ denotes the coherent errors), and $\mathcal{L}(\rho)$ is
the noise Lindblad operator corresponding to the decoherent coupling with the
environment with noise strength $\lambda$. We denote the ideal and noisy
process by $\mathcal{E}_{I}$ and $\mathcal{E}_{N}$ respectively. Given a small
time step $\delta t$ with a time interval $\delta$, both the ideal and noisy
case of the evolution can be written as

$\rho_{i}(t+\delta)=\mathcal{E}_{i}(t)\rho_{i}(t),$ (14)

with $i=I,N$. Similar to the idea of PEC in Eq. 9, we can find a recovery
process $\mathcal{E}_{Q}$ satisfying

$\displaystyle\mathcal{E}_{I}=\mathcal{E}_{Q}\mathcal{E}_{N}+\mathcal{O}(\delta
t^{2})\approx\mathcal{E}_{Q}\mathcal{E}_{N},$ (15)
$\displaystyle\mathcal{E}_{Q}=\sum_{k}q_{k}\mathcal{B}_{k}=Q\sum_{k}\text{sgn}(q_{i})\frac{|q_{i}|}{Q},$

with $Q=\sum_{k}|q_{k}|$. Given the full evolution time $T$, we can apply the
Monte Carlo sampling to realize the continuous recovery $\mathcal{E}_{Q}$ with
a time interval $\delta$ at the mitigation cost of
$\Upsilon_{\text{SEM}}=Q^{\frac{T}{\delta t}}$ [262, 131].

###  III.2 Data-driven approaches

We now introduce the data-driven approaches for error mitigation, including
zero noise extrapolation (ZNE) and error mitigation with learning process such
as Clifford data regression (CDR). These methods utilize different types of
data from circuits with various error rates (e.g., ZNE), or from near-Clifford
circuits under similar circuit structure (e.g., CDR), and they can be
naturally combined with each other to better utilize data resources [263,
264].

####  III.2.1 Zero-Noise Extrapolation

ZNE is a typical error mitigation scheme using classical post-processing. It
uses data collected at different error rates to fit the function of
expectation values with respect to the error rates, and then extrapolate to
the zero noise limit [244, 130, 245, 243, 263, 265].

Fitting methods of ZNE. The chosen fitting method is essential to the
performance of ZNE. A basic fitting method is the Richardson extrapolation,
which utilizes the polynomial relationship between noisy and ideal expectation
values yielded by Tailor expansion when the error rate is low. This estimation
is slightly coarse. More refined models, such as exponential and poly-
exponential fitting functions, can be introduced in the extrapolation method
to improve the actual mitigation performance in specific cases. We discuss
these fitting methods in detail.

Richardson extrapolation. Richardson extrapolation was the first introduced
extrapolation method [244, 245] and later experimentally demonstrated using
superconducting qubits [110].

Suppose the quantum circuit outputs a noisy $N_{q}$-qubit quantum state
$\rho_{\varepsilon}$, where $\varepsilon$ denotes the noise parameter which
characterizes the error rate in the circuit. Then the expectation value of the
target observable $O$ under error rate $\varepsilon$ can be regarded as a
function towards $\varepsilon$: $\langle
O\rangle(\varepsilon)=\text{Tr}(\rho_{\varepsilon}O)$. In order to evaluate
the error-free expectation value $\langle O\rangle(0)$, we choose a set of
error rates $\\{\lambda_{0}\varepsilon,...,\lambda_{n}\varepsilon\\}$ with
$n+1$ different coefficients $\\{\lambda_{i}\\}_{i=0,...,n}$($\lambda_{0}=1$)
, and then run circuits under these error rates to obtain the set of $n+1$
expectation values $\\{\langle O\rangle(\lambda_{i}\varepsilon)\\}_{i}$. It is
known to be difficult to reduce the error rate in the circuit, we thus usually
boost it with amplified coefficients
$1=\lambda_{0}<\lambda_{1}<...<\lambda_{n}$ using some noise scaling methods.
Linear spacing of the amplified coefficients (e.g.
$\\{\lambda_{i}=i+1\\}_{i=0,...,n}$) is commonly used [245, 249] and some
specific spacing methods have been explored for better performance [266].

The expectation value for $\lambda=0$ (the error-free case) can be estimated
using

$\langle O\rangle_{\text{Rid}}=\sum_{i=0}^{n}\gamma_{i}\langle
O\rangle(\lambda_{i}\varepsilon),$ (16)

where the fitting coefficients $\\{\gamma_{i}\\}$ are chosen to satisfy both
$\sum_{i=0}^{n}\gamma_{i}=1$ and $\sum_{i=0}^{n}\gamma_{i}\lambda_{i}^{j}=0\
\text{for}\ j=1,...,n$, and the solution gives

$\gamma_{i}=\prod_{i\neq j}\frac{\lambda_{j}}{\lambda_{j}-\lambda_{i}}.$ (17)

If the error rate $\varepsilon$ is sufficiently low, we can expand the
function $\langle O\rangle(\varepsilon)$ according to the Taylor expansion

$\langle O\rangle(\varepsilon)=\langle
O\rangle(0)+\sum_{j=1}^{n}a_{j}\varepsilon^{j}+\mathcal{O}(\varepsilon^{n+1})$
(18)

with some constant coefficients $a_{j}$. Then the Richardson mitigated value
$\langle O\rangle_{\text{Rid}}$ can be rewritten using Eq. 16 and 18

$\displaystyle\langle O\rangle_{\text{Rid}}$
$\displaystyle=\sum_{i=0}^{n}\gamma_{i}\left(\langle
O\rangle(0)+\sum_{j=1}^{n}a_{j}\lambda^{j}\varepsilon^{j}+\mathcal{O}(\varepsilon^{n+1})\right)$
$\displaystyle=\langle
O\rangle(0)+\sum_{j=1}^{n}\left(\sum_{i=0}^{n}\gamma_{i}\lambda_{i}^{j}\right)\varepsilon^{j}+\mathcal{O}(\varepsilon^{n+1})$
$\displaystyle=\langle O\rangle(0)+\mathcal{O}(\varepsilon^{n+1}),$

Here we use $\sum_{i=0}^{n}\gamma_{i}=1$ and
$\sum_{i=0}^{n}\gamma_{i}\lambda_{i}^{j}=0\ \text{for}\ j=1,...,n$. We can see
that the error rate is suppressed from the original $\varepsilon$ to the order
$\mathcal{O}(\varepsilon^{n+1})$, under the weak noise assumption. However,
the assumption of a valid Taylor expansion for low error rate may be
inaccurate in the large circuit limit, and the lack of specific information on
the noise channel sets limitation to the efficacy of Richardson extrapolation.

As for the sampling overhead of Richardson extrapolation, we need to take into
account the variance increase introduced in Eq. 16, which reads [131]

$\text{Var}(\langle
O\rangle_{\text{Rid}})=\sum_{i=0}^{n}\gamma_{i}^{2}\text{Var}(\langle
O\rangle_{\lambda_{i}\varepsilon}).$ (19)

The variance of estimating $\langle O\rangle_{\text{Rid}}$ is about
$\sum_{i=0}^{n}\gamma_{i}^{2}$ larger than evaluating $\langle O\rangle$
without error mitigation, so it requires around
$\Upsilon_{\text{Rid}}=\sum_{i=0}^{n}\gamma_{i}^{2}$ more samples as the
sampling overhead of Richardson extrapolation.

Besides, an optimized protocol [266] for the implementation of Richardson
extrapolation has been proposed to further explore the relevant parameters of
Richardson extrapolation.

Exponential Extrapolation. The error rate $\varepsilon$ used in Richardson
extrapolation quantifies the local noise strength, but in the context of NISQ
error mitigation, it is more natural to consider the mean circuit error count
[243] $\mu$, where $\mu\approx N\varepsilon$ and $N$ denotes the number of
gates in the circuit. In the large circuit limit when $N\gg 1$, the number of
errors happened in the circuit (denoted by $k$) follows the Poisson
distribution with probability [243]

$p_{k}=e^{-\mu}\frac{\mu^{k}}{k!}.$ (20)

Denoting the expectation value with $k$ errors occurred as $\langle
O_{k}\rangle$, the expectation value with mean circuit error count $\mu$ is

$\langle O_{\mu}\rangle=\sum_{k=0}^{\infty}p_{k}\langle
O_{k}\rangle=e^{-\mu}\sum_{k=0}^{\infty}\frac{\mu^{k}}{k!}\langle
O_{k}\rangle,$ (21)

where the factor $e^{-\mu}$ implies an exponential decay with $\mu$ of the
expectation value. By simply assuming an exponential function

$\langle O\rangle(\mu)=e^{-f\mu}\langle O\rangle(0)$ (22)

with a parameter $f$ which denotes the observable decay rate, we can then
obtain an two-points exponential extrapolation result with mean circuit error
count $\mu$ and $\lambda\mu$ ($\lambda>1$)

$\langle O\rangle_{exp}={\left(\frac{\langle O\rangle(\lambda\mu)}{\langle
O\rangle^{\lambda}(\mu)}\right)}^{\frac{1}{\lambda-1}},$ (23)

Endo et al. [130] first considered the exponential decay curve for error
extrapolation, where the advantage of exponential extrapolation over
Richardson extrapolation has been numerically demonstrated [130, 265]. And Cai
[243] provided a general multi-exponential extrapolation framework

$\langle O\rangle(\mu)=\sum_{i=1}^{d}A_{i}(e^{-f_{i}})^{\mu},\text{and}\
\sum_{i=1}^{d}A_{i}=\langle O\rangle(0),$ (24)

which achieves a much lower estimation bias in numerical simulation.

Poly-exponential extrapolation. The poly-exponential extrapolation [265, 267]
assumes a more general model that the exponential decay with the mean circuit
error count $\mu$ has a polynomial expansion which can be fitted to the
function

$\langle O\rangle(\mu)=e^{\sum_{i=1}^{d}f_{i}\lambda^{i}}\langle O\rangle(0).$
(25)

Then the single-exponential extrapolation modeled by Eq. 22 is a particular
case of the more general poly-exponential extrapolation.

Noise scaling methods. The implementation of ZNE, especially Richardson
extrapolation, depends non-trivially on the ability to boost the error rate
$\varepsilon$ in the circuit by different amplified
coefficients$\\{\lambda_{i}\\}_{i}$ [266]. Next we will introduce a variety of
methods to amplify the error rate in a controlled manner.

Identity insertion [130, 268] is a hardware-agnostic approach which replaces a
particular gate or layer $\mathcal{C}$ by
$\mathcal{C}\\{\mathcal{C}^{\dagger}\mathcal{C}\\}^{n}$ for non-negative
integer $n$. The circuit after insertion is logically equivalent to the
original circuit but the error rate is amplified by a factor $2n+1$ as the
circuit depth increases. There are many variant identity insertion methods
proposed [269, 270] for ZNE, which explore the trade-off between the inserted
gate number and the required measurement to achieve the same (or higher)
accuracy. And more general approaches using unitary folding can obtain an
arbitrary real scaling factor [265]. The limitation of identify insertion
method is that the amplified error rate may arbitrarily deviate from the
desired one [271]. Given a self-adjoint noisy gate
$\tilde{\mathcal{U}}=\mathcal{E}\circ\mathcal{U}$, after identity insertion of
$\tilde{\mathcal{U}}^{\dagger}\tilde{\mathcal{U}}=\tilde{\mathcal{U}}\tilde{\mathcal{U}}$
we have
$\tilde{\mathcal{U}}\tilde{\mathcal{U}}^{\dagger}\tilde{\mathcal{U}}=\mathcal{E}(\mathcal{U}\mathcal{E}\mathcal{U})\mathcal{E}\circ\mathcal{U}=\tilde{\mathcal{E}}\circ\mathcal{U}$.
The actually amplified noise
$\tilde{\mathcal{E}}=\mathcal{E}(\mathcal{U}\mathcal{E}\mathcal{U})\mathcal{E}$
is gate-dependent for general noise channel $\mathcal{E}$, leading to an
unpredictable amplified error rate in the circuit.

Stretching the control pulses in time [245, 110, 271] is also used to
implement circuits under different error rate. Given an open quantum system
described by the time-dependent multi-qubit Hamiltonian $H(t)$, the time
evolution is

$\frac{\partial}{\partial
t}\rho_{\varepsilon}(t)=-i[H(t),\rho_{\varepsilon}(t)]+\varepsilon\mathcal{L}(\rho_{\varepsilon}(t)),$
(26)

where $\varepsilon$ refers to the error rate in the system, and the noise term
$\mathcal{L}$ is a Lindblad operator. Then we rescale the Hamiltonian $H(t)$
to $\frac{1}{\lambda}H(\frac{t}{\lambda})$ with a stretched factor $\lambda>1$
and obtain a new state $\rho_{\varepsilon}^{\prime}(t)$. Under the assumption
that the generator $\mathcal{L}$ is invariant after time rescaling, and also
independent from the way the Hamiltonian $H(t)$ is encoded, we can increase
the evolution time by a factor $\lambda$. Then the error rate is boosted from
$\varepsilon$ to $\lambda\varepsilon$ according to
$\rho_{\varepsilon}^{\prime}(ct)=\rho_{\lambda\varepsilon}(t)$. This method
has been experimentally demonstrated on superconducting processors [110, 271]
using up to 26 qubits [271].

Note that two methods above do not require the exact form of the noise
channel, and there are some approaches that utilize the specific structure of
certain noise models. Parameter noise scaling [265] is designed for the pulse
error caused by imperfect control or finite precision of the physical
parameters in the parametric quantum gates. Considering a quantum gate
$G(\vec{\theta})=\text{exp}(-i\sum_{j=1}^{k}\theta_{j}H_{j})$ which is
parameterized by $k$ classical control parameters
$\vec{\theta}=(\theta_{1},\theta_{2},...,\theta_{k})$ with Hamiltonian
operators $H_{1},H_{2},...,H_{k}$, the actually implemented gate suffered from
pulse error is modeled by $G(\vec{\theta^{{}^{\prime}}})$ with

$\theta^{{}^{\prime}}_{j}=\theta_{j}+\hat{\varepsilon}_{j},\ \text{for}\
j=1,2,...,k,$ (27)

where $\hat{\varepsilon}$ is a random variable following Gaussian distribution
with zero mean and variance $\sigma_{j}^{2}$ which denotes the effect of the
stochastic calibration noise. Given the value of variance $\sigma_{j}^{2}$, we
can rescale the noise by a factor $\lambda>1$ using

$\theta^{re}_{j}=\theta_{j}+\hat{\delta}_{j},$ (28)

here the $\hat{\delta}_{j}$ is sampled from a zero-mean Gaussian distribution
with variance $(\lambda-1)\sigma_{j}^{2}$. Therefore, the effective parameter
after rescaling is

$\theta^{{}^{\prime}}_{j}=\theta_{j}+\sqrt{\lambda}\hat{\varepsilon}_{j},\
\text{for}\ j=1,2,...,k,$ (29)

which achieves a noise rescaling by a factor $\lambda$ without knowledge of
the Hamiltonian operators $H_{j}$.

Pauli twirling techniques [272, 273] focus on the cases that error rates of
single-qubit gates are much lower than those of two-qubit gates and
measurement. By applying randomly chosen Pauli gates before and after the
Clifford two-qubit gates, arbitrary noise channel can be converted into an
effective stochastic Pauli channel [257, 258]. Then we can use additional
Pauli gates to tune the practical error rate to any target value [244], but it
requires the full knowledge of the exact Pauli noise channel. We can also
combine efficient Pauli error reconstruction methods such as Cycle Error
Reconstruction [258, 259] to amplify the error rate for ZNE [260].

Gate Trotterization [274] is a local noise scaling technique acting at the
level of individual gates. We can replace each gate $U$ of the circuit with
the product of $\lambda$ equal gates using the gate Trotterization technique

$U\to(U^{\frac{1}{\lambda}})^{\lambda},\ \lambda=0,1,2...\ .$ (30)

However, the way the $U^{\frac{1}{\lambda}}$ is compiled by the hardware
depends on $\lambda$, so the circuit depth may not increase by the expected
factor.

Leveraging other error mitigation schemes as a noise scaling tool has also
been explored for ZNE [243, 275]. Different from above tuning approaches which
produces boosted error rate, we can also utilize other error mitigation
methods such as PEC to reduce the error rate [243].

####  III.2.2 Multi-dimensional variant of ZNE

Least square fitting [276] may be regarded as a multi-dimensional version but
not reliant on the Richardson extrapolation. It utilizes a hypersurface fit
where one axis refers to the measurement results and the other axes describe
the effect of difference noise parameters. The least square fitting method has
been experimentally demonstrated on Rigetti’s 8-qubit quantum processor [276].

Considering only one noise parameter $\varepsilon$ with $k$ different noise
level $\varepsilon_{1},\varepsilon_{2},...,\varepsilon_{k}$, the polynomial
model up to $n$ order in Eq. 18 can be generalized to a linear equation

$\left(\begin{matrix}1&\varepsilon_{1}&\cdots&\varepsilon_{1}^{n}\\\
1&\varepsilon_{2}&\cdots&\varepsilon_{2}^{n}\\\ \vdots&\vdots&\ddots&\vdots\\\
1&\varepsilon_{k}&\cdots&\varepsilon_{k}^{n}\\\
\end{matrix}\right)\left(\begin{matrix}\langle O\rangle(0)\\\ \alpha_{1}\\\
\vdots\\\ \alpha_{n}\end{matrix}\right)=\left(\begin{matrix}\langle
O\rangle(\varepsilon_{1})\\\ \langle O\rangle(\varepsilon_{2})\\\ \vdots\\\
\langle O\rangle(\varepsilon_{k})\end{matrix}\right).$ (31)

We can solve this linear equation using standard least-squares methods and
obtain the expansion parameters $\left(\begin{array}[]{cccc}\langle
O\rangle(0),&\alpha_{1},&\cdots,&\alpha_{n}\end{array}\right)^{T}$, where
$\langle O\rangle(0)$ is the error-free expectation value of the observable.

Then we can extend the consideration of one noise parameter to the multiple
noise parameters. Take the spontaneous emission rate $\gamma$ and the pure
dephasing rate $\zeta$ for example, suppose we have k different noise level
for both $\gamma$ and $\zeta$ which are denoted by
$\\{\gamma_{1},\gamma_{2},...,\gamma_{k}\\}$ and
$\\{\zeta_{1},\zeta_{2},...,\zeta_{k}\\}$, respectively. The linear equation
in Eq. 31 turns to

$\left(\begin{matrix}1&\gamma_{1}&\zeta_{1}&\gamma_{1}\zeta_{1}&(\gamma_{1})^{2}&(\zeta_{1})^{2}\\\
1&\gamma_{2}&\zeta_{2}&\gamma_{2}\zeta_{2}&(\gamma_{2})^{2}&(\zeta_{2})^{2}\\\
\vdots&\vdots&\ddots&\vdots\\\
1&\gamma_{k}&\zeta_{k}&\gamma_{k}\zeta_{k}&(\gamma_{k})^{2}&(\zeta_{k})^{2}\\\
\end{matrix}\right)\left(\begin{matrix}\langle O\rangle(0)\\\ \alpha_{1}\\\
\vdots\\\ \alpha_{5}\end{matrix}\right)=\left(\begin{matrix}\langle
O\rangle(\gamma_{1},\zeta_{1})\\\ \langle O\rangle(\gamma_{2},\zeta_{2})\\\
\vdots\\\ \langle O\rangle(\gamma_{k},\zeta_{k})\end{matrix}\right).$ (32)

Here we truncate the series to the second order of $\gamma$ and $\zeta$.

####  III.2.3 Error mitigation with learning process

Instead of using the exact knowledge of the noise channel (such as PEC) or the
expectation values at different error rates (like ZNE), the learning-based
error mitigation optimizes the estimator of the observable’s expectation value
using automatical learning process [277, 278, 279, 280, 281, 282, 283].
Specifically, the learning-based error mitigation methods estimate the
observable’s expectation value of the targeted circuit using an ansatz that
typically describes the relationship between the noisy and noiseless case of
the expectation value. The training set fed into the ansatz comes from
classically simulable quantum circuits, or measurement results from relevant
quantum circuits.

Clifford data regression [277] generates training data set
$\\{O_{i}^{\text{noisy}},O_{i}^{\text{ideal}}\\}_{i}$ composed of the noisy
expectation value $O_{i}^{\text{noisy}}$ from real quantum devices and the
noiseless expectation value $O_{i}^{\text{ideal}}$ simulated on the classical
computers. To enable efficient classical simulation, the quantum circuits for
generating training set are replaced with near-Clifford circuits generated by
Markov Chain Monte Carlo sampling method. We can fit the training data set to
a linear ansatz

$O^{\text{ideal}}=f(O^{\text{noisy}},\vec{\theta})=\theta_{1}O^{\text{noisy}}+\theta_{2},$
(33)

where the parameters $\vec{\theta}=(\theta_{1},\theta_{2})$ are optimized by
minimizing the cost function

$C=\sum_{i}\left(O_{i}^{\text{ideal}}-(\theta_{1}O_{i}^{\text{noisy}}+\theta_{2})\right)^{2}.$
(34)

CDR has been experimentally demonstrated on a 16-qubit IBMQ quantum computer
and can achieve an order-of-magnitude improvement for a ground state energy
problem [277], and it also been shown to has the potential to outperform other
state-of-the-art approaches [264]. CDR can be further improved by more
efficient training set construction methods and by applying symmetries in
certain quantum systems [284].

Learning-based probabilistic error cancellation [278] does not depend on the
full tomography of the noise channel as the original PEC does. For a $L$-layer
quantum circuit consisting of single-qubit unitary gates
$\mathbf{R}=(R_{1},R_{2}...,R_{N_{q}(L+1)})$ and multi-qubits Clifford gate
layers $\mathbf{U}=(U_{1},U_{2},...,U_{L})$, the learning-based PEC applies
the Pauli gates $\mathbf{P}=(P_{1},P_{2},...,P_{2N_{q}(L+1)})$ as the error
mitigating gates before and after each single-qubit unitary gates. It can then
establish a circuit configuration where the multi-qubit Clifford gates are
fixed, the other gates ($\mathbf{R}$ and $\mathbf{P}$) are served as
variables, and $\mathbf{P}=\mathbf{I}$ indicates that all error mitigating
gates are identity gates. Here we denote the noisy and error-free expectation
value by $f(\mathbf{R},\mathbf{P})$ and
$f^{\text{ef}}(\mathbf{R},\mathbf{P})$. According to the PEC, the error-
mitigated expectation value $f^{\text{em}}(\mathbf{R},\mathbf{I})$ is

$f^{\text{em}}(\mathbf{R},\mathbf{I})=\sum_{\mathbf{R}}q(\mathbf{P})f(\mathbf{R},\mathbf{P}),$
(35)

where $q(\mathbf{P})$ is the the combination coefficients, i.e.,
quasiprobabilities. By minimizing the cost function,

$C(\mathbf{R})=\frac{1}{|\mathbb{S}|}\sum_{\mathbf{R}\in\mathbb{S}}|f^{\text{em}}(\mathbf{R},\mathbf{I})-f^{\text{ef}}(\mathbf{R},\mathbf{I})|^{2},$
(36)

where a subset of single-qubit Clifford gates $\mathbb{S}$ is chosen as the
training set so that the error-free expectation value
$f^{\text{ef}}(\mathbf{R},\mathbf{I})$ can be efficiently simulated, we can
obtain a optimal distribution $q(\mathbf{P})$ for the training set
$\mathbb{S}$ which can also be applied for the non-Clifford single-qubit
gates. This scheme requires that all the single-qubit gates in the
configuration are error-free and do not rely on the exact error model of the
multi-qubit Clifford gates. For practical reasons, since the spaces of
$\mathbf{R}$ and $\mathbf{P}$ grow exponentially with the system size $N_{q}$,
we need to truncate the spaces of the training set and error mitigating gates,
or resort to variational optimization approaches [278].

Deep learning method [279, 285, 280] trains a deep neural network to model a
noise channel utilizing the ”black box” nature of the neural network. The
training of the neural network requires the measurement results with both
noise and ideal cases. However, since classical simulations may not be
possible for large-scale non-Clifford quantum circuits, the training set can
only consist of measurement results of specific quantum circuits whose ideal
measurement outcomes are known [279, 285].

###  III.3 Measurement error mitigation

Measurement (or readout) error mitigation (MEM) schemes are designed to
improve the accuracy of the measurement results obtained from noisy quantum
devices.

####  III.3.1 MEM under classical noise assumption

Under the assumption of a classical noise model, MEM is usually applied via
classical post-processing, where measurement error is modeled by a stochastic
and invertible response matrix $\Lambda$ [286, 287, 288].

Classical noise models. The ideal quantum measurement in the computational
basis can be written in terms of positive operator valued measurement (POVM).
For a $n$-qubit system $\rho$, suppose the POVM operator is
$\Pi_{\bm{x}}=|\bm{x}\rangle\langle\bm{x}|$ satisfying $\Pi_{\bm{x}}\geq 0$
for $\forall{\bm{x}}$ and $\sum_{\bm{x}}\Pi_{\bm{x}}=\mathbb{I}$, where
$\bm{x}\in\mathbb{Z}_{2}^{N_{q}}$ refers to the POVM outcome. The probability
distribution of the measurement outcome is represented by a vector $\vec{p}$,
where the $\bm{x}$ term of $\vec{p}$ is the probability of obtaining the
outcome $\bm{x}$, given by $p(\bm{x}|\rho)=\text{Tr}(\rho\Pi_{\bm{x}})$
according to the Born’s rule.

If the POVM elements have no non-trivial off-diagonal terms, we can treat the
measurement noise channel as a classical noise channel, where the
transformation between ideal and noisy measurement probability distributions
can be written using a response matrix $\Lambda$:

$\vec{p}_{ideal}=\Lambda^{-1}\vec{p}_{noisy}.$ (37)

Here we use $\vec{p}_{ideal}$ and $\vec{p}_{noisy}$ to represent the ideal and
noisy probability distributions. And the element $\Lambda_{\bm{x},\bm{y}}$ of
the response matrix $\Lambda$ is defined as

$\Lambda_{\bm{x},\bm{y}}=\langle\bm{x}|\rho|\bm{y}\rangle,\quad\bm{x},\bm{y}\in\mathbb{Z}_{2}^{N_{q}},$
(38)

which can be estimated directly by preparing the computational state and
measuring it in the computational bases, or by using some tomographic means,
such as quantum detector tomography (QDT) [289, 286, 290]. The main idea of
QDT is to estimate an unknown set of fixed noisy POVM operators
$\\{\hat{\Pi}_{\bm{x}}\\}_{\bm{x}}$ with a set of well-known states
$\\{\rho_{i}\\}_{i}$. Once the noisy POVM operators are reconstructed, we can
extract the response matrix elements using

$\hat{\Pi}_{\bm{x}}=\sum_{\bm{y}}\Lambda_{\bm{x},\bm{y}}\Pi_{\bm{y}},\quad\forall\bm{x},$
(39)

where $\Pi_{\bm{y}}$ denotes the error-free POVM operator.

In practice, however, for non-classical noise, the estimated measurement
probability distribution obtained from
$\vec{p}_{est}=\Lambda^{-1}\vec{p}_{noisy}$ may be unphysical due to the
negative terms of the vector $\vec{p}_{est}$. Moreover, the statistical
uncertainties in the response matirx can be amplified through simply
inversion, similar to the challenges faced in high-energy physics. Unfolding
methods are thus introduced to readout error mitigation [291, 292], which show
robustness to some failure cases of matrix inversion and least-squares.

Simplified classical models. Dealing with the full response matrix of size
$2^{N_{q}}\times 2^{N_{q}}$ is not scalable beyond hundreds of qubits.
Therefore, in order to avoid the exponential consumption of measurements, many
scalable approaches [287, 293, 294] have been proposed to simplify the
response matrix model.

Tensor Product Noise (TPN) model. The simplest model assumes that the noise
acts independently on each qubit and is called the Tensor Product Noise (TPN)
model [287, 295, 296]. After such a simplification, the response matrix can be
described in terms of the tensor product of $N_{q}$ single-qubit response
matrix as

$\Lambda^{\text{TPN}}=\prod_{i=1}^{N_{q}}\left(\begin{matrix}1-p_{i}&q_{i}\\\
p_{i}&1-q_{i}\\\ \end{matrix}\right),$ (40)

where $p_{i}$ and $q_{i}$ denote the single-qubit readout error probability of
$1\to 0$ and $0\to 1$, respectively. We can see that the number of
measurements required to construct the TPN response matrix is reduced from
$\mathcal{O}(2^{N_{q}})$ to $\mathcal{O}(N_{q})$. However, in the TPN model,
the multi-qubit readout error correlations are neglected.

Continuous Time Markov Processes (CPTP) noise model. Bravyi et al. [287] have
extended the TPN model to the CTMP noise model, which takes the correlated
errors into account. It defines a matrix exponential $\Lambda=e^{G}$ to
describe the response matrix, and the generator
$G=\sum_{i=1}^{2{N_{q}^{2}}}r_{i}G_{i}$ models the readout error. Here
$r_{i}\geq 0$ refers to error rate and $G_{i}$ represents the single-qubit or
two-qubit readout error. For example, the generator describing falsely
measuring $|00\rangle$ to $|11\rangle$ on two qubits is $|11\rangle\langle
00|-|00\rangle\langle 00|$. The MEM with the CPTP model has been
experimentally demonstrated on the superconducting device [287].

Subspace Reduction. When the measured probability distribution only contains a
few principal bit strings with high probability, we can mitigate readout
errors in a renormalized subspace defined by the observed probability
distribution [293, 296]. The subspace reduction of the full $2^{N_{q}}$
dimensional response matrix efficiently circumvents the original exponential
overhead and avoids matrix inversion by using the matrix-free iterative
methods [293].

Readout symmetrizing. Due to the observation that mismeasuring $|1\rangle$ to
$|0\rangle$ is more frequent than $|0\rangle$ to $|1\rangle$, Several studies
symmetrizing the readout with targeted Pauli $\hat{X}$ gates are proposed [49,
297, 298, 299].

Readout balancing. We can exploit the readout asymmetry in the practical
cases, which aims to rebalance the measurement outcomes by minimizing the
expected number of qubits in the $|1\rangle$ state [297].

Bit-flip averaging. Bit-flip averaging [298] reduces the calibration overhead
by an exponential factor without making any specific assumptions about the
classical noise model for MEM. It applies Pauli $\hat{X}$ gates on the
randomly chosen qubits before measurements and offsets the effect of $\hat{X}$
gates using classical post-processing. This leads to a symmetrized effective
response matrix $\hat{\Lambda}$ which contains only $\mathcal{O}(2^{N_{q}})$
free parameters compared to $\mathcal{O}(2^{2N_{q}})$ in the full matrix.
Another independent work [299] utilizes the same bit-flipping protocol. It
transforms the bias between the ideal and noisy expectation value to a
multiplicative factor, and eliminates it by deviding the noisy expectation by
an estimated factor obtained in the calibration procedure.

####  III.3.2 Quantum noise model for MEM

Almost all of the readout mitigation approaches discussed above are applied
under the classical noise model assumption. However, quantum measurements
inevitably suffer from the quantum coherent noise [286, 289, 300]. With the
presence of quantum noise, the noisy POVM operators have non-trivial off-
diagonal values, so the classical assumption no longer holds. Fitting the
difference between two specific measurement statistics to the Fourier series
can effectively detect the presence of quantum noise [300]. Meanwhile, many
approaches, such as $IZ$ dephasing [286, 300], can be applied to eliminate the
quantum noise, so that the effective measurement device has only classical
noise. Then we can utilize those MEM methods mentioned above to deal with the
classical noise.

####  III.3.3 MEM combined with other techniques

MEM using QEC. On devices whose readout errors dominate over the entangling
gate errors, one can combine readout error mitigation with quantum error
correction to actively reduce readout errors on a shot-by-shot basis, called
active readout error mitigation [301].

Combination with readout compression. The compression readout compresses the
large-scale quantum state into one qubit and recovers the state amplitude
populations from the one-qubit measurement results [233]. Thus, after appling
the compression readout technique, only single-qubit measurements are
performed, so this method is free of correlated measurement errors and easy to
combine with other MEM methods.

Neural network model for MEM. We can also utilize a trained artificial neural
network to characterize an non-linear model [302] for measurement noise, which
captures both classical and quantum noise. The neural network model for MEM is
experimentally performed on the IBM’s 5-qubit quantum device [302].

###  III.4 Purification-based QEM scheme

Recently, a range of purification-based QEM schemes [242, 303, 304, 305, 306,
307, 308, 309, 310] have been proposed to extract the ideal state component
from the noisy state, without any requirements on the priori knowledge of the
quantum noise channel. These schemes assume that the state of interest is
usually a pure state, while the noise tends to corrupt it into a mixed state.

The two typical methods are Virtual Distillation (VD) [242] and Error
Suppression by Derangement (ESD) [303], which use multiple copies of the noisy
state to reduce the state preparation error. Given an error-free state
$|\psi_{0}\rangle\langle\psi_{0}|$, the noisy output state $\rho$ can be
expressed via the spectral decomposition

$\rho=(1-\lambda)|\psi\rangle\langle\psi|+\lambda\rho_{\text{err}},$ (41)

where $\lambda\geq 0$ denotes the error rate. The pure state
$|\psi\rangle\langle\psi|$ denotes the dominant eigenvector which may not
match the ideal state $|\psi_{0}\rangle\langle\psi_{0}|$ due to the coherent
error [311]. And the mixed state $\rho_{\text{err}}$ refers to the noisy
component consisting of the weighted sum of other eigenvectors in the spectral
decomposition. Using $M$ copies of noisy state, the $M$-degreed mitigated
expectation value is

$\displaystyle\langle\mathcal{O}\rangle_{\text{em}}=\frac{\text{Tr}(\rho^{M}O)}{\text{Tr}(\rho^{M})}$
$\displaystyle=\frac{(1-\lambda)^{M}\langle\psi|O\psi\rangle+\lambda^{M}\text{Tr}(O\rho_{\text{err}}^{M})}{(1-\lambda)^{M}+\lambda^{M}\text{Tr}(\rho_{\text{err}}^{M})}$
(42) $\displaystyle\approx\langle\psi|O\psi\rangle+\mathcal{O}(\lambda^{M})$

where $\mathcal{O}(\lambda^{M})$ is the approximation error which decays
exponentially with the number of copies.

Figure 4: The illustration of two implementations for virtual distillation.

We usually don’t need to prepare the exact $\rho^{M}$ state to compute the
$\langle\mathcal{O}\rangle_{\text{em}}$. An illustration of the implementation
of VD and ESD is given in Fig. 4. They both utilize the cyclic shift operator
(also called derangement operator) $S^{(M)}$ on the system with size $MN_{q}$.
Using the property of $S^{(M)}$

$S^{(M)}|\psi_{1}\rangle|\psi_{2}\rangle\dots|\psi_{M}\rangle=|\psi_{M}\rangle|\psi_{1}\rangle\dots|\psi_{M-1}\rangle,$
(43)

we have

$\frac{\text{Tr}(O\rho^{M})}{\text{Tr}(\rho^{M})}=\frac{\text{Tr}(O^{i}S^{(M)}\rho^{\otimes
M})}{\text{Tr}(S^{(M)}\rho^{\otimes M})},$ (44)

where $O^{i}$ refers to the observable $O$ acting on an arbitrary subsystem
$i$. The VD shown in Fig. 4(a) applies a layer of decomposition gates $B_{i}$
to compute Eq. 4, with a joint measurement on $M$ copies. While the ESD given
in Fig. 4(b) adds an ancilla qubit and uses controlled-$S^{M}$ gates to
estimate the expected value $\text{Tr}(O\rho^{M})$ (and for
$\text{Tr}(\rho^{M})$ setting $O=I$) by

$\text{Tr}(O\rho^{M})=2P_{0}-1.$ (45)

Here we denote the probability of measuring the ancilla qubit in the
$|0\rangle$ state by $P_{0}$.

The VD (or ESD) scheme is subject to some limitations. First, the coherent
dismatch between the target pure state $|\psi_{0}\rangle\langle\psi_{0}|$ and
the dominant eigenvector $|\psi\rangle\langle\psi|$ of the noisy state
contributes to a noise floor regardless of the increasing of $M$. Second, the
increased overhead of qubits and controlled gates with large $M$ is hard to
afford for the NISQ devices. Taking these limitations into account, many
alternative protocols have been explored.

####  III.4.1 Resource-efficient variants for purification-based schemes

Realization with classical shadows. Classical shadows [312, 313, 314] are
protocols used to efficiently predict many different properties, in particular
linear functions. Considering a trade-off between qubits (and controlled
gates) and measurements overhead, many methods implement the VD(or ESD) scheme
with classical shadows [307, 308].

Dual-state purification. This protocol [309] purifies states using
$\frac{\rho\bar{\rho}+\bar{\rho\rho}}{2}$, and the $\bar{\rho}$ is the dual
state of $\rho$. The ideal expectation value of observable $O$ is estimated as

$\langle\mathcal{O}\rangle_{\text{em}}=\frac{\text{Tr}(\frac{\rho\bar{\rho}+\bar{\rho\rho}}{2}O)}{\text{Tr}(\frac{\rho\bar{\rho}+\bar{\rho\rho}}{2})}.$
(46)

It works when the noisy state $\rho$ and its dual state $\bar{\rho}$ share the
same dominant eigenvector in their spectral decompositions. Dual-state
purification is efficient in terms of qubit overhead since its implementation
requires at most one ancilla qubit.

Combination with active qubit resets. The active qubit reset technique is
enabled by major quantum computing architectures including superconducting
qubit and trapped-ion devices. With the use of active qubit resets, we can
achieve a similar error suppression using 2$N$ \+ 1 qubits compared to $MN$ \+
1 qubits with the original VD [310], which shows a space-time trade-off in the
computational resources.

####  III.4.2 Generalized purification methods

In order to deal with the coherent dismatch, we can consider a general
polynomial function
$f(\rho,M)=c_{0}+c_{1}\rho+c_{2}\rho^{2}+...+c_{M}\rho^{M}$ [315, 307].
Recently Xiong et al. [306] provides a more general framework called
permutation filters for those schemes using permutations. It uses a filter

$y_{\text{filter}}=\frac{\text{Tr}(O\mathcal{F}_{\vec{\alpha}}(\rho,M))}{\text{Tr}(\mathcal{F}_{\vec{\alpha}}(\rho,M))}$
(47)

to substitute the polynomial
$\frac{\text{Tr}(\rho^{M}O)}{\text{Tr}(\rho^{M})}$, where $\mathcal{F}$
denotes some complex function and the parameter optimization for
$\vec{\alpha}$ is proven to be an invex problem, which can always converge to
the global optimum.

###  III.5 Quantum subspace expansion

Quantum subspace expansion (QSE) is first proposed to explore the excited
states, and it can also contribute to error mitigation in VQE with additional
classical resources [95, 316, 111, 317]. The QSE effectively mitigates
coherent errors due to the imperfect variational optimization. Given a noisy
ground state $\rho$ prepared on the quantum device and the target Hamiltonian
$H$, we consider a variational subspace spanned by $\rho$, as

$\displaystyle\\{\rho_{\text{sub}}=\sum_{i,j=1}^{M}c_{i}c_{j}^{*}P_{i}\rho
P_{j}|c_{i}\in\mathcal{C},P_{i}\in\\{I,X,Y,Z\\}^{\otimes N_{q}}$ (48)
$\displaystyle\text{and}\ \text{Tr}(\rho_{\text{sub}})=1\\},$

where the parameters $\vec{c}=(c_{1},...,c_{M})$ denote the expansion
coefficients and $M$ refers to the number of expansion terms. Then, to obtain
the error-mitigated value, we need to solve the following minimization problem

$\displaystyle\text{min}_{\vec{c}}\\{\text{Tr}(\rho_{\text{sub}}H)\\}$ (49)
$\displaystyle\text{such that Tr}(\rho_{\text{sub}})=1.$

The spectrum of the Hamiltonian within the classically expanded subspace can
be calculated as the solution to a generalised eigenvalue problem

$\hat{H}\vec{c}=E\hat{B}\vec{c}.$ (50)

Here, $E$ is a diagonal matrix with elements representing the eigenenergies.
$\hat{H}_{ij}=\text{Tr}(P_{i}\rho P_{j}H)$ and
$\hat{B_{ij}}=\text{Tr}(P_{i}\rho P_{j})$ both can be efficiently estimated
via quantum computers. With the optimized parameter obtained from
minimization, we can compute the error-mitigated expectation value using

$\langle O\rangle_{em}=\sum_{i,j=1}^{M}c_{i}c_{j}^{*}\text{Tr}(P_{i}\rho
P_{j}H).$ (51)

Note that the efficiency of the QSE method depends on the number of expansion
terms $M$, thus it may not effective in suppressing stochastic errors due to
the requirement on exponential expansion terms to project the noisy state to
the error-free subspace [131, 95, 315]. While when combined with certain
properties of the target state such as symmetry [318], the QSE can handle
stochastic errors more efficiently because the complexity of constructing the
projection subspace is reduced. Beside, a generalised framework of QSE [315]
has extended the expansion operators from Pauli operators to more general
operators relevant to the noisy state.

###  III.6 Verification methods

Verification methods [316, 318, 319] exploit the knowledge of inherent
symmetry within the quantum system, such as the spin symmetry in quantum many-
body physics. This methods focus on the errors that place state outside the
symmetry-preserving subspace.

The symmetry commonly used is Pauli symmetry $\hat{S}$, which can be
effectively estimated. Suppose the Pauli symmetry of the system with size $n$
is $\hat{S}=\prod\limits_{i=1}^{M}\hat{S}_{i}$,
$\hat{S}_{i}\in\mathbb{P}^{\otimes N_{q}}$, where $M$ is the number of non-
trivial terms in $\hat{S}$, and $\mathbb{P}^{n}$ refers to the $n$-qubit Pauli
group. The verification of the symmetry can be achieved by dropping the
circuit runs which fail the verification or by using the post-processing
appraoches, as we will discuss below.

####  III.6.1 Symmetry verification

Figure 5: Symmetry Verification.

Symmetry verification [316, 318] discards the circuit runs failed to pass the
verification, which costs additional measurement overhead. It utilizes an
ancilla qubit interacting with each qubit in the system register, which is of
the form shown in Fig. 5(a). The gates $\\{B_{i}\\}_{i=1,2}$ are the basis
transformation gates that map the eigenstate of symmetry $\hat{S}$ with
eigenvalue $s$ to $|0\rangle$ state. If the ancilla qubit reads 1, we discard
this circuit run. Note that the circuit can only detect odd number of
error(s). In the low-cost version [318] of symmetry verification shown in Fig.
5(b), it shuffles the ancilla qubit along the system register, which needs
only local CNOT and SWAP two-qubits gates.

Both cases require $\mathcal{O}(M)$ circuit depth to ensure the entanglement
between the ancilla qubit and each register qubit individually, which is
general intractable in quantum circuit. Moreover, the verification circuit
applied on the noisy state may introduce extra errors, which reduces the
reliability of verification. An alternative way to perform symmetry
verification will be detailed below.

####  III.6.2 Symmetry expansion

We can also implement verifications via post-processing approach [318].
Usually if the target quantum state is the eigenstate of $\hat{S}$ (or in the
eigen-subspace) with eigenvalue $s$, we can construct projector valued
measurement to project the output state $\rho$ to the $\hat{S}=s$ subspace.
Suppose the corresponding projector is ${\hat{\Pi}_{s}}$ and we have
$\hat{\Pi}_{s}\hat{\Pi}_{s}=\hat{\Pi}_{s}$. The effective density matrix after
projection is

$\rho_{s}=\frac{\hat{\Pi}_{s}\rho\hat{\Pi}_{s}}{\text{Tr}(\hat{\Pi}_{s}\rho\hat{\Pi}_{s})}=\frac{\hat{\Pi}_{s}\rho\hat{\Pi}_{s}}{\text{Tr}(\hat{\Pi}_{s}\rho)}.$
(52)

If we will measure an observable $O$ which commutes with the symmetry
$\hat{S}$ and $[O,\hat{S}]=0$ and $[O,\hat{\Pi}_{s}]=0$. Then the expectation
value under the projected state is

$\langle
O_{em}\rangle=\text{Tr}(O\rho_{s})=\frac{\text{Tr}(O\hat{\Pi}_{s}\rho)}{\text{Tr}(\hat{\Pi}_{s}\rho)}.$
(53)

Note that the probability of passing the verification is just the probability
of projecting the noisy state to the $\hat{S}=s$ subspace, which is
$\text{Tr}(\hat{\Pi}_{s}\rho)$. Hence we need
$\Upsilon_{\text{SE}}=\frac{1}{1-\text{Tr}(\hat{\Pi}_{s}\rho)}$ more
measurement shots as the sampling overhead for mitigation.

We can combine the post-processing approaches with QSE to further improve the
efficiency of mitigation, which called s-QSE [318]. If we choose
$\hat{\Pi}_{s}=\frac{I+s\hat{S}}{2}$, we will have

$\langle
O_{\text{em}}\rangle=\frac{\text{Tr}(O\rho)+s\text{Tr}(O\hat{\Pi}_{s}\rho)}{1+s\text{Tr}(\hat{\Pi}_{s}\rho)}.$
(54)

Then the original minimization problem (as we mentioned in the Section III.5)
is reformulated to a generalized eigenvalue problem, and the terms
$\text{Tr}(O\rho)$, $\text{Tr}(O\hat{\Pi}_{s}\rho)$ and
$\text{Tr}(\hat{\Pi}_{s}\rho)$ can be efficiently estimated via quantum
devices. The s-QSE method has been experimentally demonstrated to mitigate
errors in the VQE of $\text{H}_{2}$ with two transmon qubits [320].

Cai has extended the post-processing approaches to a general framework called
symmetry expansion [319], which encompasses a broader range of symmetry-based
error mitigation methods. For example, VD can be seen as a special case of
symmetry expansion using permutation symmetry. Moreover, the particle number
is also applied as symmetry in the post-processing approach[321].

####  III.6.3 Other verification methods

Verified Phase Estimation. It applies phase estimation to estimate expectation
values while effectively post-selects for the system register to be in the
starting state [322]. It can also be adapted to the case without the use of
control qubits, which simplifies the control circuits.

Pauli Check Sandwiching. The Pauli Check Sandwiching scheme [323] applies
multiple pairs of Pauli checks to detect the occurrence of errors, and obtains
the mitigated results using post-selecting. Each pair of Pauli checks uses one
ancilla qubit to detect a component of the error operator.

###  III.7 Mitigating the coherent errors

Coherent errors refer to the imperfect or unwanted unitary rotations acting in
the circuits, which can be modeled as

$U(\vec{\theta})=e^{-\frac{i}{2}\vec{\theta}\cdot\vec{\sigma}},$ (55)

where $\vec{\theta}=(\theta_{1},...,\theta_{4^{N_{q}}})$ quantifies the
strength of the coherent error on the $4^{N_{q}}$ Pauli bases. Coherent errors
map the noiseless pure states to another pure states, since unitary operators
maintain the quantum coherence of the states. While they are purity-
preserving, it can pose a threat to the reliable multi-qubit quantum
computation. Now we will introduce several typical methods for mitigating
coherent errors [1, 324].

####  III.7.1 Randomized compiling

Randomized compiling (RC) [257, 258] is designed to tailor coherent errors
into stochastic Pauli errors, and we can combine RC techniques with other QEM
schemes [281, 261, 271].

When the two-qubit Clifford gate errors dominant over other types of error, we
can mitigate the coherent error in the two-qubit gates using RC. After
sandwiching each two-qubit gate between randomly sampled Pauli gates and
compiling those twirling Pauli gates into the original single-qubit gates
(which can be implemented on the classical computers in advance), the newly
generated circuits are logically equivalent to the bare circuits with the same
circuit depth. And the averaged results over many logically-equivalent
circuits is exactly the desired result with tailored noise.

####  III.7.2 Hidden inverses

Hidden inverses (HI) [325, 326, 327] is first introduced to mitigate certain
coherent errors (such as over-rotations and phase misalignment) in the
trapped-ion quantum computer [325], and then it has been extended to implement
on the superconducting hardwares [327]. The HI method relies on the self-
adjoint unitary operators satisfying $U=U^{\dagger}$ or self-inverse unitary
operators with $U=U^{-1}$, such as $H$ and $CX$. Though these gates represent
the same operation in the error-free case, they may suffer from noise with
different strength in practice, due to the different compiling ways for $U$
and $U^{\dagger}(U^{-1})$. Then we can mitigate coherent errors by a local
optimization which determines to construct the same gate from the original
elementary gate sequence or the sequence of the inverted (conjugate) gate.

###  III.8 Combination with quantum error suppression techniques

####  III.8.1 Combination with QEC

Individual error reduction. The individual error reduction method [328] make
sufficient use of the limited QEC techniques which have been well demonstrated
on several qubits. It reduces error on each qubit respectively and obtains the
mitigated expectation value of the observable by post-processing. Consider the
Lindblad master equation describing the noise process

$\frac{\partial\rho}{\partial
t}=\mathcal{L}(\rho)=\sum_{i}\mathcal{L}_{i}(\rho),$ (56)

where the $\mathcal{L}_{i}$ denotes the Lindblad operator acting on the single
qubit. We can give a solution for this equation with the duration $\tau$ of
the noise process to the first order approximation

$\displaystyle\rho(t+\tau)=\mathcal{E}_{\tau}(\rho(t)),$ (57)
$\displaystyle\mathcal{E}_{\tau}(\rho)=e^{[\tau\mathcal{L}(\rho)]}\approx
1+\tau\sum{i}\mathcal{L}_{i}(\rho),$

where $\mathcal{E}_{\tau}$ is the Lindblad evolution operator. Suppose the
error on the $i$-th qubit is reduced by a known factor $g_{i}$ via QEC, then
the corresponding Lindblad evolution operator $\mathcal{E}_{\tau}^{i}$ turns
to

$\mathcal{E}_{\tau}^{i}(\rho)\approx 1+\tau\sum_{j\neq
i}\mathcal{L}_{j}(\rho)+\tau(1-g_{i})\mathcal{L}_{i}(\rho).$ (58)

Hence, after QEC applied on the $i$-th qubit, the expectation value of
observable $O$ is $\langle O\rangle_{i}=\text{Tr}(\rho^{i}O)$ with the output
density matrix $\rho^{i}$. We can estimate the ideal case of the expectation
value of observable $O$ by a linear combination of different $\langle
O\rangle_{i}$:

$\langle O\rangle_{\text{Ind}}=\text{Tr}(\rho^{\text{Ind}}O)=\langle
O\rangle-\sum_{i}\frac{1}{g_{i}}(\langle O\rangle-\langle O\rangle_{i}).$ (59)

Note that the individual error reduction has suppressed errors to the first
order and can be combined with other error mitigation methods for higher-order
error suppression. However, it relies on the accurate estimation of the factor
$g_{i}$ which results in extra cost for characterizations. And additional
quantum resources (qubits or gates) are needed to implement QEC on a single
qubit.

Code space projection. The code space projection method [329, 307] decodes
errors on logical qubits via post-processing, without additional qubits and
operations for syndrome measurements. To encode $k$ logical qubits with $N$
physical qubits, we use a stabilizer code $[[N,k]]$ defined by a stabilizer
group $S=\langle S_{1},S_{2},...,S_{N-k}\rangle$ with generators
$\\{S_{i}\\}_{i}$. And the code subspace is determined by the projection
operator

$\Pi=\prod_{i=1}^{N-k}\frac{\mathbb{I}+S_{i}}{2}=\frac{1}{2^{N-k}}\sum_{j}M_{j},$
(60)

where $M_{j}\in S$ refers to the element in the group. Given a logical
observable $O$ which can be decomposed by Pauli operators
$O=\sum_{m}\gamma_{m}P_{m}$, this method can reduce errors which take the
physical state outside the code space. Then we estimate the ideal expectation
value using

$\langle
O\rangle_{\text{CS}}=\frac{1}{c2^{N-k}}\sum_{j,m}\gamma_{m}\text{Tr}(\rho
M_{j}P_{m}).$ (61)

Stochastic and deterministic subspace expansion schemes have been proposed
[329] for estimation of Eq. 61. We can also realize the code space projection
using classical shadow techniques [307].

####  III.8.2 Dynamic Decoupling

Dynamical decoupling (DD) [330, 331, 332, 333, 334, 335] methods are designed
to suppress decoherence caused by system-environment interaction. The main
idea of DD is to apply specific pulse sequence called DD sequence to the idle
qubits, which keeps the overall logic of the circuit unchanged. Introducing
the DD techniques into the QEM schemes can improve the circuit performance
while it dosen’t result in additional sampling overhead [271, 336].

####  III.8.3 Single-qubit gate scheduling

Rather than adding additional gates to the circuit as DD does, single-qubit
gate scheduling [337] optimizes circuits by scheduling within idle windows
which are periods of idle qubit waiting for the next operation. As Late As
Possible (ALAP) scheduling is a typical scheduling technique used as a default
approach in IBM Qiskit [338] to execute the single-qubit gates at the end of
the idle windows. We can also tune the positions of single-qubit gates within
idle windows by circuit slicing and inverting [337]. A variational QEM method
designed to perform DD and single-gate scheduling within the VQA framework has
been proposed recently [336], which avoids specific configuration selections
for DD sequences and gate positions.

###  III.9 Framework for QEM

Recently, much effort has been devoted to exploring the unified framework for
QEM schemes and the general bounds for sampling overhead [249, 251, 339, 250,
340]. Takagi et al. [251] described the QEM process as the concatenation of a
quantum-classical channel composed of quantum operations and independent POVM
measurements, and a classical-classical channel for the classical post-
processing. They also introduced the maximum bias to benchmark the worst-case
performance of the QEM strategies for an arbitrary state and observable, and
then derive a fundamental bound on the sampling overhead based on the maximum
bias. Other metric such as extraction rate was proposed to characterize the
cost-effectiveness of QEM schemes within the linear quantum error mitigation
[249] framework, where the process of error mitigation is regarded as
extracting the effective state out of the unmitigated state. Statistics
principles are also applied into a general QEM formalism [339] where an
approximate quadratic reduction with gate number $N$ has drawn in error
increase.

##  IV Quantum Circuit Compilation

Compilation in the classical computing community is the process of converting
a high-level programming language into a machine language that the computer
can understand and execute fluently. Similarly, in the field of quantum
computation, quantum circuit compilation (QCC) is used to transform a quantum
algorithm always in a mathematical form into its corresponding quantum circuit
which can be executed on a real quantum device.

Usually, QCC includes three stages as follows (see Fig. 6).

  1. 1.

Decomposion (also referred to as the synthesis). Usually, a mathematically
designed quantum algorithm can be represented by several $n$-qubit unitary
operators, such as the Grover iteration in Grover’s algorithm, the controlled
exponentiation in Shor’s algorithm, the Hamiltonian simulation in HHL
algorithm and so on. However, real quantum devices, especially NISQ devices,
usually provide only some elementary quantum gates, such as single-qubit gates
and CNOT gates. Thus, it is one of the most fundamental problems that how to
decompose a specific quantum algorithm or unitary operator into as few the
provided elementary quantum gates as possible for implementation.

  2. 2.

Optimization. Despite the rapid development of quantum hardware, the qubit
lifetime of some quantum systems is still not satisfactory and the quantum
noise is also inevitable [341, 7]. The optimization aims to further reduce the
depth and size of quantum circuits with the help of auxiliary qubits after the
decomposition process to alleviate the problem of short qubit lifetime. Thus,
the optimization seems like a kind of space-time trade-off, i.e., using
ancillary qubits (space) to efficiently reduce circuit size or depth (time).

  3. 3.

Mapping (also referred to as the qubit routing, qubit allocation, quantum
circuit transformation). The final stage is to transform the logical quantum
circuit to the physical quantum circuit which can be operated on a real
quantum device. Actually after decomposition and optimization, there are only
elementary quantum gates in the logical circuit, and the size and depth of
which are sufficiently optimized. For the quantum hardware with fully-
connected architecture, any logical circuit can be directly implemented.
However, for the quantum hardware with limited connectivity architecture [341,
7], we should find the correspondence between logical qubits and physical
qubits to complete the transformation.

Figure 6: Quantum circuit compilation. Quantum circuit compilation translates
the abstract quantum algorithm into a compiled circuit optimized for a
specific quantum hardware platform, enabling efficient and highly reliable
execution through the use of decomposition, optimization, and mapping methods.

Although QCC has been separated into the above three stages, sometimes there
is no exact order or boundaries. Alternatively, one can perform all the three
stages simultaneously to complete the QCC, such as the works in Refs. [342,
343, 344]. QCC is extremely important as a bridge between quantum algorithms
and quantum hardware devices, especially for NISQ devices. There have been a
series of results on QCC and we will present them according to the different
approaches of the three stages mentioned above.

###  IV.1 Decomposition

The first result in the aspect of Decomposition was published in 1995, where
all $n$-qubit unitary operators can be expressed as compositions of a finite
number $O(n^{3}4^{n})$ of single-qubit gates and CNOT gates exactly [345].
Soon, this upper bound was lowered to $O(n4^{n})$ by Knill [346]. Nine years
later, the upper bound for the exact decomposition of $n$-qubit unitary
operator has been improved to $O(4^{n})$ [347], which coincides with the
theoretical lower bound [348]. In 2004, M$\ddot{o}$tt$\ddot{o}$nen et al. came
up with Cosine-Sine Decomposition (CSD), i.e.,

$F_{n}^{l}(R_{a})=\left(\begin{array}[]{ccc}R_{a}(\alpha_{1})&&\\\ &\ddots&\\\
&&R_{a}(\alpha_{2^{n-1}})\\\ \end{array}\right).$ (62)

to express an arbitrary $n$-qubit unitary operator in terms of $n-1$-qubit
unitary operators [349], which further decreased the upper bound to
$4^{n}-2^{n+1}$ CNOT gates. Next in 2005, Shende et al. applied the NQ matrix
decomposition to implement an arbitrary $n$-qubit operator only using
$\frac{1}{2}\times 4^{n}-3\times 2^{n-1}+1$ CNOT gates [350], where the
$n$-qubit operator can be implemented by a circuit containing three uniformly
controlled rotations and four ($n-1$)-qubit operators. And then in 2006,
M$\ddot{o}$tt$\ddot{o}$nen et al. applied CSD recursively to yield a synthesis
algorithm, in which the number of CNOT gates is
$\frac{23}{48}4^{n}-\frac{3}{2}2^{n}+\frac{1}{3}$, which firstly decreased the
upper bound less than $\frac{1}{2}4^{n}$ CNOT gates [351, 352].

The analysis of lower bound was first came up with by Barenco et al. in [345]
as a conjecture, i.e., $\frac{1}{9}4^{n}-\frac{1}{3}n-\frac{1}{9}$, based on
dimension counting. In 2004 Shende et al. raised the lower bound to
$\frac{1}{4}(4^{n}-3n-1)$ by the technique of parameter counting [348], which
also implies a depth lower bound of $\Omega(4^{n}/n)$.

We can see that there is a gap on the optimal depth within the range of
$[{\Omega}(4^{n}/n),O(4^{n})]$ for general $n$-qubit circuit optimization
without ancillary qubits. Most recently, based on Gray code and unary encoding
[353], Sun et al. have proved a fact that any $n$-qubit unitary operator can
be implemented by a quantum circuit of size

$O(4^{n})$ (63)

and depth

$O\big{(}n2^{n}+\frac{4^{n}}{m+n}\big{)}$ (64)

with $m\leq 2^{n}$ ancillary qubits [354]. They used CSD repeatedly and
factored an arbitrary unitary operator into a sequence of uniformly controlled
gates. Namely, they presented that any $n$-qubit unitary operator $U$ can be
decomposed as
$U=V^{n}_{n}(0)\cdot\prod_{i=1}^{2^{n-1}-1}V^{n}_{n-\zeta(i)}(i)\cdot
V^{n}_{n}(2^{n-1}),$ where $\zeta(n)$ is the Ruler function defined as
$\zeta(n)=\max\\{k:2^{k-1}|n\\}$, $V^{n}_{k}$ denotes an $n$-qubit uniformly
controlled gate whose index of target qubit is $k$, and different $i$ in
$V_{k}^{n}(i)$ denote different forms of $n$-qubit uniformly controlled gates
despite the same target qubit $k$. Combined with their implementation of any
$n$-qubit uniformly controlled gate, i.e., any $n$-qubit $V_{k}^{n}(i)$ can be
implemented by a circuit of size $O(2^{n})$ and depth
$O\big{(}n+\frac{2^{n}}{n+m}\big{)}$, they derived the above result. What we
should note is that their size and depth are simultaneously optimized to
asymptotically optimal based on the lower bound in [348]. In specific without
any ancillary qubits, the depth can be optimized to $O(4^{n}/n)$, which is
exactly the lower bound.

Besides the most general decomposition on $n$-qubit unitary, there have
already been some results about the special cases. When the elements of
unitary belong to the ring $\mathbb{Z}[\frac{1}{\sqrt{2}},i]$, the single-
qubit unitary can be constructed using $H$ and $T$ gates only [355], in which
the similar $n$-qubit case was conjectured to be implemented by the Clifford
and $T$ gates. This conjecture has been proved soon afterwards. An $n$-qubit
unitary operator has an exact representation over the Clifford+$T$ gate set
assisted by only one ancilla if and only if its entries are in the ring
$\mathbb{Z}[\frac{1}{\sqrt{2}},i]$, and the total gate count is $O(3^{2n}nk)$,
where $k$ is the denominator exponent of the unitary [356]. Besides the above
results on exact decomposition of unitary, there have already been some works
on approximate decomposition. Based on Solovay-Kitaev theorem, Dawson and
Nielsen designed an approximate algorithm, which runs in
$O(\log^{2.71}\frac{1}{\epsilon})$ time and produces as output a sequence of
$O(\log^{3.97}\frac{1}{\epsilon})$ quantum gates including only Hadamard,
controlled-not, and $\pi/8$ gates [357].

###  IV.2 Optimization

As seen in the rapid development of quantum hardware, the number of physical
qubits may increase much faster than the qubit lifetime. Thus, when optimizing
quantum circuits, we can consider how to reduce the circuit depth using
auxiliary qubits.

In 2000, Cleve showed the quantum Fourier transform can be approximated to the
depth of $O(\log n+\log\log(1/\epsilon))$ with sufficient ancillary qubits
[358]. Next in 2001, Moore et al. proved that a circuit of any size on $n$
qubits composed entirely of CNOT gates (CNOT circuit) can be parallelized to
$O(\log n)$ depth with $O(n^{2})$ ancilla qubits [359]. Similarly for the
circuit of controlled-Pauli and $H$ gate, they derived the same result. And
for the Clifford$+T$ circuit, Selinger in 2013 made use of four ancillary
qubits to represent Toffoli gate with a T-depth of 1 and a total depth of 7
[360]. Amy et al. described an algorithm to reduce the $T$-depth for
Clifford$+T$ circuit [361], and it has worst case runtime of cubic in the
number of $T$ gates, qubits and Hadamard gates. Next let us turn our eyes back
to the optimization of CNOT circuit because of its significance in quantum
computation [345, 362]. In 2020, Jiang et al. established an asymptotically
optimal space-depth tradeoff for the design of CNOT circuits, that is, any
$n$-qubit CNOT circuit can be parallelized to

$O(\max\\{\log n,\frac{n^{2}}{(n+m)\log(n+m)}\\})$

depth with $m$ ancillary qubits, where $m$ being any natural number [363].
Firstly, any $n$-qubit CNOT circuit can be regarded as an invertible matrix M
over the finite field $\mathbb{F}_{2}^{n\times n}$. Take a three-qubit CNOT
circuit shown in Fig. 7 as an example, the left circuit can be represented as
the right matrix.

Figure 7: The matirx representation of a CNOT circuit.

Based on the above representation, the CNOT circuit optimization problem will
be transformed into a corresponding parallel Gaussian elimination problem
mathematically. Combined with a standard technique in reversible computation,
they firstly constructed two CNOT circuits $\mathcal{C}_{1}$,
$\mathcal{C}_{2}$ with $2n$ and $3sn$ qubits respectively, where $1\leq s\leq
O(n/\log^{2}n)$. Specifically, for
$\textbf{x},\textbf{j}\in\mathbb{F}_{2}^{n}$,

$\mathcal{C}_{1}|{\textbf{x}}\rangle|{\textbf{j}}\rangle|{0}\rangle^{\otimes
3sn}=|{\textbf{x}}\rangle|{\textbf{j}\oplus\textbf{Mx}}\rangle|{0}\rangle^{\otimes
3sn},$
$\mathcal{C}_{2}|{\textbf{x}}\rangle|{\textbf{j}}\rangle|{0}\rangle^{\otimes
3sn}=|{\textbf{x}}\rangle|{\textbf{j}\oplus\textbf{M}^{-1}\textbf{x}}\rangle|{0}\rangle^{\otimes
3sn}.$

Then they applied $\mathcal{C}_{1},\mathcal{C}_{2}$ to transform the initial
state $|{\textbf{x}}\rangle|{0}\rangle^{\otimes n}|{0}\rangle^{\otimes 3sn}$
to the final state $|{\textbf{Mx}}\rangle$ after permuting the first and
second $n$ qubits. This boun d is tight and this algorithm can be directly
extended to stabilizer circuits with the reason that any stabilizer circuit
has a canonical form “H-C-P-C-P-C-H-P-C-P-C”, where H and P are one layer of
Hadamard gates and Phase gates respectively [364].

Near-term quantum devices are limited not only by the qubit lifetime, but also
by the connectivity between qubits. All the above optimization results focus
on the optimal depth, without considering the constraints on connectivity.
Thus, it is necessary to consider the quantum circuit optimization problem for
the quantum processors with sparsely connected structure. We should note that
the optimization here usually refers to reducing the circuit depth with the
help of ancillary qubits, which is quite different from the mapping introduced
in the next section, which focuses on reducing the number of SWAP gates added
to realize the logical circuit physically and keeping the number of qubits
constant.

Actually there have already existed some results on reducing the circuit size
under the limited connectivity architecture [342, 343, 344]. In 2019, based on
Gaussian elimination, Kissinger and Griend presented an algorithm which can
extract any $n$-qubit CNOT circuit to $2n^{2}$-size equivalent CNOT circuit if
the architecture contains a Hamiltonian path [342]. But if there is no
Hamiltonian in the architecture, their optimized size is $O(n^{3})$. Next in
Ref. [343], Nash et al. studied the CNOT circuit optimization on any connected
graph, and they proposed an algorithm which gives a $4n^{2}$-size equivalent
CNOT circuit for any $n$-qubit CNOT circuit. This bound has been reduced to
$2n^{2}$ [344]. Furthermore for size optimization, Wu et al. came up with an
algorithm, inspired by the result in [365], which can reduce the size of any
$n$-qubit CNOT circuit to

$O(n^{2}/\log\delta)$

on any connected graph with $\delta$ being the minimum degree, and this bound
is optimal when the limited connectivity architecture is a regular graph
[344]. In addition, they also optimized the depth of CNOT circuit with the
assistance of ancillary qubits. Specifically, their algorithm can reduce the
depth of any given $n$-qubit CNOT circuit to

$O(n^{2}/\min\\{m_{1},m_{2}\\}),$

with $3n\leq m_{1}m_{2}\leq n^{2}$ ancillary qubits and the limited
connectivity architecture is $m_{1}\times m_{2}$ grid graph [344]. And this
bound is asymptotically optimal when $m_{1}\times m_{2}=n$. Also this result
can be similarly generalized to any constant dimensional grid graph.

###  IV.3 Mapping

The study of mapping problem firstly focused on its complexity. In 2009, it
has been conjectured that mapping a non-LNN (Linear Nearest Neighbor) circuit
into an LNN circuit is NP-complete [366], because the number of possible
combinations of all permutations for all two-qubit gates in circuit is
$((n-1)!)^{k}$, where $n$ is the number of qubits and $k$ is the number of
two-qubit gates. Obviously, the number $((n-1)!)^{k}$ will become huge with
$n$ and $k$ increasing. In 2014, Shafaei et al. modeled the logical quantum
circuit as an interaction graph and the structure of quantum device as a
connectivity graph. Thus the mapping problem is a standard graph embedding
problem with connectivity and interaction graphs as the host and guest graphs,
respectively. And the objective is then to minimize the total distance between
adjacent nodes of the interaction graph. For a 2D grid connectivity graph, the
mapping problem is NP-complete [367]. In 2018, Siraichi et al. proved the
mapping problem, similar to the classical register allocation, is NP-complete
in general [350], which is also derived in [368]. Due to the NP-completeness
of mapping problem, researchers have taken more attention on how to reduce the
SWAP gates needed approximately. We will introduce and analyze the relevant
results in three aspects based on different approaches to solve the mapping
problem.

Mathematical programming. It is a most natural idea to reformulate the mapping
problem on some specific physical architectures into a corresponding
optimisation problem and to solve it with state-of-the-art tools and solvers
[369, 370, 367, 371, 372, 373]. In 2009, Hirata et al. considered a sorting of
the initial qubit order according to a function they defined at first, then
sorted this qubit order to the final one. Thus the objective is to minimize
the sum of SWAP gates of the first and second step [366]. Next in 2013,
Shafaei et al. regarded the mapping problem as improving locality of a given
quantum circuit, i.e., the minimum linear arrangement (MinLA) problem in graph
theory, and they showed the effectiveness of this approach for quantum Fourier
transformation and reversible benchmarks through experimental results [370].
Also they focused on 2D grid architecture [367] and refactored the mapping
problem to Mixed Integer Programming problem. In 2015, Lye et al. formulated
the mapping problem on multi-dimensional quantum architectures as Pseudo-
Boolean Optimization (PBO) problem, and used a state-of-the-art PBO solver
[374] to achieve a minimal number of SWAP gates [369]. In 2019, Wille et al.
formulated the mapping problem as a symbolic optimization problem that is
solved using reasoning engines like Boolean satisfiability solvers, and they
provided a method that maps small-size logical circuits to IBM’s QX
architectures with a minimal number of SWAP and $H$ gates [371].

Heuristic algorithms. To adapt real quantum devices with complex connectivity
architectures, a series of heuristic algorithms [341, 375, 376, 377, 378, 379,
380, 381] on the mapping problem were developed. Actually for mapping problem,
the most important thing is to determine how to insert the SWAP gate
effectively and efficiently. Here “effectively” means mapping the logical
quantum circuit to physical quantum circuit and running that circuit on a
quantum hardware device successfully, while “efficently” means adding as few
SWAP gates as possible during the mapping procedure. These two aspects
correspond to the two steps in the heuristic algorithm design. The first step
is to determine the candidate SWAP-gate sets at present, and the second step
is to design the evaluation function that can evaluate the effectiveness of
the present SWAP gate.

In 2011, Saeedi et al. proposed two heuristic methods to determine which
qubits should be reordered, a template matching optimization and an exact
synthesis approach [382], which, from experimental results, improved the
quantum cost by more than $50\%$ on average compared to the naive method. In
2016, Wille et al. considered the effect of inserting SWAP gate on the
following two-qubit gates and came up with a look-ahead scheme, which reduced
the number of SWAP gates of $56\%$ in the best experimental evaluation [383].
This look-ahead scheme also leads to the improvement for 2D architectures. In
2019, Zulehner, the first-place winner of IBM Q Award, devised a much better
algorithm in which quantum gates firstly would be partitioned into layers such
that each layer contains only gates acted on distinct sets of qubit. Then they
employed the $A^{*}$ algorithm [384] to search and determine how to insert the
SWAP gate for the current layer, which can help to decrease the circuit depth
simultaneously. But this algorithm is not suitable to large-scale logical
circuit because of the space and time cost of $A*$ algorithm. In order to
improve the efficiency of the searching algorithm, in 2019, Li et al. proposed
a SWAP-based BidiREctional heuristic search algorithm [375]. Instead of
searching for a mapping in the full space, they only searched for SWAPs
associated with the qubits and designed a heuristic cost function to help find
the SWAP that can reduce the sum of distances between each qubit pairs in the
front layer. Based on this idea, in 2020, Zhou et al. calculated the costs of
SWAP gates on all front layers heuristically and chose the one with the
minimum cost [376], which would further reduce the complexity in time
polynomial in all parameters. Most recently, they proposed a Monte Carlo Tree
Search framework to tackle the mapping problem [377], which enabled the search
process to go much deeper.

Machine learning techniques. Machine learning techniques have also been
exploited to provide a more precise evaluation tool for mapping problem [379,
380, 381], as they can provide a more precise assessment of the effects of the
present SWAP gate. In 2020, Pozzi et al. framed the mapping problem as a
reinforcement learning problem, then applied the combinatorial optimization
techniques to search for the action leading to the highest-quality under their
defined quality function, and finally scheduled SWAP gates [378]. In Ref.
[379], graph neural networks were used as an improved architecture to help
guide the tree search used to find an optimal set of SWAP gates, then they
used an array of mutex locks to represent the current parallelization state to
help to reduce the depth of the output circuits. In Ref. [380], Zhou et al.
took both short- and long-term rewards into consideration to design a scoring
mechanism, and proposed a search algorithm which is polynomial in all relevant
parameters.

##  V benchmarking Protocols

Benchmarking is a technique to assess the performance and capabilities of
quantum devices, and is therefore critical to accelerating progress in quantum
computing. Design a “good” benchmarking protocol for quantum devices is
complicated for the following reasons: 1) There are different physical ways to
implement quantum computing, such as ion traps, photons, superconducting, etc,
and each has its own strengths and weaknesses. It is difficult for us to
design a comprehensive protocol that can evaluate the performance of each
physical system well. 2) Different types of applications may have different
forms of assessment and core metrics. 3) When the benchmarking protocol is
applied to large-scale systems, significant resource consumption may be
required.

To date, a number of benchmarking protocols have been proposed and widely used
in experiments. These benchmarking protocols can be broadly divided into two
major categories, gate-level benchmarking and circuit-level benchmarking,
which will be introduced in this section.

Figure 8: Benchmarking methods for different levels. These methods can be
divided into three levels, basic performance benckmarking, gate-level
benckmarking and circuit-level benchmarking. At the lowest level, the basic
performance of quantum devices, such as $T_{1}$, $T_{2}$, and chip
connectivity, is calibrated, while in the gate-level benckmarking, each
individual quantum operations, such as single-qubit gate and two-qubit gate,
are optimized to enable the implementation of quantum circuits. In the
circuit-level benchmarking, the performance of the quantum system is tested by
executing complete quantum circuits or particular applications. Higher-level
benchmarks reflect a integrated performance, while lower-level benchmarks are
the basic steps performed prior to the utilization of the processor.

###  V.1 Gate-level benchmarking

####  V.1.1 Randomized benchmarking (RB)

Quantum process tomography (QPT) is a standard approach to fully
characterizing quantum process, but it is sensitive to the state preparation
and measurement (SPAM) errors and infeasible for large systems. RB is a
technique for estimating the average fidelity of a set of quantum gates
without relying on accurate SPAM [385, 386].

RB is originally proposed for random unitary gates [385, 387, 388], and now is
developed as a widely used experimental technique for characterizing the
average error of quantum operations. Typically, RB is executed with gates from
the Clifford group, which is referred to Clifford randomized benchmarking
(CRB) [386, 389, 33, 390]. A CRB protocol consists of the following steps:

Step 1: Generate RB sequences. First randomly sample a sequence of gates of a
fixed length $m$ from the Clifford group $G$ on $n$-qubits, and then a global
inversion gate is added to return the qubits to the initial state in the ideal
case without noise. For each length $m$, we choose $K_{m}$ RB sequences.

Step 2: Execute the RB sequences on the noisy quantum devices. Measure the
output state to estimate its overlap with the initial state, which is named as
survival probability $\alpha_{m}$. Repeat $K_{m}$ sequences to obtain the a
single average survival probability $\overline{\alpha_{m}}$.

Step 3: Fit the results. Repeat Steps 2 for various sequence lengths $m$ to
produce a list of average survival probabilities
$\\{\overline{\alpha_{m}}\\}_{m}$, and then fit
$\\{\overline{\alpha_{m}}\\}_{m}$ to the a single exponential decay model:

$\overline{\alpha_{m}}\ =A_{0}p^{m}+B_{0}$ (65)

where $A_{0}$ and $B_{0}$ absorb state preparation and measurement errors as
well as an _edge effect_ from the error on the final gate. The quality
parameter $p$ determines the average error-rate $r$ according to the relation

$r=\frac{{{2^{n}}-1}}{{{2^{n}}}}(1-p).$ (66)

The RB protocol has been extensively developed in recent years. It has been
extended for other finite groups [391, 392, 393, 394, 395, 259, 396, 397,
398], and it has been also extensively experimentally studied [33, 399, 400,
401, 402]. More recently, a general framework for randomized benchmarking has
been discussed [403].

####  V.1.2 Cross-entropy benchmarking (XEB)

Linear XEB is a sampling-based approach to characterize an arbitrary quantum
gate [404, 405, 406], and is experimentally friendly for large-scale quantum
devices with dozens of qubits. In particular, XEB provides a measure to
approximate the fidelity of random quantum circuit (RQC) in the quantum
computational advantage experiments with limited bitstring samples [7, 8, 9].

For single-qubit gates XEB, the steps are:

Step 1: Generate XEB sequences. First randomly sample a sequence of gates of a
fixed length $m$ from the single-qubit gateset formed by the eight $\pi/2$
rotations along the following axes in the Bloch sphere representation: $\pm
X$, $\pm Y$, and $\pm(X\pm Y)$ on a single-qubit, and then end with a final
random single-qubit gate and measurement. For each length $m$, we choose
$K_{m}$ XEB sequences.

Step 2: Execute the XEB sequences on the noisy quantum devices. By comparing
the measured state probabilities with the ideal ones, we can calculate the
relative cross-entropy differences between probability distributions [405],

$\alpha{\rm{=}}\frac{{{H_{\text{uniform},\text{ideal}}}-{H_{\text{experiment},\text{ideal}}}}}{{{H_{\text{uniform},\text{ideal}}}-{H_{\text{ideal}}}}}$
(67)

where $H=-\sum\nolimits_{i}{{p_{i}}\log({q_{i}})}$ is the cross entropy
between two probability distributions $\\{p_{i}\\}$ and $\\{q_{i}\\}$, and
$H_{\text{uniform},\text{ideal}}$ is the cross-entropy between the uniform and
ideal distributions, $H_{\text{experiment},\text{ideal}}$ is the cross-entropy
between the experimentally measured and ideal distributions, and
$H_{\text{ideal}}$ is the self-entropy of the ideal distribution. Repeat
$K_{m}$ sequences to obtain the a single average survival probability
$\overline{\alpha_{m}}$.

Step 3: Fit the results. Repeat Steps 2 for various sequence lengths $m$
produces a list of average survival probabilities
$\\{\overline{\alpha_{m}}\\}_{m}$, and then fit
$\\{\overline{\alpha_{m}}\\}_{m}$ to the a single exponential decay model:

$\overline{\alpha_{m}}\ =A_{0}p^{m}+B_{0}$ (68)

where $A_{0}$ and $B_{0}$ absorb state preparation and measurement errors. The
decay parameter $p$ can be converted to the average error $r$ and Pauli error
$r_{P}$,

$r=\frac{{N-1}}{N}(1-p)$ (69)

${r_{P}}=\frac{{N-1}}{N}(r)$ (70)

with $N=2^{n}$ the dimensionality of the system having $n$ qubits.

The procedure of characterizing two-qubit gates using XEB is similar to that
of single-qubit gates, except that the circuit sequence generated in Step 1 is
different. Specifically, each cycle of the XEB circuit consisting of a layer
of random single-qubit gates and the two-qubit gate, with a final round of
random single-qubit gates appended at the end, where single-qubit gate set is
still $\pm X$, $\pm Y$, and $\pm(X\pm Y)$. Finally, divide the Pauli fidelity
of the cycle by the Pauli fidelity of the two single-qubit gates to obtain the
Pauli fidelity of the two two-qubit gates.

XEB is a powerful calibration tool for near-term quantum devices. Some works
have attempted to analyze how to spoof linear XEB, an important topic related
to the quantum computational advantage experiments [407, 408]. In addition,
Chen et al. proposed the Clifford XEB [409], which replaces the random
circuits with Clifford circuits, to enable the characterization of large-scale
quantum systems, since Clifford circuits can be simulated in polynomial time.

###  V.2 Circuit-level benchmarking

####  V.2.1 Random Quantum Circuit (RQC) Sampling

RQC sampling is outstanding candidate to demonstrate quantum computational
advantages, and also a tool to characterize the overall performance of the
quantum processor. The steps of RQC sampling are:

Step 1: Generate RQC. Each RQC is composed of $m$ cycles, and each cycle is
composed of a single-qubit gate layer and a two-qubit gate layer. In the
single-qubit gate layer, single-qubit gates are applied on all qubits and
chosen randomly from the set of $\\{\sqrt{X},\sqrt{Y},\sqrt{W}\\}$, where
$\sqrt{X}{\rm{=}}{R_{X}}(\pi/2)$, $\sqrt{Y}{\rm{=}}{R_{Y}}(\pi/2)$, and