Randomized benchmarking with gate-dependent noise

TL;DR

This paper provides a rigorous analysis of randomized benchmarking under gate-dependent noise, showing it still exhibits exponential decay and accurately reflects average fidelity, even with strong gate dependence.

Contribution

It proves that randomized benchmarking's exponential decay behavior persists under arbitrary gate-dependent noise, clarifying its operational meaning and validating its empirical effectiveness.

Findings

Decay parameter quantifies average fidelity of gate-dependent noise

Exponential decay converges to that of gate-independent noise

Analysis validated with numerical simulations of strong gate-dependent noise

Abstract

We analyze randomized benchmarking for arbitrary gate-dependent noise and prove that the exact impact of gate-dependent noise can be described by a single perturbation term that decays exponentially with the sequence length. That is, the exact behavior of randomized benchmarking under general gate-dependent noise converges exponentially to a true exponential decay of exactly the same form as that predicted by previous analysis for gate-independent noise. Moreover, we show that the operational meaning of the decay parameter for gate-dependent noise is essentially unchanged, that is, we show that it quantifies the average fidelity of the noise between ideal gates. We numerically demonstrate that our analysis is valid for strongly gate-dependent noise models. We also show why alternative analyses do not provide a rigorous justification for the empirical success of randomized benchmarking with gate-dependent noise.

Figures

Equations

$$f(\tilde{\mathcal{G}},\mathcal{G})=\int d\psi \mathrm{Tr}[\mathcal{G}(\psi )\tilde{\mathcal{G}}(\psi )].$$

$$Ap(\mathcal{E}{)}^m+B+{ϵ}_m$$

$$p(\mathcal{E})$$

$$|{ϵ}_m|$$

$$x_{b:a}$$

$$\Vert \mathcal{L}{\Vert }_{1\to 1}=\underset{M:\Vert M{\Vert }_1=1}{\sup }\Vert \mathcal{L}(M){\Vert }_1,$$

$$|A⟩⟩=\sum _j\mathrm{Tr}(B_j^{†}A)e_j$$

$$\mathcal{C}=\sum _j|\mathcal{C}(B_j)⟩⟩⟨⟨B_j|,$$

$${\mathbb{E}}_{\vec{G}\in {\mathbb{G}}^m}(Q_{\vec{G}})=A{p}^m+B$$

$${\mathcal{D}}_p(\rho )=p\rho +\frac{1-p}{d}{I}_d.$$

$${\mathcal{D}}_{p,t}(\rho )=\frac{t}{d}{I}_d+p(\rho -\frac{1}{d}{I}_d),$$

$${\mathbb{E}}_{G\in \mathbb{G}}({\mathcal{G}}^{†}\mathcal{C}\mathcal{G})={\mathcal{D}}_{p,t},$$

$$t(\mathcal{C})$$

$$p(\mathcal{C})$$

$${\mathbb{E}}_{G\in \mathbb{G}}(\tilde{\mathcal{G}}\mathcal{L}{\mathcal{G}}^{†})$$

$${\mathbb{E}}_{G\in \mathbb{G}}({\mathcal{G}}^{†}\mathcal{R}\tilde{\mathcal{G}})$$

$${\mathbb{E}}_{G\in \mathbb{G}}(\mathcal{G}\mathcal{R}\mathcal{L}{\mathcal{G}}^{†})$$

$${\mathbb{E}}_{G\in \mathbb{G}}({\mathcal{G}}^{†}\mathcal{R}{\mathcal{L}}_G\mathcal{G})={\mathcal{D}}_{p,t}.$$

$$p=p({\mathbb{E}}_{G\in \mathbb{G}}[{\mathcal{G}}^{†}\mathcal{R}{\mathcal{L}}_G\mathcal{G}])={\mathbb{E}}_{G\in \mathbb{G}}(p[\mathcal{R}{\mathcal{L}}_G])$$

$$t=t({\mathbb{E}}_{G\in \mathbb{G}}[{\mathcal{G}}^{†}\mathcal{R}{\mathcal{L}}_G\mathcal{G}])={\mathbb{E}}_{G\in \mathbb{G}}(t[\mathcal{R}{\mathcal{L}}_G])$$

$$\underset{\lambda \in \mathrm{spec}(M)}{\min }|\eta -\lambda |\le \Vert \delta M\Vert ,$$

$${\mathbb{E}}_{G\in \mathbb{G}}(\tilde{\mathcal{G}})|L⟩⟩$$

$${\mathbb{E}}_{G\in \mathbb{G}}(\tilde{\mathcal{G}}{\mathcal{L}}^{\prime }{\mathcal{G}}_{\mathrm{u}}^{†})$$

$$⟨⟨R|{\mathbb{E}}_{G\in \mathbb{G}}(\tilde{\mathcal{G}})$$

$${\mathbb{E}}_{G\in \mathbb{G}}({\mathcal{G}}_{\mathrm{u}}^{†}{\mathcal{R}}^{\prime }\tilde{\mathcal{G}})$$

$$\mathrm{vec}(\sum _j{v}_j{e}_j^{†})=\sum _j{e}_j\otimes v_j.$$

$$\mathrm{vec}(ABC)=(C^T\otimes A)\mathrm{vec}(B),$$

$${\mathbb{E}}_{G\in \mathbb{G}}({\mathcal{G}}_{\mathrm{u}}\otimes \tilde{\mathcal{G}})\mathrm{vec}({\mathcal{L}}^{\prime })$$

$${\mathbb{E}}_{G\in \mathbb{G}}(\tilde{\mathcal{G}}\otimes {\mathcal{G}}_{\mathrm{u}}{)}^T\mathrm{vec}({\mathcal{R}}^{\prime })$$

$${\mathbb{E}}_{G\in \mathbb{G}}(\mathcal{G}\mathcal{R}{\mathcal{G}}^{†})={\mathcal{D}}_{p,t}.$$

$$\Vert {\tilde{\mathcal{G}}}^{\mathcal{I}}-\mathcal{G}\mathcal{R}\Vert \le \Vert {\tilde{\mathcal{G}}}^{\mathcal{I}}-\mathcal{G}{\mathcal{D}}_{p,t}\Vert +\frac{\Vert {\mathbb{E}}_{G\in \mathbb{G}}({\mathcal{G}}^{†}{\tilde{\mathcal{G}}}^{\mathcal{I}})-{\mathcal{D}}_{p,t}\Vert }{1-{\mathbb{E}}_{G\in \mathbb{G}}(\Vert {\tilde{\mathcal{G}}}^{\mathcal{I}}-\mathcal{G}{\mathcal{D}}_{p,t}\Vert )\Vert {\mathcal{D}}_{1/p,1/t}\Vert }.$$

$$\Vert {\tilde{\mathcal{G}}}^{\mathcal{I}}-\mathcal{G}\mathcal{R}\Vert$$

$$\le \Vert {\tilde{\mathcal{G}}}^{\mathcal{I}}-\mathcal{G}{\mathcal{D}}_{p,t}\Vert +\Vert {\mathcal{D}}_{p,t}-\mathcal{R}\Vert$$

$${\mathcal{D}}_{p,t}^2+{\mathcal{D}}_{p,t}{\mathcal{R}}_1$$

$$={\mathcal{D}}_{p,t}^2+{\mathbb{E}}_{G\in \mathbb{G}}({\mathcal{G}}^{†}\mathcal{R}{\tilde{\mathcal{G}}}_1)$$

$$\Vert {\mathcal{R}}_1\Vert$$

Full text

1 Introduction

The development of practical, large-scale devices that process quantum information relies upon techniques for efficiently characterizing the quality of control operations on the quantum level. Fully characterizing quantum processes [1, 2] requires resources that scale exponentially with the number of qubits despite improvements such as compressed sensing and direct fidelity estimation [3, 4, 5, 6, 7]. However, some figures of merit can be efficiently estimated without fully characterizing quantum processes. The canonical such figure is the fidelity of an experimental implementation $\tilde{\mathcal{G}}$ of an ideal unitary channel $\mathcal{G}(\rho)=G\rho G^{\dagger}$ averaged over input states distributed according to the Haar measure $d\psi$,

[TABLE]

At present, the only efficient partial characterizations are randomized benchmarking (RB) [8, 9, 10, 11, 12] and variants thereof [13, 14, 15, 16, 17, 18, 19]. The RB protocol of Ref. [12] provides an estimate of the average fidelity over a group of operations from how the fidelity of sequences of operations decays as more operations are applied. RB protocols scale efficiently with the number of qubits and are robust to state-preparation and measurement (SPAM) errors. Consequently, RB protocols have become an important baseline for the validation and verification of quantum operations [20, 21] and have recently been used to efficiently optimize over experimental control procedures [22].

However, rigorous analyses of RB protocols make several key assumptions, namely, that the noise is Markovian and time- and gate-independent, so that RB decay curves provide somewhat heuristic estimates of noise parameters (as do all noise characterization protocols). Numerical simulations of RB experiments have shown that the standard single-exponential decay curve holds (approximately) for a variety of physically relevant noise models with gate-dependent and non-Markovian noise, although the ‘true’ and estimated fidelity often differ by a factor of approximately two [23, 25]. However, it has recently been noted that gate-dependent noise fluctuations can dominate the gate-independent decay curve, resulting in a different decay parameter than expected when compared to a standard definition of the average noise [26].

In this paper, we prove that, under the assumption of Markovian noise, the full effect of gate-dependent fluctuations can be described by a single perturbation term that decays exponentially with the sequence length $m$. We also prove that the RB decay parameter is the fidelity of a suitably defined average noise process $\mathcal{E}$. Informally, our main result (4) is that for gate-dependent, but trace-perserving and Markovian noise, the average survival probability over all RB sequences of length $m$ is

[TABLE]

for some constants $A$ and $B$, where

[TABLE]

$\delta_{1},\delta_{2}$ quantify the gate-dependence of the noise and $\delta_{2}$ is small for good implementations of gates by 3. We prove all theorems without restricting to trace-preserving noise, as we anticipate that near-term implementations of RB on encoded qubits [27] will involve significant loss from the encoded space due to imperfect error-correction procedures. We also demonstrate that the alternative analyses of Ref. [12, 28, 24, 26] do not rigorously justify fitting RB decays to a single exponential in section 7.

This paper is structured as follows. We introduce notation in section 2 and review the RB protocol of Ref. [12] in section 3. We then discuss an ambiguity in the decomposition of noise in section 4 and show how it can be exploited to cancel out the majority of gate-dependent fluctuations in section 5. We illustrate the accuracy of our new analysis numerically in section 6 and show why alternative analyses of RB for gate-dependent noise do not justify fitting RB experiments to single exponential decays in section 7. We conclude by discussing some implications of our results for implementing and interpreting RB.

2 Notation

We use the following notation throughout this paper.

• •

All operators except density operators are denoted by Roman font (e.g., $A$). The normalized identity matrix is $\hat{I}_{d}=I_{d}/\sqrt{d}$.

• •

Ideal channels are denoted by calligraphic font (e.g., $\mathcal{C}$), where the ideal channel $\mathcal{U}$ corresponding to a unitary operator $U$ is $\mathcal{U}(\rho)=U\rho U^{\dagger}$.

• •

The unital component of a trace-preserving map is the linear map defined by $\mathcal{G}_{\mathrm{u}}(A)=\mathcal{G}(A-\operatorname{Tr}(A)/d)$.

• •

A noisy implementation of an ideal channel is denoted with an overset $\tilde{}$ (e.g., $\tilde{\mathcal{C}}$ denotes the noisy implementation of $\mathcal{C}$).

• •

Channel composition is denoted by (noncommutative) multiplication (i.e., $\mathcal{C}\mathcal{B}(A)=\mathcal{C}[\mathcal{B}(A)]$).

• •

Noncommutative products of subscripted objects are denoted by the shorthand

[TABLE]

• •

Groups and sets will be denoted by blackboard font, e.g., $\mathbb{G}$.

• •

The uniform average over a set $\mathbb{X}$ is denoted by $\mathbb{E}_{x\in\mathbb{X}}[f(x)]=|X|^{-1}\sum_{x\in\mathbb{X}}f(x)$.

• •

The trace and operator norms of a matrix are $\lVert M\rVert_{1}=\operatorname{Tr}\sqrt{M^{\dagger}M}$ and $\lVert M\rVert=\max_{v}\lVert Mv\rVert/\lVert v\rVert$ respectively.

• •

The induced trace norm of a linear map $\mathcal{L}$ is

[TABLE]

where the maximization is over positive-semidefinite matrices.

We use the following matrix representation of quantum channels (variously known as the Liouville, natural and Pauli transfer matrix representation) wherever a concrete representation is required. Let $\{e_{1},\ldots,e_{n}\}$ be the canonical orthonormal unit basis of $\mathbb{C}^{n}$ and $\{B_{1},\ldots,B_{d^{2}}\}$ be a trace-orthonormal basis of $\mathbb{C}^{d\times d}$, that is, $\operatorname{Tr}(B_{j}^{\dagger}B_{k})=\delta_{j,k}$. Defining the linear map $|*\rangle\!\rangle:\mathbb{C}^{d\times d}\to\mathbb{C}^{d^{2}}$ by

[TABLE]

and the adjoint map $\langle\!\langle A|=|A\rangle\!\rangle^{\dagger}$, we have $\langle\!\langle A|B\rangle\!\rangle=\operatorname{Tr}(A^{\dagger}B)$. A channel $\mathcal{C}$ maps $\rho$ to $\mathcal{C}(\rho)$, which can be represented by the matrix

[TABLE]

where we abuse notation slightly by using the same notation for the abstract channel and its matrix representation. For the remainder of the paper, we assume the operator basis is Hermitian and that $B_{1}=\hat{I}_{d}$, so that the matrix representation of all Hermiticity-preserving maps, including all completely positive and trace-preserving (CPTP) maps, are real-valued.

3 Randomized benchmarking protocol

The RB protocol for a group of operations $\mathbb{G}$ that is also a unitary two-design [11] (e.g., the $n$-qubit Clifford group) is as follows.

Choose a positive integer $m$. 2. 2.

Choose a sequence $\vec{G}=(G_{1},\ldots,G_{m})\in\mathbb{G}^{m}$ of $m$ elements of $\mathbb{G}$ uniformly at random and set $G_{m+1}=(G_{m:1})^{\dagger}$. 3. 3.

Estimate the expectation value (also known as the survival probability) $Q_{\vec{G}}=\operatorname{Tr}[Q\tilde{\mathcal{G}}_{m+1:1}(\rho)]$ of an observable $Q$ after preparing a $d$-level system in a state $\rho$ and applying the sequence of $m+1$ operations $G_{1},G_{2},\ldots,G_{m+1}$. 4. 4.

Repeat steps 2–3 to estimate the average over random sequences $\mathbb{E}_{\vec{G}\in\mathbb{G}^{m}}(Q_{\vec{G}})$ to a desired precision. 5. 5.

Repeat steps 1–4 for different values of $m$ and fit to the decay model

[TABLE]

to estimate $p$ [12].

Despite the ubiquitous usage of RB, best practice choices for the sequence lengths and the number of sequences at each length are still unknown, though see Refs. [23, 30, 31] for some guidance.

A unitary two-design is any set of unitary channels $\mathbb{G}$ such that uniformly sampling $\mathbb{G}$ reproduces the first and second moments of the Haar measure. The canonical example of a unitary two-design is the $n$-qubit Clifford group. The unitary two-design condition can also be stated in terms of how a general channel is ‘twirled’ by the group, where the twirl of a channel $\mathcal{C}$ over a group $\mathbb{G}$ is $\mathbb{E}_{G\in\mathbb{G}}(\mathcal{G}^{\dagger}\mathcal{C}\mathcal{G})$. A unitary two-design is any group $\mathbb{G}$ such that any CPTP channel is ‘twirled’ into a depolarizing channel [11, 32]

[TABLE]

Unitary two-designs also twirl a general linear map $\mathcal{C}$ into a depolarizing map composed with a channel that uniformly decreases the trace, that is, into a channel of the form

[TABLE]

where the values of $p$ and $t$ can be computed from the following lemma [30]. Note that $p$ and $t$ are functions from linear maps to the real numbers, however, we will often omit the argument when it is clear from the context.

Lemma 1**.**

For any unitary two-design $\mathbb{G}$ and channel $\mathcal{C}$,

[TABLE]

where

[TABLE]

In 4, we prove that eq. 8 holds, up to a small and exponentially decaying perturbation, for gate-dependent noise that is time-independent and Markovian, generalizing the original analysis to more realistic noise.

4 Representing gate-dependent noise

Let $\tilde{\mathcal{G}}$ be an experimental implementation of an ideal unitary channel $\mathcal{G}(\rho)=G\rho G^{\dagger}$ with Markovian noise. It is convenient to decompose $\tilde{\mathcal{G}}$ into an ideal gate and a fictitious noise process, for example, $\tilde{\mathcal{G}}=\mathcal{G}\mathcal{E}$. However, we can also write $\tilde{\mathcal{G}}=\mathcal{L}_{G}\mathcal{G}\mathcal{R}_{G}$ for any suitable $\mathcal{L}_{G}$ and $\mathcal{R}_{G}$. We now identify an approximate decomposition that will greatly simplify the analysis of RB with gate-dependent noise.

Theorem 2**.**

Let $\mathbb{G}$ be a unitary two-design, $\{\tilde{\mathcal{G}}:G\in\mathbb{G}\}$ be a corresponding set of Hermiticity-preserving maps, and $p$ and $t$ be the largest eigenvalues of $\mathbb{E}_{G\in\mathbb{G}}(\mathcal{G}_{\mathrm{u}}\otimes\tilde{\mathcal{G}})$ and $\mathbb{E}_{G\in\mathbb{G}}(\tilde{\mathcal{G}})$ in absolute value respectively, where $\mathcal{G}_{\mathrm{u}}(A)=\mathcal{G}(A-\operatorname{Tr}(A)/d)$. There exist linear maps $\mathcal{L}$ and $\mathcal{R}$ such that

[TABLE]

Remark*.*

When the noise is independent of the gate, that is, $\tilde{\mathcal{G}}=\mathcal{L}\mathcal{G}\mathcal{R}$ for all $G\in\mathbb{G}$ for some $\mathcal{L}$ and $\mathcal{R}$, then $\mathcal{L}$ and $\mathcal{R}$ are solutions to eq. 13, as can be seen by substituting $\tilde{\mathcal{G}}=\mathcal{L}\mathcal{G}\mathcal{R}$ into eq. 13 and using 1.

For gate-dependent noise such that $\mathcal{R}$ is invertible, we can set $\tilde{\mathcal{G}}=\mathcal{L}_{G}\mathcal{G}\mathcal{R}$, so that the noise between two ideal gates $\mathcal{G}$ and $\mathcal{H}$ is $\mathcal{R}\mathcal{L}_{G}$. Substituting $\tilde{\mathcal{G}}=\mathcal{L}_{G}\mathcal{G}\mathcal{R}$ into eq. 13b and multiplying by $\mathcal{R}^{-1}$ gives

[TABLE]

As the Haar measure is unitarily invariant and $p$ and $t$ are convex functions,

[TABLE]

give the fidelity and loss of trace of the average noise between ideal gates to the identity via 1. If $\mathcal{R}$ is not invertible, we can introduce an arbitrarily small perturbation into the $\tilde{\mathcal{G}}$ to make $\mathcal{R}$ invertible and so Remark will hold to an arbitrary precision.

When $\tilde{\mathcal{G}}=\mathcal{G}\mathcal{D}_{p,t}$, $\mathbb{E}_{G\in\mathbb{G}}(\mathcal{G}_{\mathrm{u}}\otimes\tilde{\mathcal{G}})$ and $\mathbb{E}_{G\in\mathbb{G}}(\tilde{\mathcal{G}})$ each have only one nonzero eigenvalue, namely, $p$ and $t$ respectively [15, 16]. By the Bauer-Fike theorem [33], all the eigenvalues $\eta$ of a matrix $M+\delta M$ satisfy

[TABLE]

where $\operatorname{spec}(M)$ is the spectrum of $M$. As the operator norm is submultiplicative, obeys the triangle inequality, and $\lVert A\otimes B\rVert=\lVert A\rVert\lVert B\rVert$, the largest eigenvalues of $\mathbb{E}_{G\in\mathbb{G}}(\mathcal{G}_{\mathrm{u}}\otimes\tilde{\mathcal{G}})$ and $\mathbb{E}_{G\in\mathbb{G}}(\tilde{\mathcal{G}})$ will differ from $p^{\prime}$ and $t^{\prime}$ by at most $\mathbb{E}_{G\in\mathbb{G}}\lVert\tilde{\mathcal{G}}-\mathcal{G}\mathcal{D}_{p^{\prime},t^{\prime}}\rVert$ for any scalars $p^{\prime}$ and $t^{\prime}$. Therefore the largest eigenvalues in absolute value will generally be unique and hence must be real as they are the eigenvalues of real matrices.

Proof.

Equations 13a and 13b are essentially a pair of eigenvalue equations, except that $\mathbb{G}$ acts irreducibly on $|\hat{I}_{d}\rangle\!\rangle$ and its orthogonal complement. To utilize this structure, let $\mathcal{L}=|L\rangle\!\rangle\!\langle\!\langle\hat{I}_{d}|+\mathcal{L}^{\prime}$ and $\mathcal{R}=|\hat{I}_{d}\rangle\!\rangle\!\langle\!\langle R|+\mathcal{R}^{\prime}$ with $\mathcal{L}^{\prime}|\hat{I}_{d}\rangle\!\rangle=0$ and $\langle\!\langle\hat{I}_{d}|\mathcal{R}^{\prime}=0$ without loss of generality. With these definitions and as $\mathcal{G}|\hat{I}_{d}\rangle\!\rangle=|\hat{I}_{d}\rangle\!\rangle$ for any unitary (or unital) channel, we can separate each equation in eqs. 13a and 13b into two independent equations

[TABLE]

where restricting to $\mathcal{G}_{\mathrm{u}}$ enforces the orthogonality constraints on $\mathcal{L}^{\prime}|\hat{I}_{d}\rangle\!\rangle=0$ and $\langle\!\langle\hat{I}_{d}|\mathcal{R}^{\prime}=0$. Equations 17c and 17a are left- and right-eigenvector problems of a real matrix and so $t$ can be set to any eigenvalue of $\mathbb{E}_{G\in\mathbb{G}}(\tilde{\mathcal{G}})$ and $L$ and $R$ to the corresponding non-trivial eigenvectors (recall that the set of left- and right-eigenvalues of a real matrix are identical). We choose the largest eigenvalue as it will result in the smallest $\Delta_{G}=\tilde{\mathcal{G}}-\mathcal{L}\mathcal{G}\mathcal{R}$.

To turn eqs. 17d and 17b into eigenvalue equations, we can use the vectorization map that maps a matrix to a vector by stacking the columns vertically, that is,

[TABLE]

The vectorization map satisfies the identity

[TABLE]

Applying eq. 19 to eqs. 17d and 17b and using $\mathcal{G}^{T}=\mathcal{G}^{\dagger}$ (as the matrix basis is Hermitian) gives

[TABLE]

and so $p$ can be set to any eigenvalue of $\mathbb{E}_{G\in\mathbb{G}}(\mathcal{G}_{\mathrm{u}}\otimes\tilde{\mathcal{G}})$ and $\mathcal{L}^{\prime}$ and $\mathcal{R}^{\prime}$ to the corresponding non-trivial eigenvectors (or equivalently, to any eigenvalue of $\mathbb{E}_{G\in\mathbb{G}}(\tilde{\mathcal{G}}\otimes\mathcal{G}_{\mathrm{u}})$ as the two are related by a unitary transformation and so have the same spectrum).

The maps $\mathcal{L}$ and $\mathcal{R}$ identified are solutions to eigenvector equations and so are only determined up to a normalization constant. As $\mathcal{D}_{q,u}$ commutes with $\mathcal{G}$ for all scalars $q$ and $v$ and all $\mathcal{G}\in\mathbb{G}$, we can set $\mathcal{L}\to\mathcal{L}D_{q,u}$ and $\mathcal{R}\to\mathcal{D}_{r,v}\mathcal{R}$ to satisfy eq. 13c. ∎

While applying a gauge transformation $\tilde{\mathcal{G}}\to\mathcal{S}^{-1}\tilde{\mathcal{G}}\mathcal{S}$ changes $\mathcal{L}$ and $\mathcal{R}$ to $\mathcal{L}\to\mathcal{S}\mathcal{L}$ and $\mathcal{R}\to\mathcal{R}\mathcal{S}^{-1}$ respectively, the noise between ideal gates, namely, $\mathcal{R}\mathcal{L}$, is gauge invariant and so the present analysis is gauge-invariant unlike the original analysis of RB as shown in Ref. [26].

We now prove that the solutions $\mathcal{L}$ and $\mathcal{R}$ to eq. 13 can also be chosen to ensure that $\Delta_{G}=\tilde{\mathcal{G}}-\mathcal{L}\mathcal{G}\mathcal{R}$ is small when $\tilde{\mathcal{G}}\approx\mathcal{G}\mathcal{D}_{p,t}$. To prove this, it is convenient to transform to the gauge $\tilde{\mathcal{G}}^{\mathcal{I}}$ such that $\mathcal{L}=\mathcal{I}$ (i.e., set $\mathcal{S}=\mathcal{L}^{-1}$). In this gauge, eq. 13c becomes

[TABLE]

In particular, $\Delta_{G}$ is on the order of $\lVert\tilde{\mathcal{G}}-\mathcal{G}\mathcal{D}_{p,t}\rVert\in O(\sqrt{1-p})$.

Theorem 3**.**

Let $\mathbb{G}$ be a unitary two-design, $\{\tilde{\mathcal{G}}:G\in\mathbb{G}\}$ be a corresponding set of Hermiticity-preserving maps, and $p$ and $t$ be the largest eigenvalues of $\mathbb{E}_{G\in\mathbb{G}}(\mathcal{G}_{\mathrm{u}}\otimes\tilde{\mathcal{G}})$ and $\mathbb{E}_{G\in\mathbb{G}}(\tilde{\mathcal{G}})$ in absolute value respectively, where $\mathcal{G}_{\mathrm{u}}(A)=\mathcal{G}(A-\operatorname{Tr}(A)/d)$. There exists a linear map $\mathcal{R}$ satisfying eq. 13 with $\mathcal{L}=\mathcal{I}$ such that

[TABLE]

Proof.

By the triangle inequality, submultiplicativity, and unitary invariance of the operator norm,

[TABLE]

for all $G\in\mathbb{G}$.

Section 4 holds for any $\mathcal{R}$. Now let $\mathcal{R}$ be a solution to eq. 13 for $\tilde{\mathcal{G}}^{\mathcal{I}}$ and expand $\mathcal{R}=\mathcal{D}_{p,t}+\mathcal{R}_{1}$ and $\tilde{\mathcal{G}}^{\mathcal{I}}=\mathcal{G}\mathcal{D}_{p,t}+\tilde{\mathcal{G}}_{1}$. Substituting these expansions into eq. 13b and using eq. 21 gives

[TABLE]

Canceling the common term in section 4, multiplying both sides by $\mathcal{D}_{1/p,1/t}$ from the left, taking the operator norm and using the triangle inequality and submultiplicativity of the operator norm gives

[TABLE]

Rearranging and using the unitary invariance of the operator norm gives

[TABLE]

∎

5 Analyzing RB with arbitrarily gate-dependent noise

Experimental implementations of RB almost always produce relatively clean exponential decays. We now explain this empirical success for general gate-dependent noise by proving that the average survival probability over RB sequences of length $m$ is equal to the model of Ref. [28] for gate-independent noise with a single perturbation term for a suitably defined average noise model. Moreover, the decay constants correspond to the average loss of trace and fidelity of the noise between ideal gates via Remark. The perturbation term will generally be negligible for noise that is close to the identity by 3 and for $m\geq 3$ (for numerical values for a single qubit, see fig. 1). Even under this condition, more significant (and potentially oscillating) perturbations are possible for $m\approx 1$, potentially explaining experimentally observed fluctuations as seen in, e.g., [20, Fig. 3b]. While we have not investigated how the perturbation terms scale with dimension, note that we only need $\delta_{2}\lessapprox 1/2$ in order for the decay to be essentially a single exponential for $m\geq 10$, etc, as $\delta_{2}$ is exponentially suppressed.

Theorem 4**.**

For any unitary two-design $\mathbb{G}$ and set of Hermiticitiy preserving linear maps $\{\tilde{\mathcal{G}}:G\in\mathbb{G}\}$, the average survival probability over all randomized benchmarking sequences of length $m$ is

[TABLE]

for some constants $A$ and $B$, where $p$ and $t$ give the fidelity and loss of trace of the average noise between ideal gates via Remark. Moreover, the perturbation term $\epsilon_{m}$ satisfies

[TABLE]

for some $\delta_{1}$ and $\delta_{2}$ that quantify the amount of gate dependence.

Remark*.*

Equation 27 reduces to the standard single-exponential model plus an exponentially decreasing perturbation when $\tilde{\mathcal{G}}$ is a CPTP map for all $G$ as then $\mathbb{E}_{G\in\mathbb{G}}(\tilde{\mathcal{G}}\otimes\mathcal{G})$ and $\mathbb{E}_{G\in\mathbb{G}}(\mathcal{G}\otimes\tilde{\mathcal{G}})$ are both CPTP maps and so all eigenvalues are in the unit disc and the eigenvalue with eigenvector close to $|\hat{I}_{d}\rangle\!\rangle$ is 1 [34].

The constants $A$ and $B$ and the perturbation term are

[TABLE]

with $\Delta_{j}=\Delta_{G_{j}}=\tilde{\mathcal{G}}_{j}-\mathcal{L}\mathcal{G}_{j}\mathcal{R}$ and where $\mathcal{L}$ and $\mathcal{R}$ are solutions to eq. 13. The quantities bounding the perturbation term are

[TABLE]

which are in turn bounded in 3 in the gauge where $\mathcal{L}=\mathcal{I}$. Note that $\lVert\rho\rVert_{1},\lVert Q\rVert_{1}\leq 1$ in a standard gauge, however, the gauge transformation required to set $\mathcal{L}=\mathcal{I}$ may change these values and introduce negative eigenvalues to $\rho$ so that the maximization for the induced norm in $\delta_{1}$ is not restricted to positive semi-definite inputs.

Proof.

The average map over all randomized benchmarking sequences $\vec{G}$ of length $m$ is

[TABLE]

Let $\mathcal{L}$ and $\mathcal{R}$ satisfy

[TABLE]

and define

[TABLE]

Then for any integer $1\leq j\leq m+1$ we can write

[TABLE]

where $\Delta_{0:1}=I_{d}=\tilde{\mathcal{G}}_{m+1:m+2}$ by convention. Moreover, as $G_{m+1:1}=I_{d}$ and any set of $m$ elements of $\vec{G}$ are statistically independent and uniformly distributed, we can apply eq. 32a recursively and use the fact that $\mathcal{D}_{p,t}$ commutes with unital channels to obtain

[TABLE]

for $j=1$. Similarly, we can use 1 and 32 for $j\geq 2$ to obtain

[TABLE]

Therefore we can apply sections 5 and 5 and then recursively apply sections 5 and 5 to obtain

[TABLE]

The average survival probability over all random sequences is then

[TABLE]

where

[TABLE]

We now bound the perturbation term. By the triangle inequality and submultiplicativity,

[TABLE]

∎

6 Numerics

In fig. 1 we illustrate the reliability of our analysis by comparing the estimates obtained by fitting the simulated RB experiments described below under independent random unitary noise on a set of generating gates to the decay parameters predicted by 4 and by plotting the size of the perturbation terms $\delta_{1}$ and $\delta_{2}$ calculated in the gauge where $\mathcal{L}=\mathcal{I}$. Mathematica code is available at https://github.com/jjwallman/numerics.

For simplicity, we simulate RB experiments using the (projective) group

[TABLE]

where $I,X,Y,Z$ are the standard single-qubit Pauli matrices and

[TABLE]

which is a subgroup of the single-qubit Clifford group that is also a unitary two-design. For our simulations, the noisy implementation of $G=T^{t}P\in\mathbb{G}$ is $\tilde{\mathcal{G}}=\mathcal{G}^{\prime}$ where

[TABLE]

where $U$ and $V$ are independently sampled from the set $\mathbb{U}_{r}\subset\mathrm{SU}(d)$ satisfying $1-f(\mathcal{U},\mathcal{I})=r$, that is, the set of unitary channels of fixed infidelity $r$. We then calculate $p(\mathcal{R}\mathcal{L})$ from 2 and $\delta_{1},\delta_{2}$ from Remark. We set the state to $|0\rangle$ and the measurement to be a computational basis measurement and do not include SPAM errors for simplicity.

The set $\mathbb{U}_{r}$ can be sampled by sampling $V\in\mathrm{SU}(d)$ uniformly from the Haar measure and setting

[TABLE]

which gives the unique unitarily invariant distribution over $\mathbb{U}_{r}$ for some $r$ as a function of $\theta$. In order for the infidelity of the resulting gates to be $r$, we require [32]

[TABLE]

or $\theta=\arcsin(\sqrt{3r/2})$.

A single RB experiment for a fixed value of $m$ then consisted of the following.

Randomly choose $\vec{a}\in\{0,1,2\}^{m}$ and $\vec{b}\in\{1,2,3,4\}^{m}$. 2. 2.

Set $G_{j}=T^{a_{j}}P_{b_{j}}P_{b_{j-1}}T^{-a_{j-1}}=T^{a_{j}-a_{j-1}}P_{\tau_{j}}$ where $(P_{0},P_{1},P_{2},P_{3})=(I,X,Y,Z)$, $a_{m+1}=a_{0}$, $b_{m+1}=b_{0}=1$, and $P_{\tau_{j}}=T^{a_{j-1}}P_{b_{j}}P_{b_{j-1}}T^{-a_{j-1}}$. 3. 3.

Calculate $\lvert\langle 0\rvert G^{\prime}_{m+1:1}\lvert 0\rangle\rvert^{2}$

We repeat the above steps for 100 random sequences for each value of $m\in\{4,8,16,...,2048\}$ and fit the averages for each sequence to the model $\mathbb{E}_{\vec{\mathcal{C}}}=0.5p_{est}^{m}+0.5$ using Mathematica’s NonlinearModelFit function weighted by the variance over the 100 sequences at each $m$, where we have set $A=B=0.5$ because of the choice of SPAM and because the noise is unital.

7 Alternative analysis of RB with gate-dependent noise

We now show how the analyses of gate-dependent noise of Refs. [12, 28, 24, 26] are too loose to justify fitting RB decay curves to a single exponential in practical regimes. We will first show that the higher-order terms in the original analysis of RB give a systematic uncertainty of the infidelity that dominates the infidelity, even when the fidelity is dominated by stochastic errors. We will then discuss subsequent analysis and how it also results in large systematic errors.

7.1 Initial RB Analysis

The initial method of analyzing RB under gate-dependent noise was to perform a perturbation expansion to perturb about a perfect twirl [12, 28] of the average noise,

[TABLE]

Substituting $\tilde{\mathcal{G}}=\mathcal{G}\mathcal{E}+\mathcal{G}\Delta_{G}$ into the average map over all sequences,

[TABLE]

and using the triangle inequality and submultiplicativity, the $k$th order terms in $\Delta_{G}$ contribute at most

[TABLE]

to the average survival probability over all sequences of length $m$, $\langle\!\langle Q\rvert\mathbb{E}_{\vec{G}\in\mathbb{G}^{m}}(\tilde{\mathcal{G}}_{m+1:1})\lvert\rho\rangle\!\rangle$. Ref. [28] also obtained similar conditions for time-dependent noise and derived the first-order correction. As stated in Ref. [12], higher-order terms will be negligible provided $m\gamma\ll 1$, so that there is a tension between the amount of gate-dependent noise and the values of $m$ required to obtain a reasonable fit. However, there exist noise models such that $m\gamma$ is not negligible [26].

The analysis of Ref. [28] also holds for the more general gate-dependent noise model $\mathcal{L}_{G}\mathcal{G}\mathcal{R}_{G}$, where $\mathcal{L}$ and $\mathcal{R}_{G}$ are gate-independent and gate-dependent noise processes acting from the left and right (that is, after and before a ideal gate) respectively, so that some optimization of $\gamma$ is possible.

However, we now show that the higher-order terms cannot be neglected when fitting data in experimental regimes. Suppose the average survival probability over all sequences of some fixed length $m$ is measured to be $f_{m}$. From the zeroth-order model with $p=1-2r$,

[TABLE]

ignoring any finite measurement statistics. However, by eq. 48, the first-order perturbations add some $\delta\in[-m\gamma,m\gamma]$ to $f_{m}$, so that the true value of $r$ is

[TABLE]

where we have assumed $m\gg 1$ for the approximation, so that the first-order perturbation terms add a systematic uncertainty of $\gamma/2(f_{m}-B)$ to the estimate of $r$.

Therefore the systematic ucertainty will dominate the estimate of the infidelity unless $\gamma\approx r(\mathcal{E})$. To illustrate how impractical this requirement is, and thus the value of the tighter analysis in section 5, consider the following hypothetical implementation of the single-qubit Clifford group (that is, the 24 elements of SU(2) that permute the single-qubit Paulis $X$, $Y$ and $Z$ up to an overall sign). Supoose that half the elements are implemented with only depolarizing noise, that is, $\tilde{\mathcal{G}}=\mathcal{G}\mathcal{D}_{\nu}$, and the other half are implemented with an additional rotation around the $z$-axis of the Bloch sphere by some angle $\theta\in(0,\tfrac{\pi}{2})$, that is, $\tilde{\mathcal{G}}=\mathcal{G}\mathcal{D}_{\nu}\mathcal{Z}_{\theta}$ where

[TABLE]

For this noise model, the average noise is $\mathcal{E}=\mathcal{D}_{\nu}(I_{4}+\mathcal{Z}_{\theta})/2$ by eq. 46 and so

[TABLE]

using $1-r(\mathcal{E})=f(\mathcal{E},\mathcal{I})=1/2+\sum_{\sigma=X,Y,Z}\operatorname{Tr}[\sigma\mathcal{E}(\sigma)]/12$ [36]. The maximization in the induced trace norm will be achieved by any state in the $xy$ plane (e.g., $(|0\rangle+|1\rangle)/\sqrt{2}$), giving

[TABLE]

For the estimate based on the zeroth-order model to be valid to within a factor of two, we require $\gamma\leq r(\mathcal{E})$, that is, $\lvert\theta\rvert\lesssim 1-\nu$. For example, the systematic uncertainty in the estimate of $r$ is $9r$ (i.e., under-reporting the error rate $r$ by an order of magnitude) when $\nu=0.99$ and gate-dependent errors account for 10% of the total infidelity ($\theta=0.09$ in the above model) or when $\nu=0.999$ and the gate-dependent errors account for just 1% of the total infidelity ($\theta=0.009$ in the above model).

One might expect that the above analysis is pessimistic and that the contributions from gate-dependent perturbations at different sequence lengths would average out. However, Ref. [26] presented an explicit model wherein the contributions from gate-dependent perturbations are systematic.

7.2 Subsequent analyses of RB

Ref. [26] demonstrated that the RB decay rate could differ from the average gate infidelity of the average noise defined in eq. 46 by orders of magnitude. Ref. [26] also presented a significantly improved analysis of the original RB protocol, which can be viewed as an application of the method of [24] to the original RB protocol. However, this analys has three crucial drawbacks.

First, no interpretation of the decay parameter was provided, leading to frequent discussions of “the RB number” $r_{RB}$ at conferences. In the original version, Ref. [26] states that “it is not clear that $r$ (a function of the RB decay parameter) should be called (the) “fidelity”, except in the special case of gate-independent error maps,” although this statement has since been revised based on the analysis in the present paper.

Second, the analysis of Ref. [26] is vulnerable to the criticism it makes of Ref. [28], namely, that the infidelity $r(\mathcal{E})$ of a fixed decomposition of the noise (defined in eq. 46) can differ from $r_{RB}$ by orders of magnitude. This is problematic for the analysis of Ref. [26] (and also for that of Ref. [24]) because they bound a perturbation from an exponential of the form $(1-r_{RB})^{m}$ (ignoring small dimensional factors) by the diamond norm with respect to a fixed decomposition of the noise, giving a perturbation on the order of $\sqrt{r(\mathcal{E})}\gg r(\mathcal{E})\gg r_{RB}$. While some gauge optimization of their bound is possible (especially using a gauge motivated by the present analysis), it is unclear whether such optimization can close the gap between $r(\mathcal{E})$ and $r_{RB}$. Fortunately, the present paper obviates the potential problem by providing a much sharper bound.

Finally, even if gauge optimization allows the perturbation in their method to be bounded by $\sqrt{r_{RB}}$, this still renders their analysis useless in two key regimes, namely, for “large” error rates or in the linear part of the decay.

For “large” error rates where $r_{RB}\approx 0.01$, which accounts for most RB experiments to date, the best bound that can be put on the perturbation without knowing the full experimental noise process is approximately $\sqrt{r_{RB}}=0.1$. With such a large perturbation, any sequence lengths that have an average survival probability over $0.9$ or under $0.6$ are effectively useless assuming SPAM values of $A=B=0.5$ (and the range of values that contain any information decreases with $A$).

In the linear part of the decay, such as in the numerics of Ref. [26], the systematic uncertainty due to the perturbation dominates the decay. For concreteness, consider the linear regime for a single qubit,

[TABLE]

where $p(\mathcal{R}\mathcal{L})$ is the true RB decay parameter from 4. If $p(\mathcal{R}\mathcal{L})$ is estimated via the slope between $m=m_{1}$ and $m=m_{2}$, then the unknown perturbation contributes a systematic uncertainty of approximately $\sqrt{r(\mathcal{E})}/(m_{1}-m_{2})$ to the estimate of $p$, so that essentially no information is gained unless $m_{1}-m_{2}\gg\sqrt{r(\mathcal{E})}$. For the current record value of $1-p\approx 10^{-6}$, this translates to requiring $m_{2}-m_{1}>1000$ assuming $1-p(\mathcal{R}\mathcal{L})\approx r(\mathcal{E})$. However, the assumption that $1-p(\mathcal{R}\mathcal{L})\approx r(\mathcal{E})$ is unjustified as shown by their analysis, so that the required spacing of sequence lengths in order to overcome the systematic uncertainty due to their bound is unknown, rendering their analysis impractical. Consequently, the analysis of Ref. [26] does not explain the empirical success of RB to date or justify fitting to a single exponential without further analysis.

8 Conclusion

We have proven that randomized benchmarking experiments result in a single exponential with a negligible perturbation under gate-dependent noise. Moreover, we proved that the decay parameter is the fidelity of the average noise between idealized gates. More specifically, we can write $\tilde{\mathcal{G}}=\mathcal{L}\mathcal{G}\mathcal{R}_{G}$ where $\mathcal{L}$ is a solution to eq. 13. Then for trace-preserving noise the RB decay parameter is $\mathbb{E}_{G\in\mathbb{G}}[p(\mathcal{R}_{G}\mathcal{L})]$, that is, the average survival probability between ideal gates. As $\tilde{\mathcal{G}}\to\mathcal{S}^{-1}\tilde{\mathcal{G}}\mathcal{S}$ maps $\mathcal{L}\to\mathcal{S}\mathcal{L}$ and $\mathcal{R}\to\mathcal{R}\mathcal{S}^{-1}$ in 2, the decay parameter $p=\mathbb{E}_{G\in\mathbb{G}}(\mathcal{R}_{G}\mathcal{L})$ is manifestly gauge invariant, answering the recent criticism of Ref. [26].

By rigorously reducing the potential impact of gate-dependent errors, the present work enables randomized benchmarking to more reliably diagnose the presence of non-Markovian noise. More precisely, there are two key assumptions that may not hold well enough for the standard exponential decay model to be valid, namely, Markovianity and time-independence. However, the only effect that time-dependent noise can have is to alter the rate at which the survival probability decays [8, 30]. Consequently, any “revivals” in RB experiments, that is, statistically significant increases in the average survival probability as $m$ increases [23], can be attributed directly to the noise being non-Markovian.

An intriguing observation from the present analysis is that gate-dependent fluctuations may produce a significant deviation from the exponential fit for short sequences. This deviation is possible because the first few values of $m$ are the most sensitive to gate-dependent fluctuations, which may be as large as $0.1$ for $m=1$ and $0.01$ for $m=2$. Consequently, these sequence lengths should be omitted when fitting to eq. 8, but can also indicate the presence of gate-dependent noise.

Alternative protocols based on randomized benchmarking have been developed that also assume gate-independent noise [13, 14, 15, 16, 17, 18, 19]. We leave applying the current analysis to these protocols to future work.

One main open problem of the current work is that the solutions $\mathcal{L}$ and $\mathcal{R}$ from 2 are not guaranteed to be CPTP maps, even under rescaling.

The author acknowledges helpful discussions with Robin Blume-Kohout, Arnaud Carignan-Dugas, Joseph Emerson, Steve Flammia, Chris Granade, Richard Kueng, Robin Harper, Jonas Helsen, Thomas Monz, Jiaan Qi and Philipp Schindler. This research was supported by the U.S. Army Research Office through grant W911NF-14-1-0103. This research was undertaken thanks in part to funding from the Canada First Research Excellence Fund.