Character randomized benchmarking for non-multiplicity-free groups with applications to subspace, leakage, and matchgate randomized benchmarking

Jahan Claes Quantum Artificial Intelligence Laboratory (QuAIL), NASA Ames Research Center, Moffett Field, CA 94035, USA USRA (RIACS), Mountain View CA 94043, USA Department of Physics and Institute for Condensed Matter Theory, University of Illinois at Urbana-Champaign, Urbana, IL 61801, USA Eleanor Rieffel Quantum Artificial Intelligence Laboratory (QuAIL), NASA Ames Research Center, Moffett Field, CA 94035, USA Zhihui Wang zhihui.wang@nasa.gov Quantum Artificial Intelligence Laboratory (QuAIL), NASA Ames Research Center, Moffett Field, CA 94035, USA USRA (RIACS), Mountain View CA 94043, USA

(June 21, 2025)

Abstract

Randomized benchmarking (RB) is a powerful method for determining the error rate of experimental quantum gates. Traditional RB, however, is restricted to gatesets, such as the Clifford group, that form a unitary 2-design. The recently introduced character RB can benchmark more general gates using techniques from representation theory; up to now, however, this method has only been applied to “multiplicity-free” groups, a mathematical restriction on these groups. In this paper, we extend the original character RB derivation to explicitly treat non-multiplicity-free groups, and derive several applications. First, we derive a rigorous version of the recently introduced subspace RB, which seeks to characterize a set of one- and two-qubit gates that are symmetric under SWAP. Second, we develop a new leakage RB protocol that applies to more general groups of gates. Finally, we derive a scalable RB protocol for the matchgate group, a group that like the Clifford group is non-universal but becomes universal with the addition of one additional gate. This example provides one of the few examples of a scalable non-Clifford RB protocol. In all three cases, compared to existing theories, our method requires similar resources, but either provides a more accurate estimate of gate fidelity, or applies to a more general group of gates. In conclusion, we discuss the potential, and challenges, of using non-multiplicity-free character RB to develop new classes of scalable RB protocols and methods of characterizing specific gates.

I Introduction

Advances in accurate and scalable methods for characterizing the performance of quantum gates are critical for the realization of large-scale reliable quantum computers. Quantum process tomography can, in theory, completely characterize an unknown quantum channel [1, 2, 3, 4], but requires resources that scale exponentially in the number of qubits [4]. In addition, any tomographic approach will also include the effect of state preparation and measurement (SPAM) errors, which may be of the same order as the gate error that is being characterized.

Randomized benchmarking (RB) [5, 6, 7, 8] provides a method to scalably characterize gates that form a group $G$ with the additional mathematical property of being a “unitary 2-design” [9], most frequently the Clifford group [10, 11, 12]. Rather than completely characterizing a noise channel, RB determines the average fidelity, a standard measure of gate quality that can be related to other common measures such as entanglement and process fidelity [13, 14] and used to bound the gate error rate [15]. RB works by experimentally measuring the overall fidelity of a random circuit as a function of the number of applied gates $U\in G$ and fitting this to an exponential decay. The parameters of the decay then determine the average fidelity of a single gate. Unlike tomographic methods, RB provides an estimate for the average fidelity that is independent of SPAM errors.

Standard RB, however, is limited to groups that form a unitary 2-design and whose elements can be efficiently compiled (i.e. decomposed) into elementary gates. This limitation prevents standard RB from characterizing any set of quantum gates that are large enough to be universal for quantum computation [11, 12], and also prevents standard RB from characterizing smaller subgroups of 2-designs. There are ongoing efforts to extend RB to a larger class of gates. Interleaved RB was proposed to characterize individual Clifford group elements [16] as well as the $T$ -gates needed for universal quantum computation [17], but these methods are specific to the gates considered and only produce bounds on the fidelity. Ref. [18] developed a method to extract the fidelity of the dihedral group on one qubit, which is not a unitary 2-design and includes the $T$ gate, while [19] proposed a method of extending dihedral RB to an arbitrary number of qubits. Refs. [20, 21] extended this work by deriving decay formulas for the fidelity of random circuits of arbitrary groups, but these formulas involved fitting sums of multiple exponentials, and the decay parameters could not be related to the average fidelity. Ref. [22] introduced character RB to address these limitations, providing a method that only requires fitting a single exponential decay and directly predicts the average fidelity. However, this was only explored for “multiplicity-free” groups, a mathematical limitation on the group’s representations (see below).

In this work, we provide a generalization of character RB that applies to groups with multiplicity, which we underpin with rigorous derivations. This rigor enables us to provide conditions under which instantiations of the framework yield practical RB protocols. We illustrate our generalized approach with applications to three distinct situations of practical interest: benchmarking of gates with subspace preserving properties, characterization of leakage, and benchmarking of the matchgate group.

Our main contributions include:

•

We provide a derivation of character RB for non-multiplicity-free groups $G$ . This RB method allows us to directly predict the average fidelity of the gates in $G$ as in [22] but unlike [20, 21]. For non-multiplicity-free groups, our method potentially requires fitting a sum of multiple exponentials rather than a single exponential; however, the number of exponentials is significantly reduced compared to [20, 21].
•

As a primary motivation for this generalization, we improve the recently introduced subspace RB [23] designed to characterize gates that preserve a subspace of the full Hilbert space. Our generalization, and its rigorous derivations, has immediate application to near-term quantum processors, including to benchmarking the gates implemented on the ion-trap quantum processor benchmarked in [23]. Gates that preserve a proper subspace can never form a 2-design, and are never multiplicity-free, necessitating a generalized RB procedure. The original work on subspace RB established decay formulas for the fidelity of certain random circuits but could only give loose bounds on the average fidelity of the gates; our method, in contrast, allows us to directly estimate the average fidelity using a similar number of experiments as the original subspace RB. While we illustrate our approach for the $U_{ZZ}$ gate seen in [23], the method can be applied directly to other gates with the same SWAP symmetry as the $U_{ZZ}$ gate. It also provides grounding for benchmarking gates with other subspace-preserving symmetries, though creativity will be required to determine when and how these gates can be combined with single qubit gates to obtain a group with the properties that yield a practical character RB protocol. The rigorous derivations underlying our approach enables us to provide examples of noise under which the estimated fidelity yielded by [23] deviates substantially from the exact fidelity provided by our method.
•

We present a new protocol for leakage RB [24, 25, 26], a benchmarking protocol designed to characterize qubits that can “leak” into a non-computational section of the Hilbert space. Our approach reduces the assumptions on control in the leakage subspace required by the original leakage RB work [26]. Such control is frequently unrealistic for quantum hardware. Our approach can be applied immediately to determine certain leakage channel error rates in, for example, quantum dot architectures, though further research will need to be done to obtain a leakage RB protocol that enables the determination of more general parameters including the average fidelity on the computational subspace.
•

We introduce a new scalable RB procedure for the matchgate group [27], a class of quantum circuits that, like the Clifford group, is efficiently simulable [27, 28, 29, 30] but is very close to universal [30, 31, 32, 29, 33, 34, 35]. This procedure necessarily requires the full non-multiplicity-free character RB, and represents, along with the dihedral group [19, 22], one of the few non-Clifford groups that can be scalably benchmarked.

Non-multiplicity-free character RB is a general framework for benchmarking groups of quantum gates. It provides a method for characterizing individual gates when the gates can be combined into operations that form a group, as we illustrate in the case of subspace RB. This RB framework also expands the family of groups that can be scalably benchmarked, as we demonstrate with the matchgate group. Scalable benchmarking protocols are necessary to measure gate quality in large quantum processors, especially in the presence of non-local errors such as crosstalk. While we provide one new example of a scalable benchmarking protocol, we expect the framework of non-multiplicity-free character RB will lead researchers to develop further scalable examples. Benchmarking multiple overlapping groups (or subgroups of groups) may allow more accurate error characterization. While it remains an art to find the groups and constructions that yield practical character RB protocols, we expect the grounding that our work provides to support the discovery of practical protocols for various gate sets in a variety of quantum devices in the years to come.

Our paper is organized as follows. Section II provides mathematical background on the Liouville representation and the definition of average fidelity. Section III outlines the full non-multiplicity-free RB protocol, and proves that it correctly estimates the average fidelity of the gates. The next sections consist of applications. Section IV demonstrates how our method can be used to rigorously estimate the fidelity of gate sets that preserve subspaces, such as those studied in [23]. Section V applies our framework to formulate a leakage RB protocol with fewer assumptions than the current state-of-the-art [26]. Section VI reviews the matchgate group, and describes how our method can be used to derive a scalable RB protocol for this group. Each of our applications are accompanied by computer simulations of benchmarking experiments; all our computer simulations can be reproduced in under a day on a standard laptop. We conclude in Section VII with discussion of possible extensions of our work, including some of the challenges. We relegate technical details to appendices, including Appendix A which demonstrates that our method is robust to gate-dependent errors, and Appendix B which provides a self-contained and straightforward proof that generalizations of the Clifford group to qudits for $d$ prime form a unitary $2$ -design, which may be of independent interest.

II Mathematical Preliminaries

In this paper, we use the Liouville representation of quantum channels. In the Liouville representation, given some fixed basis $\{|i\rangle\}$ of our Hilbert space $\mathcal{H}$ , a density matrix $\rho=\sum_{ij}\rho_{ij}|i\rangle\langle j|$ is represented by a column vector $|\rho\rangle\rangle=\sum_{ij}\rho_{ij}|i\rangle\otimes|j\rangle$ , where we use a double-bracket $|\cdot\rangle\rangle$ to distinguish elements of $\mathcal{H}\otimes\mathcal{H}$ from elements of $\mathcal{H}$ . In the case of a pure state $\rho=|\psi\rangle\langle\psi|$ we will also sometimes write $|\psi\rangle\rangle$ in place of $|\rho\rangle\rangle$ . A quantum channel $\Lambda(\rho)=\sum_{i}A_{i}\rho A_{i}^{\dagger}$ is represented by a matrix $\hat{\Lambda}=\sum_{i}A_{i}\otimes A_{i}^{*}$ . In this representation, matrix multiplication corresponds to composition

\widehat{\Lambda_{1}\circ\Lambda_{2}}=\hat{\Lambda}_{1}\hat{\Lambda}_{2},

matrix-vector multiplication corresponds to applying a quantum channel

\hat{\Lambda}|\rho\rangle\rangle=|\Lambda(\rho)\rangle\rangle,

and the inner product of two vectors corresponds to the Hilbert-Schmidt inner product of the corresponding density matrices

\langle\langle\sigma|\rho\rangle\rangle=\operatorname{Tr}(\sigma^{\dagger}\rho).

In particular, if $M$ is a projector into some measurement outcome, the overlap $\langle\langle M|\rho\rangle\rangle$ gives the probability of measuring $M$ from a state $\rho$ . For a more detailed treatment of the Liouville representation, see [36].

Given a unitary group $G$ acting on our Hilbert space $\mathcal{H}$ , the natural action of $U\in G$ on density matrices is given by $U(\rho)=U\rho U^{\dagger}$ . In the Liouville representation, such an operator is represented by $\hat{U}=U\otimes U^{*}$ . The map $\phi:U\mapsto U\otimes U^{*}$ forms a representation [37] of the group $G$ on $\mathcal{H}\otimes\mathcal{H}$ that we will refer to as the natural representation of $G$ . We can also define the $\bm{G}$ -twirl of a quantum channel $\Lambda$ as

\hat{\Lambda}_{G}=\frac{1}{|G|}\sum_{U\in G}\hat{U}^{\dagger}\hat{\Lambda}\hat{U}.

(1)

where $|G|$ is the order of the group. We can also define the $G$ -twirl by compact groups by replacing the discrete average by the integral over the Haar measure. As we will see, $\Lambda_{G}$ has properties similar to the original channel $\Lambda$ , but it has a simpler structure that makes it more tractable to study.

If a noisy implementation of a gate $U$ results in applying the channel $(\Lambda\circ U)$ , we want to characterize how close the noise channel $\Lambda$ is to the identity. We will focus on one common measure of noise, the average fidelity $F_{\Lambda}$ , given by

F_{\Lambda}:=\int d\psi\langle\langle\psi|\hat{\Lambda}|\psi\rangle\rangle.

(2)

Here, $d\psi$ is the unitary-invariant Haar or Fubini-Study measure on $\mathcal{H}$ . The integrand $\langle\langle\psi|\hat{\Lambda}|\psi\rangle\rangle$ is the probability of preserving a state $|\psi\rangle$ after the noise operator $\Lambda$ has been applied. The average fidelity is then simply the average of this probability over all possible input states.

III The generalized character randomized benchmarking procedure

Let $G$ be the unitary group on $\mathcal{H}$ that we wish to benchmark. We will assume $G$ is either finite or compact, so that every unitary representation decomposes into irredicible representations. Let $\phi:G\rightarrow\mathcal{L}(\mathcal{H}\otimes\mathcal{H})$ be the natural representation of $G$ , which decomposes into irreducible representations as $\phi\simeq a_{1}\phi_{1}\oplus\cdots\oplus a_{I}\phi_{I}$ , where $a_{i}\in\mathbb{Z}^{+}$ is the multiplicity of the irrep $\phi_{i}$ . Let $\mathcal{H}\otimes\mathcal{H}\simeq\bigoplus_{i}\mathbbm{C}^{a_{i}}\otimes\mathcal{H}_{i}$ be the corresponding decomposition of Hilbert space, such that each $\phi_{i}$ acts nontrivially only on a single copy of $\mathcal{H}_{i}$ . We will make the standard RB assumption that the gate error $\Lambda$ associated with $U\in G$ is independent of $U$ , although this can be relaxed [38, 39, 40, 22](see Appendix A).

Let $\overline{G}\subseteq G$ be a subgroup of our unitary group with natural representation $\overline{\phi}\simeq\overline{a}_{1}\overline{\phi}_{1}\oplus\cdots\oplus\overline{a}_{\overline{I}}\overline{\phi}_{\overline{I}}$ and corresponding decomposition $\mathcal{H}\otimes\mathcal{H}\simeq\bigoplus_{i}\mathbbm{C}^{\overline{a}_{\overline{i}}}\otimes\overline{\mathcal{H}}_{\overline{i}}$ . We choose $\overline{G}$ such that for every $i\in\{1,...,I\}$ , there exists a corresponding $\overline{i}\in\{1,...,\overline{I}\}$ such that $\mathbbm{C}^{\overline{a}_{\overline{i}}}\otimes\overline{\mathcal{H}}_{\overline{i}}\subseteq\mathbbm{C}^{a_{i}}\otimes\mathcal{H}_{i}$ . One may satisfy this condition by choosing $\overline{G}=G$ , but we will see below that for this procedure to scale with the number of qubits we must choose $\overline{G}\subsetneq G$ . We denote the character of the irrep $\overline{\phi}_{\overline{i}}$ by $\chi_{\overline{i}}(U):=\text{Tr}\left[{\overline{\phi}}_{\overline{i}}(U)\right]$ .

Our RB procedure consists of the following steps:

1.

For each $i\in\{1,...,I\}$ , choose an initial state $|\rho_{i}\rangle\rangle$ and measurement projector $|M_{i}\rangle\rangle$ such that $|\langle\langle M_{i}|\hat{P}_{\overline{i}}|\rho_{i}\rangle\rangle|$ is large as possible (see Section III.3 below), where $\hat{P}_{\overline{i}}$ is the projector onto $\overline{\mathcal{H}}_{\overline{i}}$ .
2.

For a given $N$ , choose unitaries $U_{0}\in\overline{G}$ and $U_{1},...,U_{N}\in G$ randomly and uniformly (note elements can be repeated). In the case of a compact group rather than a finite group, choose elements according to the Haar measure. Compute $U_{N+1}=U_{1}^{\dagger}\cdots U_{N}^{\dagger}$ .
3.

Prepare the state $|\rho_{i}\rangle\rangle$ . Apply the gates $(U_{1}U_{0}),U_{2},...,U_{N+1}$ sequentially, where $({U}_{1}{U}_{0})$ is compiled as a single element of $G$ .
4.

Perform a measurement of the observable $M_{i}$ .

Repeat steps 2-4 many times, to estimate the character-weighted survival probability

S_{i}(N)=\frac{1}{|G|^{N+1}}\sum_{\begin{subarray}{c}U_{0}\in\overline{G}\\ U_{1},...,U_{N}\in G\end{subarray}}\chi^{*}_{\overline{i}}(U_{0})\text{Pr}_{U_{0},...,U_{N+1}}

(3)

for each $i$ , where $\text{Pr}_{U_{0},...,U_{N+1}}$ is the probability of measuring $|M_{i}\rangle\rangle$ after applying gates $(U_{1}U_{0}),...,U_{N+1}$ to $|\rho_{i}\rangle\rangle$ , including the effect of gate and SPAM errors.

6.

Repeat steps 2-5 for different values of $N$ .
7.

Fit each character-weighted survival probability to a function of the form

$S_{i}(N)=\sum_{j=1}^{a_{i}}C_{i,j}\lambda_{i,j}^{N}$ (4)

where the $C_{i,j}$ and $\lambda_{i,j}$ are (possibly complex) fitting parameters independent of $N$ . Note that if $\chi_{\overline{i}}$ is complex we may have $S_{i}$ complex, but if $\chi_{\overline{i}}$ is real the $C_{i,j}$ and $\lambda_{i,j}$ are restricted to be real or come in complex-conjugate pairs.
8.

Estimate the average fidelity of the gate error $\Lambda$ as

$F_{\Lambda}=\frac{\sum_{i=1}^{I}\left[\text{dim}(\mathcal{H}_{i})\sum_{j=1}^{a_{i}}\lambda_{i,j}\right]+d}{d^{2}+d}$ (5)

where $d:=2^{n}$ is the dimension of Hilbert space.

A similar RB procedure was first proposed in [22] for groups with all $a_{i}=1$ , the so-called multiplicity-free groups. In this case, each character-weighted survival probability becomes a single exponential decay. Character RB had been previously proposed for the multiplicity-free dihedral group on one qubit [18], and a related approach has been used to simplify standard RB [41].

We note if we omit the initial gate $U_{0}$ and the character-weighting $\chi_{\overline{i}}^{*}(U_{0})$ , we get the method of [19, 20, 21]; in this case, we get a single survival probability $S(N)$ that is given by $S(N)=\sum_{i,j}C_{i,j}\lambda_{i,j}^{N}$ . Determining the $\lambda_{i,j}$ then requires fitting all the parameters $C_{i,j}$ and $\lambda_{i,j}$ simultaneously, and quickly becomes infeasible for a modestly large number of parameters. We see that while both our method and the method of [19, 20, 21] involve simultaneously fitting multiple exponential decays, our method significantly reduces the number of parameters in each fit. For example, if $\phi\simeq 2\phi_{1}\oplus\phi_{2}\oplus\phi_{3}$ , our method requires fitting three functions, corresponding to $\phi_{1}$ , $\phi_{2}$ , and $\phi_{3}$ , where the first function is a sum of two exponential decays and the latter two functions are single exponential decays. In contrast, [19, 20, 21] require fitting a single exponential function that is the sum of four exponential decays, one for each copy of each irrep. In addition, the method of [19, 20, 21] cannot determine $F_{\Lambda}$ , because it is not possible to match the observed parameters $\{\lambda_{i,j}\}$ to their corresponding $\mathcal{H}_{i}$ in order to use Eq. 5.

The remainder of this section is devoted to deriving this procedure, for groups that are not necessarily multiplicity-free. Much of this is a straightforward extension of the derivation of [22], although the generalization to gate-dependent noise (Appendix A) is much less straightforward.

III.1 Deriving the decays

To derive the form of the character-weighted survival, Eq. 4, we will need two facts from representation theory.

Fact 1 (Schur’s Lemma).

Let $\phi:G\rightarrow\mathcal{L}(V)$ be a representation of a group $G$ on a vector space $V$ , which decomposes into irreducible representations as $\phi\simeq a_{1}\phi_{1}\oplus\cdots\oplus a_{I}\phi_{I}$ , where $a_{i}\in\mathbb{Z}^{+}$ are positive integers. The corresponding decomposition of $V$ is $V\simeq\bigoplus_{i}\mathbbm{C}^{a_{i}}\otimes V_{i}$ . In terms of this decomposition, any linear map $\hat{\eta}\in\mathcal{L}(V)$ satisfying $\hat{\eta}\phi(U)=\phi(U)\hat{\eta}$ for all $U\in G$ is of the form

\hat{\eta}\simeq\bigoplus_{i}\hat{Q}_{i}\otimes\hat{\mathbbm{1}}_{i}

(6)

where $\hat{Q}_{i}$ is some $a_{i}\times a_{i}$ matrix for each $i$ .

Fact 2 (Projection formula).

Let $\phi$ and $V$ be as above. Given an irrep $\phi_{i}:G\rightarrow\mathcal{L}(V_{i})$ , define the character $\chi_{i}:G\rightarrow\mathbbm{C}$ of $\phi_{i}$ as $\chi_{i}(U):=\operatorname{Tr}\left(\phi_{i}(U)\right)$ . Then we can write the projector onto $\mathbbm{C}^{a_{i}}\otimes V_{i}$ as

\hat{P}_{i}=\frac{\text{dim}(V_{i})}{|G|}\sum_{U\in G}\chi_{i}(U)^{*}\phi(U).

(7)

For proofs of both facts, see [37].

Given these results, we can prove the key property of $G$ -twirls that allows us to compute the average fidelity.

Theorem 1 (Form of $G$ -twirls).

If $G$ is any unitary group acting on $\mathcal{H}$ , let $\phi\simeq a_{1}\phi_{1}\oplus\cdots\oplus a_{I}\phi_{I}$ be the decomposition of the natural representation into irreps, and let $\mathcal{H}\otimes\mathcal{H}\simeq\bigoplus_{i}\mathbbm{C}^{a_{i}}\otimes\mathcal{H}_{i}$ be the corresponding decomposition of $\mathcal{H}\otimes\mathcal{H}$ . If $\Lambda$ is any quantum channel, the $G$ -twirl of $\Lambda$ is of the form

\hat{\Lambda}_{G}\simeq\bigoplus_{i}\hat{Q}_{i}\otimes\hat{\mathbbm{1}}_{i}\;,

(8)

where $Q_{i}$ is defined as in Fact. 6.

Proof.

We apply Eq. 1 to observe that

	$\displaystyle\hat{\Lambda}_{G}\hat{U}$	$\displaystyle=\frac{1}{\|G\|}\sum_{U^{\prime}\in G}\hat{U}^{\prime\dagger}\hat{\Lambda}\hat{U}^{\prime}\hat{U}$
		$\displaystyle=\frac{1}{\|G\|}\sum_{U^{\prime}\in G}\hat{U}\hat{U}^{\dagger}\hat{U}^{\prime\dagger}\hat{\Lambda}\hat{U^{\prime}}\hat{U}$
		$\displaystyle=\hat{U}\frac{1}{\|G\|}\sum_{(U^{\prime}U)\in G}(\hat{U}^{\prime}\hat{U})^{\dagger}\hat{\Lambda}(\hat{U}^{\prime}\hat{U})=\hat{U}\hat{\Lambda}_{G}$

for any $U\in G$ . We can then apply Fact 6. ∎

We are now ready to derive the formula for the character-weighted survival probability $S_{i}(N)$ . This proof follows the logic of [22], adapted for non-multiplicity-free groups. Our notation assumes finite groups; for compact groups, one simply replaces the discrete average over the group with an integral over the Haar measure. Writing out Eq. 3 explicitly, including the effect of preparation and measurement errors $\Lambda_{P}$ and $\Lambda_{M}$ , we have

\displaystyle{S}_{i}(N)

\displaystyle=\frac{1}{|G|^{N}|\overline{G}|}\sum_{U_{0},...,U_{N}}\underbracket{\chi_{\overline{i}}^{*}(U_{0})}_{\hat{P}_{i}}\langle\langle M_{i}|\hat{\Lambda}_{M}\hat{\Lambda}\hat{U}_{N+1}\hat{\Lambda}\hat{U}_{N}\cdots\hat{\Lambda}\hat{U}_{2}\hat{\Lambda}\hat{U}_{1}\underbracket{\hat{U}_{0}}_{\hat{P}_{\overline{i}}}\hat{\Lambda}_{P}|\rho_{i}\rangle\rangle.

The sum over $U_{0}$ gives the projection $|\overline{G}|\hat{P}_{\overline{i}}/\text{dim}(\overline{\mathcal{H}}_{\overline{i}})$ according to Eq. 7. To do the sum over $U_{1},...,U_{N}$ , we can define new group elements $D_{1},...,D_{N}$ by $D_{i}=U_{i}\cdots U_{1}$ . In terms of the $D_{i}$ , we then have $U_{i}=D_{i}D_{i-1}^{\dagger}$ , with the convention that $D_{N+1}=\mathbbm{1}$ . Note that summing over $U_{1},...,U_{N}$ is the same as summing over $D_{1},...,D_{N}$ . We therefore may write

\displaystyle{S}_{i}(N)

\displaystyle=\frac{1}{\text{dim}(\overline{\mathcal{H}}_{\overline{i}})|G|^{N}}\sum_{D_{1},...,D_{N}\in G}\langle\langle M_{i}|\hat{\Lambda}_{M}\hat{\Lambda}\underbrace{\hat{D}_{N}^{\dagger}\hat{\Lambda}\hat{D}_{N}}_{\hat{\Lambda}_{G}}\cdots\underbrace{\hat{D}_{2}^{\dagger}\hat{\Lambda}\hat{D}_{2}}_{\hat{\Lambda}_{G}}\underbrace{\hat{D}_{1}^{\dagger}\hat{\Lambda}\hat{D}_{1}}_{\hat{\Lambda}_{G}}\hat{P}_{\overline{i}}\hat{\Lambda}_{P}|\rho_{i}\rangle\rangle.

We can now easily perform the sum over the $D_{i}$ , since each sum just gives a $G$ -twirl according to Eq. 1. Performing this sum, and using Thm. 8, gives

	$\displaystyle S_{i}(N)$	$\displaystyle=\frac{1}{\text{dim}(\mathcal{\overline{H}}_{\overline{i}})}\langle\langle M_{i}\|\hat{\Lambda}_{M}\hat{\Lambda}\left(\hat{\Lambda}_{G}\right)^{N}\hat{P}_{\overline{i}}\hat{\Lambda}_{P}\|\rho_{i}\rangle\rangle$
		$\displaystyle=\frac{1}{\text{dim}(\mathcal{\overline{H}}_{\overline{i}})}\langle\langle M_{i}\|\hat{\Lambda}_{M}\hat{\Lambda}\left(\bigoplus_{i^{\prime}}\hat{Q}_{i^{\prime}}\otimes\mathbbm{1}_{i^{\prime}}\right)^{N}\hat{P}_{\overline{i}}\hat{\Lambda}_{P}\|\rho_{i}\rangle\rangle$
		$\displaystyle=\frac{1}{\text{dim}(\mathcal{\overline{H}}_{\overline{i}})}\langle\langle M_{i}\|\hat{\Lambda}_{M}\hat{\Lambda}\left(\hat{Q}_{i}^{N}\otimes\mathbbm{1}_{i}\right)\hat{P}_{\overline{i}}\hat{\Lambda}_{P}\|\rho_{i}\rangle\rangle$

where in the last line, we used the fact that the range of $\hat{P}_{\overline{i}}$ is included in $\mathbbm{C}^{a_{i}}\otimes\mathcal{H}_{i}$ . We see that the effect of the character-weighting is to produce a projector that restricts our attention to a single $i$ . If we diagonalize $\hat{Q}_{i}$ as $\hat{Q_{i}}=\sum_{j=1}^{a_{i}}|e_{i,j}\rangle\rangle\lambda_{i,j}\langle\langle\overline{e}_{i,j}|$ with $\langle\langle\overline{e}_{i,j}|{e}_{i,j^{\prime}}\rangle\rangle=\delta_{j,j^{\prime}}$ , then $\hat{Q}_{i}^{N}=\sum_{j=1}^{a_{i}}|e_{i,j}\rangle\rangle\lambda_{i,j}^{N}\langle\langle\overline{e}_{i,j}|$ , and we may write the final form of $S_{i}(N)$ as

S_{i}(N)=\sum_{j=1}^{a_{i}}\frac{\langle\langle M_{i}|\hat{\Lambda}_{M}\hat{\Lambda}\Big{(}|e_{i,j}\rangle\rangle\langle\langle\overline{e}_{i,j}|\otimes\mathbbm{1}_{i}\Big{)}\hat{P}_{\overline{i}}\hat{\Lambda}_{P}|\rho_{i}\rangle\rangle}{\text{dim}(\overline{\mathcal{H}}_{\overline{i}})}\lambda_{i,j}^{N}

which is precisely the form given in Eq. 4. Notice that the $\lambda_{i,j}$ depend only on the gate error $\Lambda$ , and not the SPAM errors $\Lambda_{P},\Lambda_{M}$ which are absorbed into the constant prefactor.

III.2 Computing the fidelity

Finally, we prove the fidelity can be estimated according to Eq. 5. This was first derived in [21], although we will adopt a simpler proof here using techniques introduced in [13, 14]. The key realization is that both the fidelity and the trace of a channel are invariant under twirling by an arbitrary group: $F_{\Lambda}=F_{\Lambda_{G}}$ and $\operatorname{Tr}(\hat{\Lambda})=\operatorname{Tr}(\hat{\Lambda}_{G})$ (see Eq. 1). In particular, if we choose $G$ to be the full unitary group it is known that the full twirl of a channel is simply a depolarizing channel [14, 13]¹¹1In our notation, this can be seen by noting that the natural representation of the full unitary group decomposes into two irreps which act on $|\mathbbm{1}\rangle\rangle$ and the orthogonal complement of $|\mathbbm{1}\rangle\rangle$ , respectively, and then applying Fact 6.:

\hat{\Lambda}_{G}:=\int dU\ \hat{U}^{\dagger}\hat{\Lambda}\hat{U}=p\mathbbm{1}+(1-p)\frac{1}{d}|\mathbbm{1}\rangle\rangle\langle\langle\mathbbm{1}|.

(9)

In terms of the parameter $p$ , we can directly compute $F_{\Lambda_{G}}=p+\frac{1-p}{d}$ . Similarly, we can also directly compute $\operatorname{Tr}(\hat{\Lambda}_{G})=pd^{2}+(1-p)$ . Combining these equations gives

F_{\Lambda}=\frac{\operatorname{Tr}(\hat{\Lambda})+d}{d^{2}+d}.

(10)

To complete the proof, we note that $\operatorname{Tr}(\hat{\Lambda})$ can be written in terms of the matrices $\hat{Q}_{i}$ in Eq. 8 as

\operatorname{Tr}(\hat{\Lambda})=\sum_{i=1}^{I}\left[\text{dim}(\mathcal{H}_{i})\operatorname{Tr}(\hat{Q}_{i})\right]=\sum_{i=1}^{I}\left[\text{dim}(\mathcal{H}_{i})\sum_{j=1}^{a_{i}}\lambda_{i,j}\right]

which, combined with Eq. 10, gives Eq. 5 as desired.

III.3 Scaling and Feasibility

We note that experimentally determining $S_{i}(N)$ requires Monte Carlo sampling of $U_{0},U_{1},...,U_{N}$ . Each term in this sample is bounded by $\max_{U_{0}\in\overline{G}}(|\chi_{\overline{i}}(U_{0})|)=\text{dim}({\mathcal{H}}_{\overline{i}})$ . Therefore, the standard deviation of the samples is bounded by $\text{dim}({\mathcal{H}}_{\overline{i}})$ , and the sample mean has uncertainty bounded by $\text{dim}({\mathcal{H}}_{\overline{i}})/\sqrt{\text{no. samples}}$ . To determine the relative uncertainty, we consider $S_{i}(N)\approx\sum_{j=1}^{a_{i}}C_{i,j}$ which is given by

	$\displaystyle\sum_{j=1}^{a_{i}}C_{i,j}$	$\displaystyle=\sum_{j=1}^{a_{i}}\frac{\langle\langle M_{i}\|\hat{\Lambda}_{M}\hat{\Lambda}\Big{(}\|e_{i,j}\rangle\rangle\langle\langle\overline{e}_{i,j}\|\otimes\mathbbm{1}_{i}\Big{)}\hat{P}_{\overline{i}}\hat{\Lambda}_{P}\|\rho_{i}\rangle\rangle}{\text{dim}(\mathcal{\overline{H}}_{\overline{i}})}$
		$\displaystyle\approx\frac{\langle\langle M_{i}\|\hat{P}_{\overline{i}}\|\rho_{i}\rangle\rangle}{\text{dim}(\mathcal{\overline{H}}_{\overline{i}})}$

where we’ve approximated $\Lambda,\Lambda_{M},\Lambda_{P}\approx\mathbbm{1}$ . The relative uncertainty in $S_{i}(N)$ is therefore bounded by

\frac{\sigma_{i}}{|S_{i}(N)|}\lesssim\frac{\text{dim}(\overline{\mathcal{H}}_{\overline{i}})^{2}}{|\langle\langle M_{i}|\hat{P}_{\overline{i}}|\rho_{i}\rangle\rangle|\sqrt{\text{no. samples}}}

We see that to efficiently benchmarking a group $G$ , we must have $I$ , $a_{i}$ , and $\dim(\overline{\mathcal{H}}_{\overline{i}})$ all small. $I$ must be small so that we only need to estimate a small number of character-weighted survival probabilities $S_{i}(N)$ , $a_{i}$ must be small so that we may fit a function with a small number of parameters, and $\dim(\overline{\mathcal{H}}_{\overline{i}})$ must be small for our Monte Carlo estimation of $S_{i}(N)$ to converge quickly. Note that for any $G$ the natural representation satisfies $\sum_{i=1}^{I}a_{i}\dim(\mathcal{H}_{i})=4^{n}$ where $n$ is the number of qubits, so that choosing $\overline{G}=G$ will not suffice if the number of qubits is large. In particular, to scalably benchmark a group, we must choose $G$ so that the number of irreps $I$ grows slowly with $n$ , the multiplicity $a_{i}$ of each irrep is bounded by a small constant, and $\overline{G}$ has corresponding irreps $\overline{\mathcal{H}}_{\overline{i}}$ whose dimension grows slowly with $n$ . These scaling considerations are similar to those discussed in [22] for multiplicity-free RB, except in our case we allow $a_{i}$ to be bounded rather than strictly $1$ .

Note that the optimal $|\rho_{i}\rangle\rangle$ with largest $|\langle\langle M_{i}|\hat{P}_{\overline{i}}|\rho_{i}\rangle\rangle|$ is necessarily a pure state, since any mixed state $|\rho_{i}\rangle\rangle=\sum_{\gamma}p_{\gamma}|\psi_{\gamma}\rangle\rangle$ has

|\langle\langle M_{i}|\hat{P}_{\overline{i}}|\rho_{i}\rangle\rangle|\leq\sum_{\gamma}p_{\gamma}|\langle\langle M_{i}|\hat{P}_{\overline{i}}|\psi_{\gamma}\rangle\rangle|\leq\max_{\gamma}|\langle\langle M_{i}|\hat{P}_{\overline{i}}|\psi_{\gamma}\rangle\rangle|.

Ref. [22] considered the case of mixed initial states, and included a protocol for sampling from a mixed state $|\rho_{i}\rangle\rangle=\sum_{\gamma}p_{\gamma}|\psi_{\gamma}\rangle\rangle$ provided one can efficiently prepare the states $\{|\psi_{\gamma}\rangle\rangle\}$ . However, we see that it suffices to take the initial state to be one of the efficiently preparable $|\psi_{\gamma}\rangle\rangle$ , which simplifies initial state preparation.

Our scaling estimates are based on the typical case; however, there are a few worst-case failure modes. First, the noise may have some symmetry that restricts $\langle\langle\overline{e}_{i,j}|\hat{P}_{\bar{i}}\approx 0$ for some $(i,j)$ . In this case, the corresponding $\lambda_{i,j}$ will not be accurately estimated by the fitting function. To remedy this, one may choose a set of projectors $\left\{\hat{P}_{\overline{i},1},...,\hat{P}_{\overline{i},k}\right\}$ such that each $\langle\langle\overline{e}_{i,j}|$ has overlap with at least one $\hat{P}_{\overline{i},\alpha}$ . This requires at most $a_{i}$ projectors. We can then define

\hat{P}_{\overline{i}}=\sum_{\alpha}\hat{P}_{\overline{i},\alpha}\quad\chi_{\overline{i}}=\sum_{\alpha}\chi_{\overline{i},\alpha}.

The modified character-weighted survival probability will require taking additional data to achieve the same relative uncertainty, since the corresponding $\dim(\overline{\mathcal{H}}_{\overline{i}})=\sum_{\alpha}\dim(\overline{\mathcal{H}}_{\overline{i},\alpha})$ will be larger, but is otherwise identical.

The fitting procedure may also have difficulty fitting multiple exponential decays [43, 44], especially if the decay rates are similar [44]. In the case of similar decays, the fit might have numerous local minima; worse, the fitting function might simply set the coefficient of one of the decays to zero and the corresponding decay rate to some arbitrary value, and fit the curve using fewer exponential decays. This can be detected during the fitting procedure, and corrected by either taking more data to more closely constrain the fit or by simply fitting fewer exponential decays. For a detailed discussion of methods used to fit multiexponential decays and their failure modes, we refer to [45, 46, 47].

IV Application: Subspace randomized benchmarking

As an application of the general character RB method, we can improve on the recently introduced subspace randomized benchmarking method [23]. Subspace RB characterizes the error associated with a group of gates $G$ that preserve a subspace of the Hilbert space. In [23], a benchmarking procedure is introduced that yields two decay parameters that are functions of the noise channel, but the procedure does not give an estimate for the average fidelity or other quantities with simple physical interpretations. The multiplicity-free character RB of [22] is not directly applicable to this situation, as we will see that any group that preserves subspaces necessarily decomposes into irreps with multiplicity. However, using our method we can easily characterize the average fidelity of such gates.

To simplify our discussion, we will focus on the particular case discussed in [23]. The system considered in [23] can implement arbitrary symmetric single qubit gates $U_{1}:=U\otimes U$ as well as the two-qubit entangling gate $U_{ZZ}:=\exp\{-i\frac{\pi}{4}Z\otimes Z\}$ . The symmetric single qubit gates have negligible error compared to the entangling gate, so the goal of the experiment is to characterize the fidelity of $U_{ZZ}$ . This is accomplished by combining the elementary gates into elements of a benchmarking group $G$ , using a fixed number of the relevant gate $U_{ZZ}$ , and then designing an RB procedure to benchmark elements of $G$ . It is straightforward to see that any $U\in G$ made up of products of $U_{1}$ and $U_{ZZ}$ operators preserves the triplet and singlet subspaces

	$\displaystyle\mathcal{H}_{T}$	$\displaystyle:=\text{span}\left\{\|00\rangle,\frac{\|01\rangle+\|10\rangle}{\sqrt{2}},\|11\rangle\right\}$
	$\displaystyle\mathcal{H}_{S}$	$\displaystyle:=\text{span}\left\{\frac{\|01\rangle-\|10\rangle}{\sqrt{2}}\right\}.$

This implies that every gate $U\in G$ decomposes as $U=U_{T}\oplus U_{S}$ , with $U_{T}$ and $U_{S}$ acting on the triplet and singlet spaces, respectively.

Our method differs from the original in several ways. Most notably, we combine the elementary gates into elements $U\in G$ such that $G$ forms a group. This requires a moderate increase in complexity of the combined gates; [23] combined their gates into unitaries involving three $U_{ZZ}$ gates, while our construction requires four. However, in return for this increased complexity, our method offers several advantages. Rather than estimate decay parameters with no clear physical interpretation, our method produces direct estimates of the average fidelity. In addition, the derivation of the form of the exponential decays in [23] required assumptions on the relative phases of $U_{T}$ and $U_{S}$ that could not actually be realized on their experimental platform. In contrast, our method yields rigorous decays thanks to the underlying group structure of $G$ .

The original subspace RB can be extended to sets of gates $G$ that preserve some arbitrary splitting of $\mathcal{H}$ into subspaces $\mathcal{H}=\mathcal{H}_{1}\oplus\mathcal{H}_{2}$ provided the set $G$ can be written as

G=\{U_{1,b_{1}}\oplus\sigma U_{2,b_{2}}:\sigma=\pm,\ (b_{1},b_{2})\in B_{1}\times B_{2}\}

where $G_{1}:=\{U_{1,b_{1}}:b_{1}\in B_{1}\}$ and $G_{2}:=\{U_{2,b_{2}}:b_{2}\in B_{2}\}$ are groups and unitary 2-designs ²²2Ref. [23] claimed it was sufficient to require $G_{2}$ to be a unitary $1$ -design, but this appears to be an error. A similar error was made in [26], from which much of [23] is derived. (see below for the definition of a $2$ -design) acting on $\mathcal{H}_{1}$ and $\mathcal{H}_{2}$ respectively (here, $B_{1}$ and $B_{2}$ are just index sets for the groups $G_{1}$ and $G_{2}$ ). However, it is difficult to construct such a $G$ in a way that is experimentally relevant; indeed, [23] could not do this for the simple case of two qubits, and we avoid attempting such a construction here. A more useful approach, which mirrors our approach below, is to construct an arbitrary group out of the elementary gates and perform character RB on whatever irreps result. This method can likely be used to benchmark other two-qubit gates that are symmetric under SWAP besides $U_{ZZ}$ , and may also prove useful for gates that preserve other subspaces.

IV.1 Constructing the benchmarking group

Ref. [23] constructed their benchmarking set $G$ using a generalization of the Clifford group [11, 11, 12] to a $d$ -level system [49]. We will follow a similar procedure, modified to ensure $G$ forms a group. For a $d$ -level system, analogues of the $X$ and $Z$ qubit operators are defined as [50]:

X|z\rangle=|z+1\rangle\qquad Z|z\rangle=\omega^{z}|z\rangle

where $\omega:=e^{\frac{2\pi i}{d}}$ and addition is performed modulo $d$ . These generalized $X$ and $Z$ operators are unitary but not Hermitian, and the set $\{X^{a}Z^{b}:a,b\in\mathbbm{Z}_{d}\}$ forms a (complex) orthogonal basis for the set of all $d\times d$ matrices. Note that for $d=2$ we recover the usual Pauli matrices.

Specializing to $d=3$ , define the generalized Pauli group as $\mathcal{P}:=\{\omega^{\eta}X^{a}Z^{b}:\eta,a,b\in\mathbbm{Z}_{d}\}$ . The fact that $\mathcal{P}$ is a group follows from the commutation relation $ZX=\omega XZ$ . The generalized Clifford group is defined to be the set of all unitaries that stabilize $\mathcal{P}$ [49]:

G_{T}=\{U:U\mathcal{P}U^{\dagger}=\mathcal{P}\}.

An element $U\in G_{T}$ is defined (up to a global phase) by its action on $X$ and $Z$ . Defining $UXU^{\dagger}=\omega^{\eta_{x}}X^{a_{x}}Z^{b_{x}}$ and $UZU^{\dagger}=\omega^{\eta_{z}}X^{a_{z}}Z^{b_{z}}$ , and noting

	$\displaystyle ZX$	$\displaystyle=\omega XZ$
	$\displaystyle UZU^{\dagger}UXU^{\dagger}$	$\displaystyle=\omega UXU^{\dagger}UZU^{\dagger}$
	$\displaystyle\omega^{\eta_{x}+\eta_{z}}X^{a_{z}}Z^{b_{z}}X^{a_{x}}Z^{b_{x}}$	$\displaystyle=\omega^{1+\eta_{x}+\eta_{z}}X^{a_{x}}Z^{b_{x}}X^{a_{z}}Z^{b_{z}}$
	$\displaystyle\omega^{a_{x}b_{z}}X^{a_{x}+a_{z}}Z^{b_{x}+b_{z}}$	$\displaystyle=\omega^{1+a_{z}b_{x}}X^{a_{x}+a_{z}}Z^{b_{x}+b_{z}}$

we see that we must have $a_{x}b_{z}-a_{z}b_{x}=_{3}1$ , where $=_{3}$ denotes equality mod $3$ . This is the only restriction on $\eta_{x},\eta_{z},a_{x},a_{z},b_{x},b_{z}$ [49], leading to a total of 216 elements of $G_{T}$ . We can find the action of $U\in G_{T}$ on a general element $X^{a}Z^{b}$ by

\displaystyle\begin{split}UX^{a}Z^{b}U^{\dagger}=&\ (UXU^{\dagger})^{a}(UZU^{\dagger})^{b}\\ =&\ \omega^{P}X^{aa_{x}+ba_{z}}Z^{ab_{x}+bb_{z}}\end{split}

where

P:=\eta_{x}a+\eta_{z}b+2(a^{2}-a)a_{x}b_{x}+2(b^{2}-b)a_{z}b_{z}+abb_{x}a_{z}.

The action of $U$ on a general density matrix then follows by linearity.

Our benchmarking group $G$ is constructed by combining the elementary symmetric gates to act as $G_{T}$ on the triplet subspace, where the three levels $|0\rangle,|1\rangle,|2\rangle$ correspond to the triplet basis $|00\rangle,\frac{|01\rangle+|10\rangle}{\sqrt{2}},|11\rangle$ . The most general composite gate is formed by alternatively applying $U_{1}$ and $U_{ZZ}$ gates to our qubits. A straightforward calculation shows that if such a circuit applies an operator $U_{T}$ to the triplet subspace, its action on the singlet subspace is necessarily given by $(-1)^{n_{z}}\omega^{\eta}\det(U_{T})^{1/3}$ , where $n_{z}$ is the number of entangling $U_{ZZ}$ gates. By varying the single-qubit unitaries $U_{1}$ , we find computationally that all elements of $G_{T}$ and all relative phases $\omega^{\eta}$ can be generated by circuits of exactly four $U_{ZZ}$ gates, as shown in Fig. 1 ³³3Ref. [23] required a shorter circuit of only three entangling gates. However, this circuit cannot implement all relative phases between the subspaces and thus does not result in a group.. In total, then, the benchmarking group is given by

G:=\{U_{T}\oplus\omega^{\eta}\det(U_{T})^{1/3}:U_{T}\in G_{T},\eta=0,1,2\}

where the first summand acts on the triplet subspace and the second acts on the singlet subspace. Note that every group element contains exactly four entangling gates, so the average fidelity of $G$ gives a useful measure of the fidelity of the entangling gate.

Refer to caption — Figure 1: The elements of the benchmarking group $G$ are constructed by composing elementary gates as shown above to implement elements of $G_{T}$ on the triplet subspace. Each group element contains exactly four entangling gates.

Subrep	Projector	$\chi_{i}(U_{T}\oplus U_{S})$
$\mathcal{H}_{T0}$	$\hat{P}_{T0}=\frac{1}{3}\|\mathbbm{1}_{T}\rangle\rangle\langle\langle\mathbbm{1}_{T}\|$	$1$
$\mathcal{H}_{S0}$	$\hat{P}_{S0}=\|\mathbbm{1}_{S}\rangle\rangle\langle\langle\mathbbm{1}_{S}\|$	$1$
$\mathcal{H}_{T\perp}$	$\hat{P}_{T\perp}=\mathbbm{1}_{T}-\hat{P}_{T0}$	$\|\operatorname{Tr}(U_{T})\|^{2}-1$
$\mathcal{H}_{TS}$	$\hat{P}_{TS}=\text{Projector onto }\mathcal{H}_{T}\otimes\mathcal{H}_{S}$	$\operatorname{Tr}(U_{T})\operatorname{Tr}(U_{S})^{*}$
$\mathcal{H}_{ST}$	$\hat{P}_{ST}=\text{Projector onto }\mathcal{H}_{S}\otimes\mathcal{H}_{T}$	$\operatorname{Tr}(U_{T})^{*}\operatorname{Tr}(U_{S})$

Table 1: Subrepresentations of the standard representation for groups that preserve the triplet and singlet subspaces, and their corresponding projectors and characters.

IV.2 Irreps of the benchmarking group

For $G$ given above, the natural representation decomposes into the irreps $\mathcal{H}_{T0}$ , $\mathcal{H}_{S0}$ , $\mathcal{H}_{T\perp}$ , $\mathcal{H}_{TS}$ , and $\mathcal{H}_{ST}$ , which are described in Table 1. These are all clearly subrepresentations of the natural representation; for proof that they are in fact irreducible, we will use the concept of a unitary $\mathbf{t}$ -design [9].

Let $S$ be a set of unitaries acting on a space $\mathcal{H}$ . A balanced polynomial of degree $t$ is a polynomial in the matrix elements of $U$ and $U^{*}$ where each term in the polynomial has degree $d<t$ in the elements of $U$ and degree $d$ in the elements of $U^{*}$ . $S$ is a unitary $t$ -design if for balanced polynomial $p(U,U^{*})$ of degree $t$ , averaging $p(U,U^{*})$ over $S$ is the same as averaging over all unitaries on $\mathcal{H}$ (weighted by the Haar measure)

\frac{1}{|S|}\sum_{U\in S}p(U,U^{*})=\int dU\ p(U,U^{*}).

A classic example is the Clifford group, which forms a unitary $3$ -design [9, 52, 53].

The group $G_{T}$ forms a unitary 2-design [54] (see Appendix B for a proof). This allows us to prove the representations in Table 1 are irreducible, using the following fact:

Fact 3 (Schur normalization).

Let $\chi$ be the character of a representation. The representation is irreducible iff

\frac{1}{|G|}\sum_{U\in G}|\chi(U)|^{2}=1.

For a proof, see [37].

The representations $\mathcal{H}_{T0}$ and $\mathcal{H}_{S0}$ are 1D, thus irreducible. For the representation $\mathcal{H}_{T\perp}$ , we have

	$\displaystyle\frac{1}{\|G\|}\sum_{U\in G}\|\chi_{T\perp}(U)\|^{2}$	$\displaystyle=\frac{1}{3\|G_{T}\|}\sum_{\begin{subarray}{c}U_{T}\in G_{T}\\ \eta=0,1,2\end{subarray}}\|\chi_{T\perp}(U_{T})\|^{2}$
		$\displaystyle=\frac{1}{\|G_{T}\|}\sum_{G_{T}}\left(\|\operatorname{Tr}(U_{T})\|^{2}-1\right)^{2}$
		$\displaystyle=\int dU_{\alpha}\ \left(\|\operatorname{Tr}(U_{\alpha})\|^{2}-1\right)^{2}$
		$\displaystyle=1$

where the second equality follows from the unitary 2-design property, and the third follows from the fact that $\mathcal{H}_{T\perp}$ is an irrep of the natural representation of the full unitary group on $\mathcal{H}_{T}$ . Finally, for $\mathcal{H}_{TS}$ and $\mathcal{H}_{ST}$ we have

	$\displaystyle\frac{1}{\|G\|}\sum_{U\in G}\|\chi_{ST}(U)\|^{2}$	$\displaystyle=\frac{1}{3\|G_{T}\|}\sum_{\begin{subarray}{c}U_{T}\in G_{T}\\ \eta=0,1,2\end{subarray}}\|\operatorname{Tr}(U_{T})\|^{2}$
		$\displaystyle=\int dU_{T}\ \|\operatorname{Tr}(U_{T})\|^{2}$
		$\displaystyle=1$

where the second equality follows from the unitary 2-design property and the third follows from the fact that the direct representation of the full unitary group on $\mathcal{H}_{T}$ is irreducible.

Note that $\mathcal{H}_{T0}$ and $\mathcal{H}_{S0}$ are two irreducible copies of the trivial representation, so that $G$ is necessarily non-multiplicity-free ⁴⁴4It follows that $G$ also cannot form a 2-design, as 2-designs are always multiplicity free; in particular, the natural representation of a 2-design decomposes into precisely two non-isomorphic irreps, acting on $|\mathbbm{1}\rangle\rangle$ and the orthogonal complement of $|\mathbbm{1}\rangle\rangle$ [9, 22].. The remaining irreps are all unique, since they have different character functions.

IV.3 Benchmarking $G$

The form of the decay curves corresponding to each irrep is given by

$\displaystyle S_{0}(N)$	$\displaystyle=C_{0}\lambda_{0}^{N}+B$	(11)
$\displaystyle S_{TS}(N)$	$\displaystyle=C_{TS}\lambda_{TS}^{N}$
$\displaystyle S_{ST}(N)$	$\displaystyle=C_{ST}\lambda_{ST}^{N}$
$\displaystyle S_{T\perp}(N)$	$\displaystyle=C_{T\perp}\lambda_{TS}^{N}.$

Note that from our general form Eq. 4 we would expect that $S_{0}(N)$ is the sum of two exponential terms, with each $\lambda_{0,j}$ corresponding to an eigenvalue of $\hat{\Lambda}_{G}$ restricted to $\mathcal{H}_{0}$ . However, we know that for trace-preserving noise $\langle\langle\mathbbm{1}|\hat{\Lambda}_{G}=\langle\langle\mathbbm{1}|$ , which implies that one of the eigenvalues is $1$ .

We define two different subgroups $\overline{G}_{1},\overline{G}_{2}\subseteq G$ for our benchmarking procedure. We will use $\overline{G}_{1}$ to construct $S_{0}(N)$ and $S_{T\perp}(N)$ , and $\overline{G}_{2}$ to construct $S_{TS}(N)$ and $S_{ST}(N)$ . We define

	$\displaystyle\overline{G}_{1}$	$\displaystyle:=\{X^{a}Z^{b}\oplus\omega^{\eta}:a,b,\eta=0,1,2\}$
	$\displaystyle\overline{G}_{2}$	$\displaystyle:=\{Z^{b}\oplus\omega^{\eta}:b,\eta=0,1,2\}.$

For $\overline{G}_{1}$ , we can define the following character functions and their corresponding projectors:

	$\displaystyle\chi_{\overline{0}}(X^{a}Z^{b}\oplus\omega^{\eta})$	$\displaystyle=1$
	$\displaystyle\hat{P}_{\overline{0}}$	$\displaystyle=\frac{1}{3}\|\mathbbm{1}_{T}\rangle\rangle\langle\langle\mathbbm{1}_{T}\|+\|\mathbbm{1}_{S}\rangle\rangle\langle\langle\mathbbm{1}_{S}\|$
	$\displaystyle\chi_{\overline{T\perp}}(X^{a}Z^{b}\oplus\omega^{\eta})$	$\displaystyle=\omega^{-a}$
	$\displaystyle\hat{P}_{\overline{T\perp}}$	$\displaystyle=\frac{1}{3}\|Z\rangle\rangle\langle\langle Z\|$

We see that $\hat{P}_{\overline{0}}$ projects into $2\overline{\mathcal{H}}_{\overline{0}}\subseteq 2\mathcal{H}_{0}$ and $\hat{P}_{\overline{T\perp}}$ projects into $\overline{\mathcal{H}}_{\overline{T\perp}}\subseteq\mathcal{H}_{T\perp}$ , as required. We also see that $\text{dim}(\overline{\mathcal{H}}_{\overline{T\perp}})=1$ , so that $S_{T\perp}(N)$ will have the best possible relative error (see Section III.3).

For $\overline{G}_{2}$ , we can define the character functions and corresponding projectors

	$\displaystyle\chi_{\overline{TS}}(Z^{b}\oplus\omega^{\eta})$	$\displaystyle=\omega^{b-\eta}$
	$\displaystyle\hat{P}_{\overline{TS}}$	$\displaystyle=\|T\rangle\|S\rangle\langle T\|\langle S\|$
	$\displaystyle\chi_{\overline{ST}}(Z^{b}\oplus\omega^{\eta})$	$\displaystyle=\omega^{-b+\eta}$
	$\displaystyle\hat{P}_{\overline{ST}}$	$\displaystyle=\|S\rangle\|T\rangle\langle S\|\langle T\|$

where $|T\rangle:=(|01\rangle+|10\rangle)/\sqrt{2}$ is the triplet state satisfying $Z|T\rangle=\omega|T\rangle$ and $|S\rangle:=(|01\rangle-|10\rangle)/\sqrt{2}$ is the singlet state. We again see that $P_{\overline{TS}}$ projects into $\overline{\mathcal{H}}_{\overline{TS}}\subseteq\mathcal{H}_{TS}$ and $\text{dim}(\overline{\mathcal{H}}_{\overline{TS}})=1$ , so that $S_{TS}(N)$ will also have the best possible relative error.

As our initial states, we choose

|\rho_{i}\rangle\rangle=\left\{\begin{array}[]{ll}|00\rangle\rangle,&i=0,T\!\!\perp\\ |01\rangle\rangle,&i=TS,ST\end{array}\right.

Here, we’ve restricted ourselves to initial states that are a mixture of $Z$ -basis product states, for ease of preparation.

As our measurement projectors, we choose

|M_{i}\rangle\rangle=\left\{\begin{array}[]{ll}|00\rangle\rangle+|11\rangle\rangle,&i=0,T\!\!\perp\\ |01\rangle\rangle,&i=TS,ST\\ \end{array}\right.

Here, we’ve restricted our measurement projectors to correspond to $Z$ measurements, for ease of measuring.

With these choices, the $S_{i}(N)$ are approximately

\displaystyle S_{i}(N)

\displaystyle\approx\frac{\langle\langle M_{i}|\hat{P}_{\overline{i}}|\rho_{i}\rangle\rangle}{\text{dim}(\mathcal{\overline{H}}_{\overline{i}})}=\left\{\begin{array}[]{ll}\frac{2}{3},&{i}={0}\\ \frac{e^{-i\pi/3}}{3},&{i}=T\!\!\perp\\ \frac{1}{4},&{i}=TS,ST\\ \end{array}\right.

Note that $\lambda_{ST}=\lambda_{TS}^{*}$ , so it is unnecessary to compute both $S_{TS}(N)$ and $S_{ST}(N)$ . Note also that $\lambda_{0}$ and $\lambda_{T\perp}$ are both necessarily real, as are $C_{0}$ and $B$ . The remaining parameters are complex. For convenience, we will rotate $S_{T\perp}(N)$ by $e^{i\pi/3}$ so that $S_{T\perp}(N)$ is approximately real.

We demonstrate our method by generating random error channels and simulating our RB procedure. To generate a random error channel $\Lambda$ on a $d$ -dimensional Hilbert space, we generate a random unitary on a $(d^{2}+d)$ dimensional Hilbert space and trace out $d^{2}$ auxiliary degrees of freedom; to adjust the fidelity, we take a convex combination of the resulting channel with the identity channel. All channels generated by this method are guaranteed to be completely positive trace-preserving (CPTP), thus valid error channels, and every CPTP channel can be generated via this method [36]. For each error channel, we take data at 15 different values of $N$ , and sample unitary operators at each value of $N$ until we have applied a total of $150,000$ unitary operators in total. For each string of unitary operators, we perform full state-vector simulation to apply the RB sequence of operators, and then generate a measurement outcome of $0$ or $1$ using the appropriate probability, and compute the character-weighted average. In Fig. 2, we show the exact value of $S_{i}(N)$ , the data we take to estimate $S_{i}(N)$ , and the fit to $S_{i}(N)$ according to Eq. 11 for a single random error channel $\Lambda$ .

From the fit data, we can estimate $F_{\Lambda}$ by applying Eq. 5:

F_{\Lambda}=\frac{1+\lambda_{0}+8\lambda_{T\perp}+3\lambda_{TS}+3\lambda_{ST}+4}{20}.

(12)

Note that the imaginary parts of $\lambda_{TS}$ and $\lambda_{ST}$ always cancel to give a real $F_{\Lambda}$ as expected. We use this formula to estimate the fidelity of our randomly generated error channels, and compare our estimate to the true fidelity in Fig. 3. We see that the true fidelity and the estimated fidelity agree within the error bars set by the uncertainty of our fits.

We can directly compare this with the original subspace RB method [23]. That method served to estimate only $\lambda_{0}$ and $\lambda_{T\perp}$ ( $t$ and $r$ in their notation), and they could only form a measure of gate fidelity using these quantities. They defined a so-called “extended sub-fidelity” $\tilde{F}_{\Lambda}$ , which they obtained by replacing $\lambda_{ST}$ and $\lambda_{TS}$ with the weighted average of the other eigenvalues: $\lambda_{ST}+\lambda_{TS}\approx 2\frac{1+\lambda_{0}+8\lambda_{T\perp}}{10}$ . Explicitly, the extended sub-fidelity is given by ⁵⁵5Our formula differs slightly from the corresponding formula in [23]. Ref. [23] considered approximating the process (also called entanglement) fidelity rather than the average fidelity; however, the average fidelity can be determined from the process fidelity[14, 13]. To be consistent with the rest of our paper, we have translated their approximation of the process fidelity into the corresponding approximation of the average fidelity.

\tilde{F}_{\lambda}=\frac{16\lambda_{T\perp}+2\lambda_{0}+7}{25}.

It is obvious that if $F_{\Lambda}\rightarrow 1$ , $\tilde{F}_{\Lambda}\rightarrow 1$ as well, but the reverse is not necessarily true. We can compare the extended sub-fidelity to the exact fidelity for the various noise sources explored in [23]. We consider intensity errors, which correspond to an overrotation $e^{-i\epsilon ZZ}$ ; optical pumping errors, which cause amplitude-damping on each qubit; inhomogenous fields, which cause phase-damping on each qubit; and SWAP errors, which interchange the qubits. The results are shown in Fig. 4. We see that while for most error sources $F_{\Lambda}\approx\tilde{F}_{\Lambda}$ , there exist worse-case errors, such as SWAP, that cannot be detected by $\tilde{F}_{\Lambda}$ . This was also noted in [23] as a limitation of their method.

Our work also improves upon the original work in the mathematical assumptions needed to derive the benchmarking decays. Ref. [23] derived their decay formulas under the assumption that their benchmarking set was of the form $\{U_{T}\oplus\sigma\phi_{U_{T}}:U_{T}\in G_{T},\sigma=\pm\}$ , where $\phi_{U_{T}}$ is some uncontrolled phase that occurs on the singlet space and $\sigma$ is a controllable phase between the singlet and triplet spaces. However, in practice they could not control $\sigma$ using a constant number of $U_{ZZ}$ gates. Instead, they implemented only $\{U_{T}\oplus\phi_{U_{T}}:U_{T}\in G_{T}\}$ and assumed the form of the decay would not change. In our work, by contrast, we have rigorously derived decay formulas for a group of gates that can be directly compiled into elementary symmetric gates using a constant number of $U_{ZZ}$ .

We note that our method does require one additional capability that was not required in the original work: in order to estimate $S_{TS}(N)$ , it is necessary to initialize and measure the $|01\rangle$ state. This requires additional experimental overhead to individually address and measure each qubit at the beginning and end of the benchmarking procedure. However, such overhead only contributes to the SPAM errors $\Lambda_{P},\Lambda_{M}$ , and does not affect our estimates of the entangling error. In any case, our method to measure $\lambda_{0}$ and $\lambda_{T\perp}$ does not require individual addressing, and can be viewed as a mathematically rigorous method to extract these parameters with no additional experimental requirements.

V Application: Leakage randomized benchmarking

We may also use our generalized character RB to improve the leakage RB introduced in [26]. In leakage RB, like subspace RB, one is given a group $G$ that preserves the splitting of the Hilbert space into subspaces $\mathcal{H}=\mathcal{H}_{1}\oplus\mathcal{H}_{2}$ . In leakage RB, however, $\mathcal{H}_{1}\oplus\mathcal{H}_{2}$ does not represent the computational Hilbert space, and the goal is not to compute the average fidelity of the group operations. Instead, $\mathcal{H}_{1}$ represents the computational space of a quantum system (e.g. the two lowest-level states that encode a qubit), while $\mathcal{H}_{2}$ represents the leakage space outside the computational space. Leakage RB determines the average probability of “leaking” from $\mathcal{H}_{1}$ to $\mathcal{H}_{2}$ or “seeping” from $\mathcal{H}_{2}$ to $\mathcal{H}_{1}$ . Noting that the probability of a state $|\rho\rangle\rangle$ being in subspace $\alpha=1,2$ is given by $\langle\langle\mathbbm{1}_{\alpha}|\rho\rangle\rangle$ , define the leakage $L$ and seepage $S$ by

	$\displaystyle L:=\int d\psi_{1}\langle\langle\mathbbm{1}_{2}\|\hat{\Lambda}\|\psi_{1}\rangle\rangle=\frac{1}{d_{1}}\langle\langle\mathbbm{1}_{2}\|\hat{\Lambda}\|\mathbbm{1}_{1}\rangle\rangle$		(13)
	$\displaystyle S:=\int d\psi_{2}\langle\langle\mathbbm{1}_{1}\|\hat{\Lambda}\|\psi_{2}\rangle\rangle=\frac{1}{d_{2}}\langle\langle\mathbbm{1}_{1}\|\hat{\Lambda}\|\mathbbm{1}_{2}\rangle\rangle.$		(14)

In addition, leakage RB determines the average fidelity restricted to the subspace $\mathcal{H}_{1}$

F_{\Lambda,1}=\int d\psi_{1}\langle\langle\psi_{1}|\hat{\Lambda}|\psi_{1}\rangle\rangle.

(15)

which is the appropriate measure of gate quality, since all computations take place in $\mathcal{H}_{1}$ . Leakage RB is relevant for any system in which qubits are encoded in the subspace of a larger Hilbert space, which includes superconducting qubits [57, 58], quantum dots, [59, 60, 61, 62, 63], and trapped ions [64, 65, 66].

The original leakage RB could only be applied to a group

G=\{U_{1,b_{1}}\oplus\sigma U_{2,b_{2}}:(b_{1},b_{2})\in B_{1}\times B_{2},\ \sigma=\pm 1\}

(16)

such that $G_{1}=\{U_{1,b_{1}}:b_{1}\in B_{1}\}$ and $G_{2}=\{U_{2,b_{2}}:b_{2}\in B_{2}\}$ form $2$ -designs on their respective subspaces ⁶⁶6Ref. [26] originally claimed it was sufficient for $G_{2}$ to be a unitary 1-design, but this appears to be an error. This is a very stringent condition, as it requires being able to independently control the computational and leakage subspaces. In many experimental implementations such control is not realistic; an experimental implementation of a gate $U_{1,b}$ on the computational subspace will naturally implement some $U_{2,b}$ on the leakage subspace. It is therefore desirable to develop a leakage RB that can be applied to more general groups.

Using our method, we can derive a leakage RB procedure that is more general than the one described in [26]. Let $G$ be a group of unitary gates that preserve the subspaces of $\mathcal{H}$ , and let $\Lambda$ be their shared error channel. To estimate $L$ and $S$ , we will require that the only trivial representations of $G$ are $|\mathbbm{1}_{1}\rangle\rangle$ and $|\mathbbm{1}_{2}\rangle\rangle$ , while to estimate $F_{\Lambda,1}$ we additionally require that the subrepresentation $\mathcal{H}_{1\perp}\subseteq\mathcal{H}_{1}\otimes\mathcal{H}_{1}$ orthogonal to $|\mathbbm{1}_{1}\rangle\rangle$ is an irrep of multiplicity $1$ .

If we write our group $G$ as

	$\displaystyle G$	$\displaystyle=\{U_{b,\sigma}:b\in B,\ \sigma=\pm 1\}$
		$\displaystyle=\{U_{1,b}\oplus\sigma U_{2,b}:b\in B,\ \sigma=\pm 1\}.$

then the first condition is satisfied provided $\{U_{1,b}:b\in B\}$ and $\{U_{2,b}:b\in B\}$ are unitary $1$ -designs, while the second condition is satisfied if provided these groups are unitary $2$ -designs with dimensions $d_{1}\neq d_{2}$ (see Appendix C for proofs). Note that our requirements are significantly weaker than the original leakage RB, as we are only assuming the ability to implement an independent phase on the leakage space.

We outline our procedure for determining $L$ , $S$ , and $F_{\Lambda,1}$ for such groups $G$ . Our procedure, like the original leakage RB, requires that SPAM errors do not mix the the subspaces $\mathcal{H}_{1}$ and $\mathcal{H}_{2}$ , or at least that such mixing is negligible compared to the gate errors. In our derivations we will assume $\hat{\Lambda}_{M}=\hat{\Lambda}_{P}=\hat{\mathbbm{1}}$ , although the generalization to errors that act only within the subspaces is trivial.

Our modified leakage RB procedure consists of the following steps:

1.

Choose an initial state $|\rho\rangle\rangle\in\mathcal{H}_{1}$ and measurement projector $|M\rangle\rangle=|\mathbbm{1}_{1}\rangle$ .
2.

For a given $N$ , choose unitaries $U_{0},U_{1},...,U_{N}\in G$ randomly and uniformly. Compute $U_{N+1}=U_{1}^{\dagger}\cdots U_{N}^{\dagger}$ .
3.

Prepare the state $|\rho\rangle\rangle$ . Apply the gates $(U_{1}U_{0}),U_{2},...,U_{N+1}$ sequentially, where $({U}_{1}{U}_{0})$ is compiled as a single element of $G$ .
4.

Perform a measurement of the observable $M$ to determine if the state is still in $\mathcal{H}_{1}$ .
5.

Repeat steps 2-4 many times, to estimate the trivial character-weighted survival probability

$S_{0}(N)=\frac{1}{|G|^{N+1}}\sum_{U_{0},...,U_{N}\in G}\text{Pr}_{U_{0},...,U_{N}}$ (17)

where $\text{Pr}_{U_{0},...,U_{N+1}}$ is the probability of remaining in $\mathcal{H}_{1}$ after applying gates $(U_{1}U_{0}),...,U_{N+1}$ to $|\rho\rangle\rangle$ .
6.

Repeat steps 2-5 for different values of $N$ .
7.

Fit the survival probability to a function of the form

$S_{0}(N)=A\lambda^{N}+B$ (18)

where $A$ , $B$ , and $\lambda$ are independent of $N$ .
8.

Estimate $L$ and $S$ as

$\displaystyle L$ $\displaystyle=(1-B)(1-\lambda)$ (19)

$\displaystyle S$ $\displaystyle=B(1-\lambda)$ (20)
9.

Use the original character RB (section III) to measure the character-weighted survival probability $S_{1\perp}$ associated to the irrep $\mathcal{H}_{1\perp}$ . Fit

$S_{1\perp}(N)=C\lambda_{1\perp}^{N}$

to estimate $\lambda_{1\perp}$ .
10.

Estimate $F_{\Lambda,1}$ as

$F_{\Lambda,1}=\frac{(d_{1}^{2}-1)\lambda_{1\perp}+(d_{1}+1)(1-L)}{d_{1}^{2}+d_{1}}.$ (21)

In the remainder of this section, we prove the correctness of this procedure and provide an example of such leakage RB.

V.1 Deriving $L$ and $S$

Written out explicitly, the zeroth character-weighed survival probability is

S_{0}(N)=\langle\langle\mathbbm{1}_{1}|\hat{\Lambda}\hat{\Lambda}_{G}^{N}\hat{P}_{0}|\rho\rangle\rangle.

where $\hat{P}_{0}$ is the projector onto the trivial irrep, and we have made the same substitutions as in Section III.1 to reduce the sum over $\{U_{0},...,U_{N}\}$ to $G$ -twirls and a projector. We know from Thm. 8 that $\hat{\Lambda}_{G}$ has a block-diagonal form $\hat{\Lambda}_{G}=\bigoplus_{i}\hat{Q}_{i}\otimes\hat{\mathbbm{1}}_{i}$ , where $i$ indexes the irreps. Because $\hat{\Lambda}_{G}$ is multiplied by the projector $\hat{P}_{0}$ in Eq. 17, we may ignore all terms except $\hat{Q}_{0}\otimes\mathbbm{1}_{0}$ . In terms of the eigendecomposition of $\hat{Q}_{0}$ , we may write $\hat{Q}_{0}\otimes\mathbbm{1}_{0}=|e_{0}\rangle\rangle\langle\langle\overline{e}_{0}|+\lambda|e_{1}\rangle\rangle\langle\langle\overline{e}_{1}|$ , so that

{S}_{0}(N)=\langle\langle\mathbbm{1}_{1}|\hat{\Lambda}|e_{0}\rangle\rangle\langle\langle\overline{e}_{0}|\rho\rangle\rangle+\langle\langle\mathbbm{1}_{1}|\hat{\Lambda}|e_{1}\rangle\rangle\langle\langle\overline{e}_{1}|\rho\rangle\rangle\lambda^{N}

where we have used the fact, noted in Section IV, that one eigenvalue of $\hat{Q}_{0}$ is always $1$ . This justifies the fit Eq. 18.

So far, we have simply repeated the steps in Section III.1 with slight modifications. However, in order to estimate $L$ and $S$ we will need to explicitly determine the eigendecomposition of $\hat{Q}_{0}\otimes\mathbbm{1}_{0}$ . We first note that the $\hat{P}_{0}$ subspace is spanned by the orthonormal vectors

\frac{1}{\sqrt{d_{1}}}|\mathbbm{1}_{1}\rangle\rangle:=|\hat{\mathbbm{1}}_{1}\rangle\rangle\qquad\frac{1}{\sqrt{d_{2}}}|\mathbbm{1}_{2}\rangle\rangle:=|\hat{\mathbbm{1}}_{2}\rangle\rangle.

Thus in terms of these basis vectors, we may write

\hat{Q}_{0}\otimes\mathbbm{1}_{0}=|\hat{\mathbbm{1}}_{\alpha}\rangle\rangle Q_{\alpha\beta}\langle\langle\hat{\mathbbm{1}}_{\beta}|

for some constants $Q_{\alpha\beta}$ . Noting that $M_{\alpha\beta}=\langle\langle\hat{\mathbbm{1}}_{\alpha}|\hat{\Lambda}_{G}|\hat{\mathbbm{1}}_{\beta}\rangle\rangle=\langle\langle\hat{\mathbbm{1}}_{\alpha}|\hat{\Lambda}|\hat{\mathbbm{1}}_{\beta}\rangle\rangle$ , we can use the definitions of $L$ and $S$ , (Eqs. 13 and 14) to determine the constants $Q_{\alpha\beta}$ :

Q_{\alpha\beta}=\left(\begin{matrix}1-L&\sqrt{\frac{d_{2}}{d_{1}}}S\\ \sqrt{\frac{d_{1}}{d_{2}}}L&1-S\end{matrix}\right)_{\alpha\beta}.

From the explicit form of $Q_{\alpha\beta}$ , we can determine the eigendecomposition of $\hat{Q}_{0}\otimes\mathbbm{1}_{0}$ via straightforward algebra [26, 23]:

	$\displaystyle\|e_{0}\rangle\rangle$	$\displaystyle=\frac{S}{\sqrt{d_{1}}(L+S)}\|\hat{\mathbbm{1}}_{1}\rangle\rangle+\frac{L}{\sqrt{d_{2}}(L+S)}\|\hat{\mathbbm{1}}_{2}\rangle\rangle$
	$\displaystyle\|\overline{e}_{0}\rangle\rangle$	$\displaystyle=\sqrt{d_{1}}\|\hat{\mathbbm{1}}_{1}\rangle\rangle+\sqrt{d_{2}}\|\hat{\mathbbm{1}}_{2}\rangle\rangle$
	$\displaystyle\|e_{1}\rangle\rangle$	$\displaystyle=\sqrt{d_{2}}\|\hat{\mathbbm{1}}_{1}\rangle\rangle-\sqrt{d_{1}}\|\hat{\mathbbm{1}}_{2}\rangle\rangle$
	$\displaystyle\|\overline{e}_{1}\rangle\rangle$	$\displaystyle=\frac{L}{\sqrt{d_{2}}(L+S)}\|\hat{\mathbbm{1}}_{1}\rangle\rangle-\frac{S}{\sqrt{d_{1}}(L+S)}\|\hat{\mathbbm{1}}_{2}\rangle\rangle$
	$\displaystyle\lambda$	$\displaystyle=1-L-S$

Putting this together, we can evaluate the zeroth character-weighted survival probability as

{S}_{0}(N)=\frac{S}{L+S}+\frac{L}{L+S}(1-L-S)^{N+1}

We then have that $B=\frac{S}{L+S}$ , which can be combined with $\lambda=1-L-S$ to immediately give Eqs. 19 and 20.

V.2 Deriving $F_{\Lambda,1}$

To establish Eq. 21, we first prove the following:

F_{\Lambda,1}=\frac{\operatorname{Tr}(\hat{\Lambda}\hat{P}_{11})+d_{1}(1-L)}{d_{1}^{2}+d_{1}}

(22)

where $\hat{P}_{11}$ is the projector onto $\mathcal{H}_{1}\otimes\mathcal{H}_{1}$ . We use a similar method as in our proof of Eq. 10. We first note that the restricted average fidelities of $\hat{\Lambda}$ and $\hat{P}_{11}\hat{\Lambda}\hat{P}_{11}:=\hat{\Lambda}_{11}$ are equal. $\hat{\Lambda}_{11}$ is an error channel restricted to the $\mathcal{H}_{1}$ subspace. We can twirl $\hat{\Lambda}_{11}$ by the full unitary group on $\mathcal{H}_{1}$ to get a depolarizing channel

(\Lambda_{11})_{G}=p\mathbbm{1}_{1}+q\frac{1}{d_{1}}|\mathbbm{1}_{1}\rangle\rangle\langle\langle\mathbbm{1}_{1}|.

Note that we have $p$ and $q$ rather than $p$ and $(1-p)$ as in Eq. 9; this is because $\hat{\Lambda}_{11}$ is not necessarily trace-preserving. We can directly compute $F_{(\Lambda_{11})_{G}}=p+\frac{q}{d_{1}}$ . Similarly, we can also directly compute $\operatorname{Tr}\left((\hat{\Lambda}_{11})_{G}\right)=pd_{1}^{2}+q$ . Finally, we can directly compute $p+q=\frac{1}{d_{1}}\langle\langle\mathbbm{1}_{1}|(\hat{\Lambda}_{11})_{G}|\mathbbm{1}_{1}\rangle\rangle=\frac{1}{d_{1}}\langle\langle\mathbbm{1}_{1}|\hat{\Lambda}|\mathbbm{1}_{1}\rangle\rangle=1-L$ . Combining these three equations gives Eq. 22.

To estimate $\operatorname{Tr}(\hat{\Lambda}\hat{P}_{11})$ , we can divide this trace up into two pieces:

\operatorname{Tr}(\hat{\Lambda}\hat{P}_{11})=\langle\langle\hat{\mathbbm{1}}_{1}|\hat{\Lambda}|\hat{\mathbbm{1}}_{1}\rangle\rangle+\operatorname{Tr}(\hat{\Lambda}\hat{P}_{1\perp})=(1-L)+\operatorname{Tr}(\hat{\Lambda}\hat{P}_{1\perp})

where $\hat{P}_{1\perp}$ is the projector onto $\mathcal{H}_{1\perp}$ . The latter trace is simply $(d_{1}^{2}-1)\lambda_{1\perp}$ . Plugging this in to Eq. 22 gives Eq. 21 as desired.

V.3 Example: Two-qubit logical encodings

Here, we illustrate the advantages of our leakage RB over the original leakage RB of [26] via a single-qubit example where [26] is not applicable.. We consider an encoding of a single logical qubit into the $S_{z}=0$ subspace of two physical qubits. This encoding is frequently used in quantum dot qubits [60, 61, 62]. The computational space $\mathcal{H}_{1}$ is spanned by

|0\rangle:=\frac{|01\rangle-|10\rangle}{\sqrt{2}},\qquad|1\rangle:=\frac{|01\rangle+|10\rangle}{\sqrt{2}}

and the leakage space $\mathcal{H}_{2}$ is spanned by

|2\rangle:=|00\rangle,\qquad|3\rangle:=|11\rangle.

Let’s assume we implement single-qubit rotations on our computational space by the operators

R_{X}=X_{C}\oplus Z_{L}\qquad R_{Z}=Z_{C}\oplus\frac{X_{L}+Z_{L}}{\sqrt{2}},

where implementing an $X$ or $Z$ rotation on the computational space naturally induces a specific rotation on the leakage space.

We will take our benchmarking group to be the group generated by these two rotations, $G=\langle R_{X},R_{Z}\rangle$ . This group has a total of 16 elements. It cannot be written as direct sum of a group acting on $\mathcal{H}_{1}$ and a group acting on $\mathcal{H}_{2}$ as in Eq. 16, so the leakage RB of [26] does not apply. However, elementary calculation shows that the natural representation of this group contains exactly two trivial irreps, spanned by $|\mathbbm{1}_{1}\rangle\rangle$ and $|\mathbbm{1}_{2}\rangle\rangle$ , and we can therefore use our procedure to estimate $L$ and $S$ .

We illustrate this method by generating random error channels and simulating the RB procedure. In Figs. 5, we show the exact value of $S_{0}(N)$ , the data we take to estimate $S_{0}(N)$ , and the fit to $S_{0}(N)$ according to Eq. 18. In Fig. 6, we repeat the same fitting procedure for a set of randomly generated error channels, and estimate $L$ and $S$ using Eq. 19. We see that the true values of $L$ and $S$ and our estimate for $L$ and $S$ agree within the error bars set by the uncertainty in our fits.

We cannot apply our method to find $F_{\Lambda,1}$ because in this example $\mathcal{H}_{2\perp}$ and $\mathcal{H}_{1\perp}$ share an irrep. This reflects the overall difficulty in applying leakage RB to physically realistic circumstances. While this work provides the most widely applicable method for leakage RB currently available, more work is needed to develop a truly general procedure.

VI Application: Matchgate RB

We can also use our method to introduce a new procedure for scalably benchmarking circuits made of matchgates. Matchgates are 2-qubit gates of the form

G(A,B)=\left(\begin{matrix}a_{11}&0&0&a_{12}\\ 0&b_{11}&b_{12}&0\\ 0&b_{21}&b_{22}&0\\ a_{21}&0&0&a_{22}\end{matrix}\right)

with $\det(A)=\det(B)$ . In other words, a matchgate acts as $A$ on the even parity subspace spanned by $\{|00\rangle,|11\rangle\}$ and as $B$ on the odd parity subspace spanned by $\{|01\rangle,|10\rangle\}$ . Without loss of generality we may assume $\det(A)=\det(B)=1$ . The set of matchgates acting on a line of nearest neighbors is efficiently simulable [27, 29, 30, 28]. However matchgates acting on next-nearest-neighbors [30] or acting on any nontrivial connectivity graph [34, 31] are universal, as are matchgates plus arbitrary one-qubit gates [32, 29], matchgates plus a single $G(A,B)$ with $\det(A)\neq\det(B)$ [33], matchgates acting on entangled input states [35], and matchgates plus adaptive measurements [35]. Implementations of arbitrary matchgates have been proposed for trapped atom systems [68] and have been experimentally demonstrated in photonic systems [69].

We will derive a benchmarking procedure that determines the average fidelity of circuits composed of matchgates using a number of experiments that scales polynomially in the number of qubits. Our method is the matchgate equivalent of traditional Clifford RB, which characterizes the average fidelity of circuits composed of Hadamard, phase, and CNOT gates, and also requires a number of experiments that scales polynomially in the number of qubits. However, we will see that benchmarking matchgate circuits requires the full machinery of non-multiplicity-free character RB.

VI.1 The matchgate group

Consider a line of $n$ qubits with nearest-neighbor connectivity. Let $G$ be the matchgate group on $n$ qubits, the group of all unitaries generated from nearest-neighbor matchgates. Naively, $G$ could contain arbitrarily long circuits of matchgates. However, one can prove that every element of $G$ can be realized using circuits of at most $4n^{3}$ nearest-neighbor matchgates [30, Thm. 5]. We will provide a simplified proof of this fact below.

Following [29, 30], our primary tool to understand $G$ will be the Jordan-Wigner transformation [70]. Define $2n$ Majorana operators $\{c_{i}\}$ as

	$\displaystyle c_{2k-1}$	$\displaystyle=Z_{1}\cdots Z_{k-1}X_{k}$
	$\displaystyle c_{2k}$	$\displaystyle=Z_{1}\cdots Z_{k-1}Y_{k}$

for $k=1,...,n$ . The $\{c_{m}\}$ are Hermitian operators satisfying $\{c_{\ell},c_{m}\}=2\delta_{\ell m}$ . Polynomials in the $\{c_{m}\}$ form a Hermitian basis for the space of all density matrices, so a unitary $U$ is defined by its action on the $\{c_{m}\}$ up to a potential phase. Because of our restriction $\det(A)=\det(B)=1$ , there is no phase freedom on the matchgates or any product of matchgates, so the action of $U\in G$ is entirely determined by its action on the $\{c_{m}\}$ . We make two claims [30]:

Claim 1.

Every $U\in G$ in the matchgate group acts on the Majorana operators as a proper rotation. In other words, there exists some $R\in SO(2n)$ such that $Uc_{\ell}U^{\dagger}=R_{\ell m}c_{m}$ .

Claim 2.

Any unitary operator $U\in U(2^{n})$ that acts on the Majorana operators as a proper rotation is in the matchgate group $G$ . In particular, such a $U$ can be decomposed into a product of at most $2n^{3}$ nearest-neighbor matchgates.

These two claims together imply that the matchgate group is isomorphic to $SO(2N)$ , and that every element of the matchgate group can be efficiently implemented in a quantum circuit. In particular, this shows that the matchgate group is a compact group, thus we can apply character RB.

VI.1.1 Proof of claims

Proof of Claim 1.

We provide a simplification of the proof in [30]. We prove that a nearest-neighbor matchgate acting on qubits $k$ and $k+1$ acts as a rotation mixing $c_{2k-1}$ , $c_{2k}$ , $c_{2k+1}$ , and $c_{2k+2}$ , and that all such rotations are realized by matchgates. It then follows that all products of matchgates also act as rotations on the Majorana operators.

Without loss of generality, we can restrict ourselves to $k=1$ , so our Majorana operators are given by

	$\displaystyle c_{1}$	$\displaystyle=X_{1}$	$\displaystyle c_{3}$	$\displaystyle=Z_{1}X_{2}$
	$\displaystyle c_{2}$	$\displaystyle=Y_{2}$	$\displaystyle c_{4}$	$\displaystyle=Z_{1}Y_{2}.$

We can write an infinitesimal matchgate as $U=\mathbbm{1}-i\epsilon M$ , where $M$ must be of the form

\alpha_{12}Z_{1}-\alpha_{13}Y_{1}X_{2}-\alpha_{14}Y_{1}Y_{2}+\alpha_{23}X_{1}X_{2}+\alpha_{24}X_{1}Y_{2}+\alpha_{34}Z_{2}

with $\alpha_{ab}\in\mathbbm{R}$ . One can directly check that $U$ satisfies

	$\displaystyle Uc_{1}U^{\dagger}$	$\displaystyle=c_{1}+2\epsilon\alpha_{12}c_{2}+2\epsilon\alpha_{13}c_{3}+2\epsilon\alpha_{14}c_{4}$
	$\displaystyle Uc_{2}U^{\dagger}$	$\displaystyle=-2\epsilon\alpha_{12}c_{1}+c_{2}+2\epsilon\alpha_{23}c_{3}+2\epsilon\alpha_{24}c_{4}$
	$\displaystyle Uc_{3}U^{\dagger}$	$\displaystyle=-2\epsilon\alpha_{13}c_{1}-2\epsilon\alpha_{23}c_{2}+c_{3}+2\epsilon\alpha_{34}c_{4}$
	$\displaystyle Uc_{4}U^{\dagger}$	$\displaystyle=-2\epsilon\alpha_{14}c_{1}-2\epsilon\alpha_{24}c_{2}-2\epsilon\alpha_{34}c_{3}+c_{4}$

so that $Uc_{i}U^{\dagger}=R_{ij}c_{j}$ with

R=\mathbbm{1}+2\epsilon\left(\begin{matrix}0&\alpha_{12}&\alpha_{13}&\alpha_{14}\\ -\alpha_{12}&0&\alpha_{23}&\alpha_{24}\\ -\alpha_{13}&-\alpha_{23}&0&\alpha_{34}\\ -\alpha_{14}&-\alpha_{24}&-\alpha_{34}&0\end{matrix}\right)

We therefore see that infinitesimal matchgates generate the whole Lie algebra $\mathfrak{so}(4)$ of real antisymmetric matrices. By exponentiating the infinitesimal matchgates, we generate the full set of matchgates; in this process, we generate the full group $SO(4)$ as well. ∎

Proof of Claim 2.

We note, following [30], that every $R\in SO(2n)$ can be decomposed into $n(2n-1)$ rotations that act as the identity on all but $2$ basis elements $c_{\ell},c_{m}$ by the Hoffman algorithm [71, 72]. In turn, a rotation mixing $c_{\ell}$ and $c_{m}$ with $\ell<m$ can be decomposed into a product of $s:=\left(\lceil\frac{m}{2}\rceil-\lceil\frac{\ell}{2}\rceil-1\right)$ rotations that exchange $(c_{\ell}\leftrightarrow c_{\ell+2})$ , $(c_{\ell+2}\leftrightarrow c_{\ell+4})$ , …, $(c_{\ell+2s-2}\leftrightarrow c_{\ell+2s})$ , followed by a rotation that mixes $c_{\ell+2s}$ and $c_{m}$ , followed by $s$ rotations that exchange $(c_{\ell+2s}\leftrightarrow c_{\ell+2s-2})$ , $(c_{\ell+2s-2}\leftrightarrow c_{\ell+2n-4})$ , …, $(c_{\ell+2}\leftrightarrow c_{\ell})$ . Each of these rotations only involve Majorana operators associated to neighboring qubits, and thus can be written as a matchgate. Thus, $R$ can be realized as the product of a total of $n(2n-1)(2s+1)<4n^{3}$ matchgates, as claimed. ∎

We note that an arbitrary rotation between two Majorana operators

\left(\begin{matrix}c_{\ell}\\ c_{m}\end{matrix}\right)\rightarrow\left(\begin{matrix}\cos(\theta)&\sin(\theta)\\ -\sin(\theta)&\cos(\theta)\end{matrix}\right)\left(\begin{matrix}c_{\ell}\\ c_{m}\end{matrix}\right)

is generated by the unitary $U=e^{\frac{\theta}{2}c_{\ell}c_{m}}$ . In the case where $\left|\lceil\frac{m}{2}\rceil-\lceil\frac{\ell}{2}\rceil\right|\leq 1$ , this $U$ is a nearest-neighbor matchgate. For example, if $\ell=3$ , $m=5$ , then we have $U=e^{-i\frac{\theta}{2}Y_{2}X_{3}}$ . Thus, the above decomposition of $R$ into $<4n^{3}$ two-Majorana rotations gives an explicit formula for the matchgates needed to construct $R$ . We provide Python code to realize the Hoffman decomposition of $R$ into elementary rotations, as well as the reduction of $R$ to a matchgate circuit, at [73].

VI.2 Irreps of the matchgate group

We want to understand how the natural representation of $G$ decomposes into irreps. This is most convenient in the basis of polynomials of $\{c_{m}\}$ . Note that $c_{m}^{2}=1$ , so our polynomials are at most degree $1$ in any given $c_{m}$ and there are $4^{N}$ such polynomials. Explicitly, an orthonormal basis of $\mathcal{H}\otimes\mathcal{H}$ is given by

$\displaystyle\frac{1}{2^{N/2}}\mathbbm{1}:=$	$\displaystyle\|\hat{\mathbbm{1}}\rangle\rangle$
$\displaystyle\frac{1}{2^{N/2}}c_{m_{1}}:=$	$\displaystyle\|m_{1}\rangle\rangle$	$\displaystyle 1\leq m_{1}\leq 2n$
$\displaystyle\frac{1}{2^{N/2}}c_{m_{1}}c_{m_{2}}:=$	$\displaystyle\|m_{1}m_{2}\rangle\rangle$	$\displaystyle 1\leq m_{1}<m_{2}\leq 2n$
	$\displaystyle\vdots$	$\displaystyle\vdots\qquad$
$\displaystyle\|m_{1}\cdots m_{2n-1}\rangle\rangle$		$\displaystyle 1\leq m_{1}<\cdots\leq 2n$
$\displaystyle\|1\cdots 2n\rangle\rangle.$

Define $\mathcal{H}_{i}:=\text{span}\{|m_{1}\cdots m_{i}\rangle\rangle\}$ to be the space spanned by degree- $i$ basis elements, for each $i=0,...,2n$ . Then $\mathcal{H}_{i}\simeq\bigwedge^{i}\mathbbm{C}^{2n}$ , the $i$ -fold wedge product of $\mathbbm{C}^{2n}$ . It’s clear that $\hat{U}$ preserves each $\mathcal{H}_{i}$ , so that each $\mathcal{H}_{i}$ is a subrepresentation. On $\mathcal{H}_{1}$ , $\hat{U}$ acts as the rotation operator $R$ associated to $U$ :

\hat{U}|i_{1}\rangle\rangle=R_{i_{1}j_{1}}|j_{1}\rangle\rangle.

On general $\mathcal{H}_{i}$ , $\hat{U}$ acts as the wedge product of the rotation operator:

\displaystyle\hat{U}|\ell_{1}\cdots\ell_{i}\rangle\rangle=\sum_{\mathclap{\begin{subarray}{c}m_{1}<\cdots<m_{i}\\ \sigma\in S^{i}\end{subarray}}}(-1)^{\sigma}R_{\ell_{1}m_{\sigma 1}}\cdots R_{\ell_{i}m_{\sigma i}}|m_{1}\cdots m_{i}\rangle\rangle.

Claim 3.

The natural representation of the matchgate group decomposes into the irreps

\mathcal{H}_{0}\oplus\mathcal{H}_{1}\oplus\cdots\oplus\mathcal{H}_{n,1}\oplus\mathcal{H}_{n,2}\oplus\cdots\oplus\mathcal{H}_{2n-1}\oplus\mathcal{H}_{2n}.

where $\mathcal{H}_{n}=\mathcal{H}_{n,1}\oplus\mathcal{H}_{n,2}$ . Explicitly, we have

	$\displaystyle\mathcal{H}_{n,1}$	$\displaystyle=\text{span}\{\|\ell_{1}\cdots\ell_{n}\rangle\rangle+i^{n}(-1)^{\sigma(\ell,m)}\|m_{1}\cdots m_{n}\rangle\rangle\}$
	$\displaystyle\mathcal{H}_{n,2}$	$\displaystyle=\text{span}\{\|\ell_{1}\cdots\ell_{n}\rangle\rangle-i^{n}(-1)^{\sigma(\ell,m)}\|m_{1}\cdots m_{n}\rangle\rangle\}$

where $\{m_{a}\}$ is the complement of $\{\ell_{a}\}$ and $\sigma(\ell,m)$ is the permutation that takes $(\ell_{1},...,\ell_{n},m_{1},...,m_{n})\mapsto(1,...,2n)$ . Note that if $n$ is even these are real representations, while for $n$ odd these representations are complex conjugates of each other. The irreps $\mathcal{H}_{i}$ and $\mathcal{H}_{2n-i}$ are isomorphic for $i\neq n$ , but no other irreps are isomorphic to each other.

Proof.

Define the Hodge star operator $*:\mathcal{H}_{i}\rightarrow\mathcal{H}_{2n-i}$ by

*|\ell_{1}\cdots\ell_{i}\rangle\rangle=(-1)^{\sigma(\ell,m)}|m_{1}\cdots m_{2n-i}\rangle\rangle

where $\{m_{a}\}$ is the complement of $\{\ell_{a}\}$ and $\sigma(\ell,m)$ is the permutation that takes $(\ell_{1},...,\ell_{i},m_{1},...,m_{2n-i})\mapsto(1,...,2n)$ . It is straightforward to show that $*$ commutes with the action of $U$ , and thus provides the isomorphism of representations $\mathcal{H}_{i}\simeq\mathcal{H}_{2n-i}$ when $i\neq n$ . We defer the proof that the $\mathcal{H}_{i}$ , $\mathcal{H}_{n,1}$ , and $\mathcal{H}_{n,2}$ are in fact irreducible to chapter 4 of [74]. ∎

VI.3 Benchmarking the matchgate group

Let $\overline{G}\subset G$ be the subgroup of the matchgate group generated by $R\in SO(2n)$ with $R$ diagonal. Such an $R$ is always of the form $R=\text{diag}\{\sigma_{1},...,\sigma_{2n}\}$ with $\sigma_{1}\sigma_{2}\cdots\sigma_{2n}=1$ . The action on a state $|m_{1}\cdots m_{i}\rangle\rangle\in\mathcal{H}_{i}$ is given by

\hat{U}|m_{1}\cdots m_{i}\rangle\rangle=\sigma_{i_{1}}\cdots\sigma_{i_{m}}|m_{1}\cdots m_{i}\rangle\rangle

and therefore the states $|i_{1}\cdots i_{m}\rangle\rangle$ are the irreps of the natural representation of $\overline{G}$ . Because of the constraint $\sigma_{1}\sigma_{2}\cdots\sigma_{2N}=1$ , each irrep has multiplicity 2, with the irrep spanned by $|m_{1}\cdots m_{i}\rangle\rangle$ isomorphic to the irrep spanned by $|\ell_{1}\cdots\ell_{2n-i}\rangle\rangle$ with $\{\ell_{a}\}$ the complement of $\{m_{a}\}$ . For each $i=0,...,n$ , we can define a character function and corresponding projector

	$\displaystyle\chi_{\overline{i}}(R)=$	$\displaystyle\sigma_{1}\cdots\sigma_{i}$
	$\displaystyle\hat{P}_{\overline{i}}=$	$\displaystyle\|1\cdots i\rangle\rangle\langle\langle 1\cdots i\|$
		$\displaystyle+\|(i+1)\cdots 2n\rangle\rangle\langle\langle(i+1)\cdots 2n\|.$

These projectors project into the multiplicty-two irreps $\mathcal{H}_{i}\oplus\mathcal{H}_{2n-i}$ for $\overline{i}=0,...,(n-1)$ , and project into the two inequivalent irreps $\mathcal{H}_{n,1}\oplus\mathcal{H}_{n,2}$ for $\overline{i}=n$ .

As our initial state, for each $i=0,...,n$ we choose

|\rho_{i}\rangle\rangle=\left\{\begin{array}[]{ll}|0\cdots+\cdots 0\rangle\rangle,&i=2k-1\\ |0\cdots 0\rangle\rangle,&i=2k.\end{array}\right.

where $k$ th qubit is in the $+$ state of the $X$ operator for $i=2k-1$ . Provided we can prepare both $X$ -basis and $Z$ -basis single qubit states, we can prepare $|\rho_{i}\rangle\rangle$ .

As our measurement projector, for each $i=0,...,n$ we choose

|M_{i}\rangle\rangle=\left\{\begin{array}[]{ll}\frac{1}{2}(X_{k}+\mathbbm{1}),&i=2k-1\\ \frac{1}{2}\left(\prod_{\alpha>n-k}Z_{\alpha}+\mathbbm{1}\right),&i=2k.\end{array}\right.

For $i=2k-1$ , this corresponds to a measurement of the $k$ th qubit in the $X$ basis, while for $i=2k$ this corresponds to a measurement of the product of the last $k$ qubits in the $Z$ basis.

With these choices, the $S_{i}(N)$ are approximately

S_{i}(N)\approx\frac{\langle\langle M_{i}|\hat{P}_{\overline{i}}|\rho_{i}\rangle\rangle}{\text{dim}(\mathcal{\overline{H}}_{\overline{i}})}=\left\{\begin{array}[]{ll}1,&i=0\\ \frac{1}{2},&1\leq i\leq n\end{array}\right.

and the relative uncertainty does not depend on the number of qubits. This is therefore a scalable method to benchmark the matchgate group.

The form of the decay is given by

S_{i}(N)=\left\{\begin{array}[]{ll}C_{0}\lambda_{0}^{N}+B,&i=0\\ C_{i,1}\lambda_{i,1}^{N}+C_{i,2}\lambda_{i,2}^{N},&1\leq i\leq n.\end{array}\right.

(23)

For each $i$ , either $\lambda_{i,1},\lambda_{i,2},C_{i,1},C_{i,2}\in\mathbbm{R}$ or $\lambda_{i,1}=\lambda_{i,2}^{*}$ and $C_{i,1}=C_{i,2}^{*}$ , since $S_{i}(N)$ is always real. For the case of $i=n$ , we know that the former case holds when $n$ is even and the latter when $n$ is odd, by Claim 3. For $1\leq i<n$ , one should assume whichever case gives the best fit. Note that in all cases, we fit at most $4$ real parameters.

As an example, we simulate a noisy implementation of the matchgate group on $n=3$ qubits. In Fig. 7, we show the exact value of $S_{i}(N)$ , the data we take to estimate $S_{i}(N)$ , and the fit to $S_{i}(N)$ according to Eq. 23 for a single random error channel $\Lambda$ . In Fig. 8, we do the same fitting procedure for a set of randomly generated error channels, and estimate their fidelity. We see that the true fidelity and the estimated fidelity agree within the error bars set by the uncertainty of our fits.

VII Conclusion and Discussions

In this work, we extended the recently introduced character RB of [22] to groups with multiplicity. Compared to earlier work on benchmarking arbitrary groups [20, 21], our method allows us to accurately determine the fidelity and fit fewer exponentials to experimental data. The generalization to non-multiplicity-free groups was essential to deriving a rigorous version of subspace RB and a scalable RB protocol for the matchgate group. This generalization also allowed us to develop an improved leakage RB protocol.

While we derived the character RB procedure in more generality than [22], our generalization still requires groups of small multiplicity, since the multiplicity of the group determines the number of exponential decays in our fit function. Robustly fitting a sum of many exponential decays is challenging, especially when the decay rates are roughly equal [43, 44]. It is likely straightforward to benchmark groups in which the trivial irrep has multiplicity three, as the corresponding decay $S_{0}(N)=A+B\lambda_{0,1}^{N}+C\lambda_{0,2}^{N}$ has only five real parameters. An irrep of multiplicity three with a real character function $\chi$ has a decay with six parameters, which may be feasible with sufficient data. A general irrep of multiplicity three, however, requires fitting nine real parameters, which is likely unfeasible for realistic amounts of data. Higher-multiplicity irreps are correspondingly more difficult. All of the groups we considered in the examples in this paper decomposed into irreps with multiplicity at most 2.

All our applications involved a group that preserved some subspace of the Hilbert space. In the case of subspace RB, the group preserved the triplet and singlet subspaces; in the case of leakage RB, the computational and leakage subspaces; and in the case of matchgate RB, the even and odd parity subspaces. Any group that preserves subspaces necessarily has multiplicity, since there is always a copy of the trivial irrep in each subspace. It is an open question whether non-multiplicity-free character RB has useful applications to groups that do not preserve subspaces but nonetheless have multiplicity.

One group related to the matchgate group that would be of immediate experimental interest is the XY group, the subgroup of the matchgate group generated by only nearest-neighbor XY mixers $U_{XY}(\theta)=\exp\left\{i\theta(X_{1}\otimes X_{2}+Y_{1}\otimes Y_{2})\right\}$ . Unlike general matchgates, XY mixers can be naturally realized on superconducting qubits [75, 76], and they are a necessary ingredient in extensions of the QAOA algorithm [77, 78, 79]. In addition, XY mixers are efficiently simulable on a line but become universal on nontrivial graphs, just like the full matchgate group[31]. However, XY mixers on $N$ qubits preserve the $(N+1)$ subspaces of definite Hamming weight; this implies that the trivial representation of the XY group must have multiplicity $(N+1)$ . Thus, our method cannot be used to scalably benchmark the XY group; even $N=2$ qubits is likely infeasible. On the other hand, [80] recently introduced a compilation of general two-qubit matchgates into products of four XY mixers and single-qubit gates. Using this decomposition, the average fidelity of the resulting two-qubit matchgates can be used as a proxy for the fidelity of the XY mixers. This method is similar to the benchmarking framework in our Sec IV, where we compile group elements into a fixed number of gates of interest (in our case, $U_{ZZ}$ ), with the modification that [80] allows the gate of interest $XY(\theta)$ to vary. It is an open question if there is a generalization of this compilation to the matchgate group on $N>2$ qubits.

While our leakage RB necessitates the fewest assumptions to date, it is still too restrictive for many experimental implementations. Most notably, our RB requires the set of gates to be a group, which may be unrealistic; often, the gates will only form a group modulo rotations in the leakage space. In experimental implementations of leakage RB, this problem is usually simply ignored and an exponential decay is posited to exist with the usual relation to the leakage rate [58, 63]. It is worth exploring whether the methods used here can be further extended to such sets of gates that are only groups in the computational subspace, modulo rotations in the leakage subspace, to provide a more rigorous foundation for leakage RB experiments.

There are two obvious directions for further applications of character RB, with or without multiplicity. First, character RB has the potential to drastically expand the family of groups that can be scalably benchmarked. This requires both finding a group $G$ that can be efficiently compiled into elementary gates whose multiplicity is bounded as the number of qubits $n$ increases, as well as finding a subgroup $\overline{G}\subseteq G$ whose irreps have slowly growing dimension. As a simple example, the subgroups of the Clifford group considered in [20] likely have a scalable protocol based on character RB, with $\overline{G}$ given by the Pauli group. Increasing the number of groups that can be scalably benchmarked gives new ways of characterizing compiled gates, especially non-Clifford gates.

Second, character RB can be used to characterize specific elementary gates by combining these gates into a group, as we did in Section IV for subspace RB. This requires finding a group that can be implemented by combining a fixed number of the gate to be characterized with known high-fidelity gates. Constructing these groups is a non-trivial task, as we have seen in the case of the $U_{ZZ}$ operator above. We leave the exploration of such applications to future work.

Note added. After the first version of this paper was posted, [80] was posted to the arXiv which also proposes a matchgate RB. Their method relies on enlarging the matchgate group with additional unitaries to avoid representations with multiplicity, but is otherwise similar to ours. As we mentioned in this paper, our character RB does not apply to the group generated by nearest-neighbor XY gates. While [80] does not propose a method to benchmark the group generated by nearest-neighbor XY mixers, they do demonstrate a method to compile two-qubit matchgate elements using a fixed number of XY mixers and additional single-qubit gates, allowing the matchgate RB to be used to characterize XY mixers, as discussed above.

Acknowledgements

JC thanks Alexandre Pyvovarov for useful discussions on representations of $\bigwedge^{i}\mathbbm{C}^{2N}$ . We are grateful for support from NASA Ames Research Center, the NASA Advanced Exploration systems (AES) program, and the NASA Transformative Aeronautic Concepts Program (TACP). We are also grateful for support from the AFRL Information Directorate under grant F4HBKC4162G001. JC was supported by the USRA Feynman Quantum Academy funded by the NAMS R&D Student Program at NASA Ames Research Center. JC and ZW are also supported by NASA Academic Mission Services, Contract No. NNA16BD14C.

References

O’Brien et al. [2004] J. L. O’Brien, G. J. Pryde, A. Gilchrist, D. F. V. James, N. K. Langford, T. C. Ralph, and A. G. White, Quantum process tomography of a controlled-NOT gate, Phys. Rev. Lett. 93, 080502 (2004).
Chuang and Nielsen [1997] I. L. Chuang and M. A. Nielsen, Prescription for experimental determination of the dynamics of a quantum black box, J. Mod. Opt 44, 2455 (1997).
Poyatos et al. [1997] J. F. Poyatos, J. I. Cirac, and P. Zoller, Complete characterization of a quantum process: the two-bit quantum gate, Phys. Rev. Lett. 78, 390 (1997).
Mohseni et al. [2008] M. Mohseni, A. T. Rezakhani, and D. Lidar, Quantum-process tomography: Resource analysis of different strategies, Phys. Rev. A 77, 032322 (2008).
Emerson et al. [2005] J. Emerson, R. Alicki, and K. Życzkowski, Scalable noise estimation with random unitary operators, J. Opt. B 7, S347 (2005).
Knill et al. [2008] E. Knill, D. Leibfried, R. Reichle, J. Britton, R. B. Blakestad, J. D. Jost, C. Langer, R. Ozeri, S. Seidelin, and D. J. Wineland, Randomized benchmarking of quantum gates, Phys. Rev. A 77, 012307 (2008).
Magesan et al. [2011] E. Magesan, J. M. Gambetta, and J. Emerson, Scalable and robust randomized benchmarking of quantum processes, Phys. Rev. Lett. 106, 180504 (2011).
Magesan et al. [2012a] E. Magesan, J. M. Gambetta, and J. Emerson, Characterizing quantum gates via randomized benchmarking, Phys. Rev. A 85, 042311 (2012a).
Dankert et al. [2009] C. Dankert, R. Cleve, J. Emerson, and E. Livine, Exact and approximate unitary 2-designs and their application to fidelity estimation, Phys. Rev. A 80, 012304 (2009).
Gottesman [1998a] D. Gottesman, Theory of fault-tolerant quantum computation, Phys. Rev. A 57, 127 (1998a).
Gottesman [1998b] D. Gottesman, The Heisenberg representation of quantum computers, quant-ph/9807006 (1998b).
Aaronson and Gottesman [2004] S. Aaronson and D. Gottesman, Improved simulation of stabilizer circuits, Phys. Rev. A 70, 052328 (2004).
Nielsen [2002] M. A. Nielsen, A simple formula for the average gate fidelity of a quantum dynamical operation, Phys. Lett. A 303, 249 (2002).
Horodecki et al. [1999] M. Horodecki, P. Horodecki, and R. Horodecki, General teleportation channel, singlet fraction, and quasidistillation, Phys. Rev. A 60, 1888 (1999).
Sanders et al. [2015] Y. R. Sanders, J. J. Wallman, and B. C. Sanders, Bounding quantum gate error rate based on reported average fidelity, New J. Phys. 18, 012002 (2015).
Magesan et al. [2012b] E. Magesan, J. M. Gambetta, B. R. Johnson, C. A. Ryan, J. M. Chow, S. T. Merkel, M. P. Da Silva, G. A. Keefe, M. B. Rothwell, T. A. Ohki, et al., Efficient measurement of quantum gate error by interleaved randomized benchmarking, Phys. Rev. Lett. 109, 080505 (2012b).
Harper and Flammia [2017] R. Harper and S. T. Flammia, Estimating the fidelity of T gates using standard interleaved randomized benchmarking, Quantum Sci. Technol. 2, 015008 (2017).
Carignan-Dugas et al. [2015] A. Carignan-Dugas, J. J. Wallman, and J. Emerson, Characterizing universal gate sets via dihedral benchmarking, Phys. Rev. A 92, 060302 (2015).
Cross et al. [2016] A. W. Cross, E. Magesan, L. S. Bishop, J. A. Smolin, and J. M. Gambetta, Scalable randomised benchmarking of non-Clifford gates, npj Quantum Inf. 2, 1 (2016).
Brown and Eastin [2018] W. G. Brown and B. Eastin, Randomized benchmarking with restricted gate sets, Phys. Rev. A 97, 062323 (2018).
França and Hashagen [2018] D. S. França and A. K. Hashagen, Approximate randomized benchmarking for finite groups, J. Phys. A 51, 395302 (2018).
Helsen et al. [2019] J. Helsen, X. Xue, L. M. K. Vandersypen, and S. Wehner, A new class of efficient randomized benchmarking protocols, npj Quantum Inf. 5, 1 (2019).
Baldwin et al. [2020] C. H. Baldwin, B. J. Bjork, J. P. Gaebler, D. Hayes, and D. Stack, Subspace benchmarking high-fidelity entangling operations with trapped ions, Phys. Rev. Res. 2, 013317 (2020).
Chasseur and Wilhelm [2015] T. Chasseur and F. K. Wilhelm, Complete randomized benchmarking protocol accounting for leakage errors, Phys. Rev. A 92, 042333 (2015).
Wallman et al. [2016] J. J. Wallman, M. Barnhill, and J. Emerson, Robust characterization of leakage errors, New J. Phys. 18, 043021 (2016).
Wood and Gambetta [2018] C. J. Wood and J. M. Gambetta, Quantification and characterization of leakage errors, Phys. Rev. A 97, 032306 (2018).
Valiant [2002] L. G. Valiant, Quantum circuits that can be simulated classically in polynomial time, SIAM J. Comput. 31, 1229 (2002).
Brod [2016] D. J. Brod, Efficient classical simulation of matchgate circuits with generalized inputs and measurements, Phys. Rev. A 93, 062332(R) (2016).
Terhal and DiVincenzo [2002] B. M. Terhal and D. P. DiVincenzo, Classical simulation of noninteracting-fermion quantum circuits, Phys. Rev. A 65, 032325 (2002).
Jozsa and Miyake [2008] R. Jozsa and A. Miyake, Matchgates and classical simulation of quantum circuits, Proc. Math. Phys. Eng. Sci. 464, 3089 (2008).
Brod and Childs [2014] D. J. Brod and A. M. Childs, The computational power of matchgates and the XY interaction on arbitrary graphs, Quantum Inf. Comput. 14, 901 (2014).
Imamoglu et al. [1999] A. Imamoglu, D. D. Awschalom, G. Burkard, D. P. DiVincenzo, D. Loss, M. Sherwin, A. Small, et al., Quantum information processing using quantum dot spins and cavity QED, Phys. Rev. Lett. 83, 4204 (1999).
Brod and Galvao [2011] D. J. Brod and E. F. Galvao, Extending matchgates into universal quantum computation, Phys. Rev. A 84, 022310 (2011).
Brod and Galvao [2012] D. J. Brod and E. F. Galvao, Geometries for universal quantum computation with matchgates, Phys. Rev. A 86, 052307 (2012).
Hebenstreit et al. [2020] M. Hebenstreit, R. Jozsa, B. Kraus, and S. Strelchuk, Computational power of matchgates with supplementary resources, arXiv:2007.08231 (2020).
Wood et al. [2015] C. J. Wood, J. D. Biamonte, and D. G. Cory, Tensor networks and graphical calculus for open quantum systems, Quantum Inf. Comput. 15, 759 (2015).
Fulton and Harris [2013] W. Fulton and J. Harris, Representation theory: a first course, Vol. 129 (Springer Science & Business Media, 2013).
Proctor et al. [2017] T. Proctor, K. Rudinger, K. Young, M. Sarovar, and R. Blume-Kohout, What randomized benchmarking actually measures, Phys. Rev. Lett. 119, 130502 (2017).
Wallman [2018] J. J. Wallman, Randomized benchmarking with gate-dependent noise, Quantum 2, 47 (2018).
Merkel et al. [2018] S. T. Merkel, E. J. Pritchett, and B. H. Fong, Randomized benchmarking as convolution: Fourier analysis of gate dependent errors, arXiv preprint arXiv:1804.05951 (2018).
Harper et al. [2019] R. Harper, I. Hincks, C. Ferrie, S. T. Flammia, and J. J. Wallman, Statistical analysis of randomized benchmarking, Phys. Rev. A 99, 052350 (2019).
Note [1] In our notation, this can be seen by noting that the natural representation of the full unitary group decomposes into two irreps which act on $|\mathbbm{1}\rangle\rangle$ and the orthogonal complement of $|\mathbbm{1}\rangle\rangle$ , respectively, and then applying Fact 6.
Bromage [1983] G. Bromage, A quantification of the hazards of fitting sums of exponentials to noisy data, Comput. Phys. Commun. 30, 229 (1983).
Clayden and Hesler [1992] N. Clayden and B. Hesler, Multiexponential analysis of relaxation decays, J. Magn. Reson. 98, 271 (1992).
Istratov and Vyvenko [1999] A. A. Istratov and O. F. Vyvenko, Exponential analysis in physical phenomena, Rev. Sci. Instrum. 70, 1233 (1999).
Holmström and Petersson [2002] K. Holmström and J. Petersson, A review of the parameter estimation problem of fitting positive exponential sums to empirical data, Appl. Math. Comput. 126, 31 (2002).
Hokanson [2013] J. Hokanson, Numerically stable and statistically efficient algorithms for large scale exponential fitting, Ph.D. thesis, Rice University (2013).
Note [2] Ref. [23] claimed it was sufficient to require $G_{2}$ to be a unitary $1$ -design, but this appears to be an error. A similar error was made in [26], from which much of [23] is derived.
Hostens et al. [2005] E. Hostens, J. Dehaene, and B. De Moor, Stabilizer states and Clifford operations for systems of arbitrary dimensions and modular arithmetic, Phys. Rev. A 71, 042315 (2005).
Gottesman [1998c] D. Gottesman, Fault-tolerant quantum computation with higher-dimensional systems, in NASA International Conference on Quantum Computing and Quantum Communications (Springer, 1998) pp. 302–313.
Note [3] Ref. [23] required a shorter circuit of only three entangling gates. However, this circuit cannot implement all relative phases between the subspaces and thus does not result in a group.
Webb [2016] Z. Webb, The Clifford group forms a unitary 3-design, Quantum Inf. Comput. 16, 1379 (2016).
Zhu [2017] H. Zhu, Multiqubit Clifford groups are unitary 3-designs, Phys. Rev. A 96, 062336 (2017).
Chau [2005] H. F. Chau, Unconditionally secure key distribution in higher dimensions by depolarization, IEEE Trans. Inf. Theory 51, 1451 (2005).
Note [4] It follows that $G$ also cannot form a 2-design, as 2-designs are always multiplicity free; in particular, the natural representation of a 2-design decomposes into precisely two non-isomorphic irreps, acting on $|\mathbbm{1}\rangle\rangle$ and the orthogonal complement of $|\mathbbm{1}\rangle\rangle$ [9, 22].
Note [5] Our formula differs slightly from the corresponding formula in [23]. Ref. [23] considered approximating the process (also called entanglement) fidelity rather than the average fidelity; however, the average fidelity can be determined from the process fidelity[14, 13]. To be consistent with the rest of our paper, we have translated their approximation of the process fidelity into the corresponding approximation of the average fidelity.
Gambetta et al. [2017] J. M. Gambetta, J. M. Chow, and M. Steffen, Building logical qubits in a superconducting quantum computing system, npj Quantum Inf. 3, 1 (2017).
Chen et al. [2016] Z. Chen, J. Kelly, C. Quintana, R. Barends, B. Campbell, Y. Chen, B. Chiaro, A. Dunsworth, A. Fowler, E. Lucero, et al., Measuring and suppressing quantum state leakage in a superconducting qubit, Phys. Rev. Lett. 116, 020501 (2016).
DiVincenzo et al. [2000] D. P. DiVincenzo, D. Bacon, J. Kempe, G. Burkard, and K. B. Whaley, Universal quantum computation with the exchange interaction, Nature 408, 339 (2000).
Petta et al. [2005] J. R. Petta, A. C. Johnson, J. M. Taylor, E. A. Laird, A. Yacoby, M. D. Lukin, C. M. Marcus, M. P. Hanson, and A. C. Gossard, Coherent manipulation of coupled electron spins in semiconductor quantum dots, Science 309, 2180 (2005).
Hanson and Burkard [2007] R. Hanson and G. Burkard, Universal set of quantum gates for double-dot spin qubits with fixed interdot coupling, Phys. Rev. Lett. 98, 050502 (2007).
Levy [2002] J. Levy, Universal quantum computation with spin-1/2 pairs and Heisenberg exchange, Phys. Rev. Lett. 89, 147902 (2002).
Andrews et al. [2019] R. W. Andrews, C. Jones, M. D. Reed, A. M. Jones, S. D. Ha, M. P. Jura, J. Kerckhoff, M. Levendorf, S. Meenehan, S. T. Merkel, et al., Quantifying error and leakage in an encoded Si/SiGe triple-dot qubit, Nat. Nanotechnol. 14, 747 (2019).
Häffner et al. [2008] H. Häffner, C. F. Roos, and R. Blatt, Quantum computing with trapped ions, Phys. Rep. 469, 155 (2008).
Hayes et al. [2020] D. Hayes, D. Stack, B. Bjork, A. C. Potter, C. H. Baldwin, and R. P. Stutz, Eliminating leakage errors in hyperfine qubits, Phys. Rev. Lett. 124, 170501 (2020).
Stricker et al. [2020] R. Stricker, D. Vodola, A. Erhard, L. Postler, M. Meth, M. Ringbauer, P. Schindler, T. Monz, M. Müller, and R. Blatt, Experimental deterministic correction of qubit loss, Nature 585, 207 (2020).
Note [6] Ref. [26] originally claimed it was sufficient for $G_{2}$ to be a unitary 1-design, but this appears to be an error.
Herrera et al. [2014] F. Herrera, Y. Cao, S. Kais, and K. B. Whaley, Infrared-dressed entanglement of cold open-shell polar molecules for universal matchgate quantum computing, New J. Phys. 16, 075001 (2014).
Ramelow et al. [2010] S. Ramelow, A. Fedrizzi, A. M. Steinberg, and A. G. White, Matchgate quantum computing and non-local process analysis, New J. Phys. 12, 083027 (2010).
Jordan and Wigner [1928] P. Jordan and E. Wigner, Über das paulische Äquivalenzverbot, Z. Physik 47, 631 (1928).
Raffenetti and Ruedenberg [1969] R. C. Raffenetti and K. Ruedenberg, Parametrization of an orthogonal matrix in terms of generalized Eulerian angles, Int. J. Quantum Chem. 4, 625 (1969).
Hoffman et al. [1972] D. K. Hoffman, R. C. Raffenetti, and K. Ruedenberg, Generalization of Euler angles to N-dimensional orthogonal matrices, J. Math. Phys. 13, 528 (1972).
Claes [2020] J. Claes, Hoffman decomposition and the matchgate group, https://github.com/jahanclaes/Hoffman-Decomposition-and-the-Matchgate-Group (2020).
Knapp [2001] A. W. Knapp, Representation theory of semisimple groups: an overview based on examples, Vol. 36 (Princeton university press, 2001).
Abrams et al. [2019] D. M. Abrams, N. Didier, B. R. Johnson, M. P. da Silva, and C. A. Ryan, Implementation of the XY interaction family with calibration of a single pulse, arXiv:1912.04424 (2019).
Foxen et al. [2020] B. Foxen, C. Neill, A. Dunsworth, P. Roushan, B. Chiaro, A. Megrant, J. Kelly, Z. Chen, K. Satzinger, R. Barends, et al., Demonstrating a continuous set of two-qubit gates for near-term quantum algorithms, Phys. Rev. Lett. 125, 120504 (2020).
Hadfield et al. [2019] S. Hadfield, Z. Wang, B. O’Gorman, E. G. Rieffel, D. Venturelli, and R. Biswas, From the quantum approximate optimization algorithm to a quantum alternating operator ansatz, Algorithms 12, 34 (2019).
Wang et al. [2020] Z. Wang, N. C. Rubin, J. M. Dominy, and E. G. Rieffel, XY mixers: Analytical and numerical results for the quantum alternating operator ansatz, Phys. Rev. A 101, 012320 (2020).
Cook et al. [2019] J. Cook, S. Eidenbenz, and A. Bärtschi, The quantum alternating operator ansatz on maximum k-vertex cover, arXiv:1910.13483 (2019).
Helsen et al. [2020] J. Helsen, S. Nezami, M. Reagor, and M. Walter, Matchgate benchmarking: Scalable benchmarking of a continuous family of many-qubit gates, arXiv preprint arXiv:2011.13048 (2020).
Moore and Russell [2015] C. Moore and A. Russell, Approximate representations, approximate homomorphisms, and low-dimensional embeddings of groups, SIAM J. Discrete Math. 29, 182 (2015).
Gowers and Hatami [2017] W. T. Gowers and O. Hatami, Inverse and stability theorems for approximate representations of finite groups, Sb. Math. 208, 1784 (2017).

Appendix A Gate-dependent errors

In this appendix, we extend the work of [39, 40, 22] on gate-dependent errors to the case of non-multiplicity-free character RB. Ref. [22] had previously generalized [39] to establish that multiplicity-free character RB is robust to gate-dependent errors. Rather than follow the method of [39, 22] we use the Fourier transform method of [40], which is more natural for groups with multiplicity. Our ultimate goal is the following theorem:

Theorem 2.

Let $G$ be a benchmarking group, and let $i$ be an irrep of the natural representation with multiplicity $a_{i}$ . Assume each gate $U\in G$ is realized as a noisy operator $\eta(U)$ , but do not assume we can write $\eta(U)=\hat{\Lambda}\hat{U}$ for some $U$ -independent noise channel $\Lambda$ . Then the character-weighted survival probability is given by

S_{i}(N)=\sum_{j=1}^{a_{i}}C_{i,j}\lambda_{i,j}^{N}+\epsilon_{N}

where $\epsilon_{N}$ is an error term satisfying $|\epsilon_{N}|<\delta_{1}\delta_{2}^{N}$ and $\delta_{1},\delta_{2}$ are both small for high-fidelity gates. Since we know that $\lambda_{i,j}\approx 1$ for high-fidelity gates, $\epsilon_{N}$ is negligible compared to $S_{i}(N)$ for moderately large $N$ .

This theorem implies we may safely use the RB protocols even in the presence of gate-dependent errors, although we will see the interpretation of the estimated fidelity is slightly modified.

In what follows, we will use the notation $\mathbbm{E}\left[\cdot\right]$ for the average $\frac{1}{|G|}\sum_{U\in G}\left(\cdot\right)$ or $\int_{G}dU\left(\cdot\right)$ to make our equations cleaner. We will also use the shorthand $d_{i}$ for $\dim(\mathcal{H}_{i})$ .

A.1 The generalized Fourier transform and its application to character RB

We first define a generalization of the Fourier transform to matrix-valued functions of a group $G$ [81, 82]. For any group $G$ we define $\widetilde{G}$ to index the irreps of $G$ , and we assume WLOG that the irreps are unitary. Given a function $\eta:G\rightarrow\mathcal{L}(\mathbbm{C}^{D})$ , for each $i\in\widetilde{G}$ we define the Fourier transform $\widetilde{\eta}(i)\in\mathcal{L}(\mathbbm{C}^{D})\otimes\mathcal{L}(\mathcal{H}_{i})$ to be

\widetilde{\eta}(i)=\mathbbm{E}\left[\eta(U)\otimes\phi_{i}^{*}(U)\right].

(24)

where $\phi_{i}:G\rightarrow\mathcal{L}(\mathcal{H}_{i})$ is the $i$ th irrep.

Given two matrix-valued functions $\eta,\xi:G\rightarrow\mathcal{L}(\mathbbm{C}^{D})$ , we can also define the convolution $(\eta*\xi)$ by

(\eta*\xi)(U_{0})=\mathbbm{E}\left[\eta(U^{\dagger})\xi(UU_{0})\right].

(25)

The generalized Fourier transform shares many properties with the usual Fourier transform; in particular, we will use the following identities [82, 40]:

$\displaystyle(\widetilde{\eta*\xi})(i)$	$\displaystyle=\widetilde{\eta}(i)\widetilde{\xi}(i)$	(26)
$\displaystyle\mathbbm{E}\left[\text{Tr}\left(\eta(U)\xi^{\dagger}(U)\right)\right]$	$\displaystyle=\sum_{i}d_{i}\text{Tr}\left(\widetilde{\eta}(i)\widetilde{\xi}^{\dagger}(i)\right)$	(27)
$\displaystyle\eta(U)$	$\displaystyle=\sum_{i}d_{i}\text{Tr}_{i}\left([\mathbbm{1}\otimes\phi_{i}^{T}(U)]\widetilde{\eta}(i)\right)$	(28)

where in the last line, $\text{Tr}_{i}\left(\cdot\right)$ is the partial trace over $\mathcal{H}_{i}$ . Eq. 26 is the analogue of the usual convolution identity for Fourier transforms, Eq. 27 is the analogue of Parseval’s identity, and Eq. 28 gives the inverse Fourier transform.

The generalized Fourier transformation is useful because it allows us to express the result of a character RB experiment in a simpler form. A character RB experiment estimates a matrix element of the operator

\hat{O}_{i}:=\mathbbm{E}\left[\eta(U_{1}^{\dagger}\cdots U_{N}^{\dagger})\eta(U_{N})\cdots\eta(U_{2})\eta(U_{1}U_{0})\chi_{\overline{i}}^{*}(U_{0})\right]

where the expectation value is over all $U_{0}\in\overline{G}$ , $U_{1},...,U_{N}\in G$ . Through the change of variables $U_{i}\rightarrow U_{i}U_{i-1}\cdots U_{1}$ for $i=1,...,N$ , we can rewrite this expression as a convolution:

	$\displaystyle\hat{O}_{i}$	$\displaystyle=\mathbbm{E}\left[\eta(U_{N}^{\dagger})\eta(U_{N}U_{N-1}^{\dagger})\cdots\eta(U_{2}U_{1}^{\dagger})\eta(U_{1}U_{0})\chi_{\overline{i}}^{*}(U_{0})\right]$
		$\displaystyle=\mathbbm{E}\left[\underbrace{(\eta\cdots\eta)}_{(N+1)\text{ times}}(U_{0})\chi_{\overline{i}}^{*}(U_{0})\right]$

Using the inverse Fourier transform (Eq. 28) we can write $(\eta*\cdots*\eta)(U_{0})$ in terms of $(\widetilde{\eta*\cdots*\eta})(i^{\prime})$ , while the convolution identity (Eq. 26) allows us to simplify $(\widetilde{\eta*\cdots*\eta})(i^{\prime})=\widetilde{\eta}(i^{\prime})^{N+1}$ . In total, we find

\displaystyle\hat{O}_{i}

\displaystyle=\sum_{i^{\prime}}d_{i^{\prime}}\text{Tr}_{i^{\prime}}\left(\left[\mathbbm{1}\otimes\mathbbm{E}\left[\chi_{\overline{i}}^{*}(U_{0})\phi_{i^{\prime}}(U_{0})\right]^{T}\right]\widetilde{\eta}(i^{\prime})^{N+1}\right).

We now use the projection formula (Fact 7) to note that $d_{\overline{i}}\mathbbm{E}\left[\chi_{\overline{i}}^{*}(U_{0})\phi_{i^{\prime}}(U_{0})\right]$ is just the projection of $\phi_{i^{\prime}}$ onto the irrep $\overline{i}$ of $\overline{G}$ . By assumption, the irrep $\phi_{\overline{i}}$ is a subrepresentation of only $\phi_{i}$ , and not a subrepresentation of any $\phi_{i^{\prime}}$ with $i^{\prime}\neq i$ . Therefore,

\displaystyle\hat{O}_{i}

\displaystyle=\frac{d_{i}}{d_{\overline{i}}}\text{Tr}_{i}\left(\left[\mathbbm{1}\otimes\hat{P}_{\overline{i}}^{T}\right]\widetilde{\eta}(i)^{N+1}\right).

We therefore see that the outcome of a character RB experiment, $S_{i}(N)$ , can be described by the Fourier transform of $\eta$ via

	$\displaystyle S_{i}(N)$	$\displaystyle=\langle\langle M_{i}\|\hat{\Lambda}_{M}\hat{O}_{i}\hat{\Lambda}_{P}\|\rho_{i}\rangle\rangle$		(29)
		$\displaystyle=\frac{d_{i}}{d_{\overline{i}}}\langle\langle M_{i}\|\hat{\Lambda}_{M}\text{Tr}_{i}\left(\left[\mathbbm{1}\otimes\hat{P}^{T}_{\overline{i}}\right]\widetilde{\eta}(i)^{N+1}\right)\hat{\Lambda}_{P}\|\rho_{i}\rangle\rangle$

and the decay of $S_{i}(N)$ is determined by the eigenvalues of $\widetilde{\eta}(i)$ .

A.2 Simplifying the decay

In the case of ideal gates $\eta_{\text{ideal}}(U)=\hat{U}$ , we have that $\widetilde{\eta}_{\text{ideal}}(i)$ is given by

\widetilde{\eta}_{\text{ideal}}(i)=\mathbbm{E}\left[\hat{U}\otimes\phi_{i}(U)\right]

This can be simplified by noting that the map $\eta_{\text{ideal}}\otimes\phi_{i}:U\mapsto\hat{U}\otimes\phi_{i}(U)$ is a representation of $G$ , and $\mathbbm{E}\left[\hat{U}\otimes\phi_{i}(U)\right]$ is the projection of this representation onto the copies of the trivial irrep (Fact 7). We can count the multiplicity of the trivial irrep in $(\eta_{\text{ideal}}\otimes\phi_{i})$ using the following fact:

Fact 4 (Schur orthonormality).

If $\chi$ is the character of an arbitrary representation $\phi$ , and $\chi_{i}$ is the character of an irrep $\phi_{i}$ , the multiplicity $a_{i}$ of $\phi_{i}$ is

a_{i}=\frac{1}{|G|}\sum_{U\in G}\chi_{i}^{*}(U)\chi(U).

For a proof, see [37].

Since the trivial irrep has $\chi_{i}(U)=1$ , we have that the multiplicity of the trivial irrep in $(\eta_{\text{ideal}}\otimes\phi_{i})$ is given by

\mathbbm{E}\left[\text{Tr}\left(\hat{U}\otimes\phi_{i}^{*}(U)\right)\right]=\mathbbm{E}\left[\chi_{i}^{*}(U)\text{Tr}\left(\hat{U}\right)\right]=a_{i}.

In other words, $\tilde{\eta}_{\text{ideal}}(i)$ is a rank- $a_{i}$ projector.

We can explicitly find the form of $\tilde{\eta}_{\text{ideal}}(i)$ by constructing $a_{i}$ trivial irreps of $(\eta_{\text{ideal}}\otimes\phi_{i})$ . Let $\{|\psi^{i}_{n}\rangle\rangle\}$ be an orthonormal basis for $\mathcal{H}_{i}$ , and let $\{|\psi^{i,j}_{n}\rangle\rangle\}$ be the corresponding basis for the $j$ th copy of $\mathcal{H}_{i}$ inside $\mathcal{H}\otimes\mathcal{H}$ . It is straightforward to show that

|\Psi^{i,j}\rangle\rangle:=\frac{1}{\sqrt{d_{i}}}\sum_{n=1}^{d_{i}}|\psi_{n}^{i,j}\rangle\rangle\otimes|\psi_{n}^{i}\rangle\rangle

spans an irrep for each $j=1,...,a_{i}$ . Therefore,

\tilde{\eta}_{\text{ideal}}(i)=\sum_{j=1}^{a_{i}}|\Psi^{i,j}\rangle\rangle\langle\langle\Psi^{i,j}|

(30)

A realistic experiment will have gates described by a function $\eta(U)$ that is some small perturbation from $\eta_{\text{ideal}}(U)$ . Perturbing $\eta_{\text{ideal}}(U)$ by a small amount will perturb $\tilde{\eta}_{\text{ideal}}(i)$ by a small amount, since the Fourier transform is a linear operation. Thus $\tilde{\eta}(i)$ is a perturbation of a rank- $a_{i}$ projector for high-fidelity gates, so that $\tilde{\eta}(i)$ has $a_{i}$ eigenvalues close to $1$ , which we will denote by $\lambda_{i,j}$ , and the remaining eigenvalues close to $0$ . This is sufficient to make $S_{i}(N)$ dominanted by $a_{i}$ exponential decays, corresponding to the $a_{i}$ largest eigenvalues (see Eq. 29). This proves Thm. 2.

A.3 Computing the average fidelity

If we define $\eta(U)=\hat{\Lambda}_{U}\hat{U}$ , with $\Lambda_{U}$ the gate-dependent error channel, then we can define an average fidelity

F_{\text{av}}=\frac{\mathbbm{E}\left[\text{Tr}(\hat{\Lambda}_{U})\right]+d}{d^{2}+d}

(31)

Comparing to Eq. 10, we see that this is simply the average of the individual fidelities $F_{\Lambda_{U}}$ .

We can express $F_{\text{av}}$ in terms of the $a_{i}$ largest eigenvalues of $\tilde{\eta}(i)$ as follows. We first note that we may write

	$\displaystyle\mathbbm{E}\left[\text{Tr}\left(\hat{\Lambda}_{U}\right)\right]$	$\displaystyle=\mathbbm{E}\left[\text{Tr}\left(\eta(U)\eta^{\dagger}_{\text{ideal}}(U)\right)\right]$
		$\displaystyle=\sum_{i}d_{i}\text{Tr}\left(\tilde{\eta}(i)\tilde{\eta}^{\dagger}_{\text{ideal}}(i)\right)$
		$\displaystyle=\sum_{i=1}^{I}\sum_{j=1}^{a_{i}}d_{i}\langle\langle\Psi^{i,j}\|\tilde{\eta}(i)\|\Psi^{i,j}\rangle\rangle$

where in the second line we used the Parseval identity (Eq. 27) to move to Fourier space, and in the third line we used the explicit form of $\tilde{\eta}_{\text{ideal}}(i)$ (Eq. 30). To first order in $\left(\tilde{\eta}(i)-\tilde{\eta}_{\text{ideal}}(i)\right)$ , we have that

\sum_{j=1}^{a_{i}}\langle\langle\Psi^{i,j}|\tilde{\eta}(i)|\Psi^{i,j}\rangle\rangle\approx\sum_{j=1}^{a_{i}}\lambda_{i,j}

Therefore, we can rewrite Eq. 31 as

F_{\text{av}}\approx\frac{\sum_{i=1}^{I}d_{i}\sum_{j=1}^{a_{i}}\lambda_{i,j}+d}{d^{2}+d}

which is the same form as Eq. 5 in the case of gate-independent noise.

Appendix B The generalized Clifford group is a unitary 2-design

In this Appendix, we prove the generalized Clifford group considered in Section IV.2 is a unitary 2-design. We will give a fully general treatment for arbitrary sets of $n$ qudits with $d>2$ prime, although we need only the case of $n=1$ , $d=3$ for our subspace benchmarking above. This result can be inferred from results proven in [54], but we give a direct proof below. We first review the construction of the generalized Clifford groups as introduced in [49].

For a $d$ -level system, define analogues of the $X$ and $Z$ qubit operators [50]:

X|z\rangle=|z+1\rangle\qquad Z|z\rangle=\omega^{z}|z\rangle

where $\omega:=e^{2\pi i/d}$ and addition is performed modulo $d$ . These generalized $X$ and $Z$ operators are unitary and satisfy $ZX=\omega XZ$ .

For a set of $n$ qudits, define the $d$ -dimensional generalization of the Pauli group as (this only holds for $d$ odd; see [49] for the definition for $d$ even):

\mathcal{P}:=\{\omega^{\eta}X_{1}^{a_{1}}Z_{1}^{b_{1}}\cdots X_{n}^{a_{n}}Z_{n}^{b_{n}}:\eta,a_{i},b_{i}\in\mathbbm{Z}_{d}\}.

We will write a general element of the Pauli group as

\omega^{\eta}X_{1}^{a_{1}}Z_{1}^{b_{1}}\cdots X_{n}^{a_{n}}Z_{n}^{b_{n}}:=\omega^{\eta}X\!\!Z(\vec{v}),\quad\vec{v}:=\left(\begin{smallmatrix}\vec{a}\\ \vec{b}\end{smallmatrix}\right).

Multiplication of general elements of the Pauli group is given by

X\!\!Z(\vec{v})X\!\!Z(\vec{w})=\omega^{\vec{v}^{T}Q\vec{w}}X\!\!Z(\vec{v}+\vec{w})

where $Q$ is defined by $Q=\left(\begin{smallmatrix}0&0\\ \mathbbm{1}&0\end{smallmatrix}\right)$ . This demonstrates that $\mathcal{P}$ is indeed a group.

The generalized Clifford group is defined to be the set of all unitaries that stabilize $\mathcal{P}$ :

G=\{U:U\mathcal{P}U^{\dagger}=\mathcal{P}\}.

An element $U\in G$ is defined (up to a global phase) by its action on $X_{i}$ and $Z_{i}$ . We define the matrix $M$ and vector $\vec{h}$ such that for each unit vector $\hat{e}_{i}\in\mathbbm{Z}_{d}^{2d}$ we have

UX\!\!Z(\hat{e}_{i})U^{\dagger}=\omega^{h_{i}}X\!\!Z(M\hat{e}_{i})

It then follows that a general element $X\!\!Z(a)$ is transformed as

\begin{array}[]{c}UX\!\!Z(\vec{v})U^{\dagger}=\omega^{\eta}X\!\!Z(M\vec{v})\\ \eta:=\left(\vec{h}-\frac{\text{diag}(M^{T}QM)}{2}\right)^{T}\vec{v}+\vec{v}^{T}\left(M^{T}QM-Q\right)\frac{\vec{v}}{2}\end{array}

(32)

Not every matrix $M$ can be realized by a unitary operator. To derive a restriction on $M$ , we consider the commutation relation (where we define $P=Q-Q^{T}$ ):

	$\displaystyle X\!\!Z(\vec{v})X\!\!Z(\vec{w})$	$\displaystyle=\omega^{\vec{v}^{T}P\vec{w}}X\!\!Z(\vec{w})X\!\!Z(\vec{v})$
	$\displaystyle UX\!\!Z(\vec{v})X\!\!Z(\vec{w})U^{\dagger}$	$\displaystyle=\omega^{\vec{v}^{T}P\vec{w}}UX\!\!Z(\vec{w})X\!\!Z(\vec{v})U^{\dagger}$
	$\displaystyle X\!\!Z(M\vec{v})X\!\!Z(M\vec{w})$	$\displaystyle=\omega^{\vec{v}^{T}P\vec{w}}X\!\!Z(M\vec{w})X\!\!Z(M\vec{v})$
	$\displaystyle\omega^{\vec{v}^{T}M^{T}PM\vec{w}}X\!\!Z(M\vec{w})X\!\!Z(M\vec{v})$	$\displaystyle=\omega^{\vec{v}^{T}P\vec{w}}X\!\!Z(M\vec{w})X\!\!Z(M\vec{v})$

where we have ignored phase factors common to both sides. We see that we must have $P=M^{T}PM$ ; such an $M$ is called a symplectic matrix. This is the only restriction on $M,h$ , as [49] demonstrated how to explicitly construct unitaries to implement any $M,h$ provided $M$ is symplectic.

To prove $G$ forms a unitary $2$ -design, we need to show (see Section IV.2 of the main text)

\frac{1}{|G|}\sum_{U\in G}p(U,U^{*})=\int dU\ p(U,U^{*})

for any balanced polynomial $p(U,U^{*})$ of degree at most $2$ in the elements of $U$ and $U^{*}$ . Any such $p(U,U^{*})$ can be written as a linear combination of terms of the form $UAU^{\dagger}BUCU^{\dagger}$ and $UDU^{\dagger}$ , where $A,B,C,D$ are matrices. We are thus reduced to proving

\frac{1}{|G|}\sum_{U\in G}UAU^{\dagger}BUCU^{\dagger}=\int dU\ UAU^{\dagger}BUCU^{\dagger}

(33)

\frac{1}{|G|}\sum_{U\in G}UDU^{\dagger}=\int dU\ UDU^{\dagger}

(34)

for arbitrary matrices $A,B,C,D$ .

In the following, we will make repeated use of an elementary identity of complex roots of unity.

Fact 5.

If $\vec{w}\in\mathbbm{Z}_{d}^{2n}\setminus\{0\}$ is any nonzero vector, then

\sum_{\vec{v}}\omega^{\vec{v}^{T}\vec{w}}=0.

B.1 Degree 1 polynomials

Let’s start by proving Eq. 34. Without loss of generality, we can assume $D=X\!\!Z(\vec{v})$ , since such matrices form a basis. The RHS of this equation is invariant under conjugation by arbitrary unitaries; thus, it must be proportional to the identity matrix. Noting that $\operatorname{Tr}(\text{RHS})=\operatorname{Tr}(D)$ and that $\operatorname{Tr}\left[X\!\!Z(\vec{v})\right]=0$ whenever $\vec{v}\neq 0$ , we find

RHS

\displaystyle=\left\{\begin{array}[]{ll}\mathbbm{1},&\vec{v}=0\\ 0,&\text{else}.\end{array}\right.

We evaluate the LHS by using Eq. 32 for the conjugation of a general Pauli element:

	LHS	$\displaystyle=\frac{1}{\|G\|}\sum_{U\in G}UX\!\!Z(\vec{v})U^{\dagger}$
		$\displaystyle=\frac{1}{\|G\|}\sum_{\begin{subarray}{c}M,\vec{h}\\ M^{T}PM=P\end{subarray}}\omega^{\eta}X\!\!Z(M\vec{v})$

We note that $\eta=\vec{h}^{T}\vec{v}+(\cdots)$ , where $(\cdots)$ denotes terms that do not depend on $\vec{h}$ . We see by Fact 5 that for fixed $M$ the sum over $\vec{h}$ gives zero unless $\vec{v}=0$ , while when $\vec{v}=0$ it is clear $\text{LHS}=\mathbbm{1}$ . This proves Eq. 34.

B.2 Degree 2 polynomials

We now turn to Eq. 33. We prove this using methods from [9], who proved the case $d=2$ . First, we note that the RHS of Eq. 33 is covariant in $B$ : sending $B\rightarrow UBU^{\dagger}$ sends $\text{RHS}\rightarrow U(\text{RHS})U^{\dagger}$ for any unitary $U$ . The only covariant linear functions of $B$ are $\frac{\operatorname{Tr}(B)\mathbbm{1}}{d^{n}}$ and $\left[B-\frac{\operatorname{Tr}(B)\mathbbm{1}}{d^{n}}\right]$ , so the RHS must be of the form [5]

\text{RHS}=q\left[B-\frac{\operatorname{Tr}(B)\mathbbm{1}}{d^{n}}\right]+p\frac{\operatorname{Tr}(B)\mathbbm{1}}{d^{n}}.

(35)

To determine $p$ we plug in $B=\mathbbm{1}$ and note that

\text{RHS}=\int dU\ UACU^{\dagger}=\frac{\operatorname{Tr}(AC)}{d^{n}}\mathbbm{1},

while simultaneously according to Eq. 35,

\text{RHS}=p\mathbbm{1}

so $p=\frac{\operatorname{Tr}(AC)}{d^{n}}$ . To determine $q$ , we consider plugging in $B=|i\rangle\langle j|$ . Denoting the result when plugging in $B=|i\rangle\langle j|$ as $(\text{RHS})_{ij}$ , we can evaluate

	$\displaystyle\sum_{i,j}\langle i\|(\text{RHS})_{ij}\|j\rangle$	$\displaystyle=\sum_{i,j}\int dU\ \langle i\|UAU^{\dagger}\|i\rangle\langle j\|UCU^{\dagger}\|j\rangle$
		$\displaystyle=\operatorname{Tr}(A)\operatorname{Tr}(C).$

On the other hand, Eq. 35 gives

\sum_{i,j}\langle i|(\text{RHS})_{ij}|j\rangle=(d^{2n}-1)q+p

so $q=\frac{d^{n}\operatorname{Tr}(A)\operatorname{Tr}(C)-\operatorname{Tr}(AC)}{d^{n}(d^{2n}-1)}$ . Thus in total, we have

\text{RHS}=\frac{d^{n}\operatorname{Tr}(A)\operatorname{Tr}(C)-\operatorname{Tr}(AC)}{d^{n}(d^{2n}-1)}\left[B-\frac{\operatorname{Tr}(B)\mathbbm{1}}{d^{n}}\right]\\ +\frac{\operatorname{Tr}(AC)\operatorname{Tr}(B)\mathbbm{1}}{d^{2n}}.

(36)

Without loss of generality, we can specialize to the case where $A=X\!\!Z(\vec{v}_{A})$ , $B=X\!\!Z(\vec{v}_{B})$ , and $C=X\!\!Z(\vec{v}_{C})$ , whence Eq. 36 gives

\text{RHS}=\left\{\begin{array}[]{ll}X\!\!Z(\vec{v}_{B}),&\vec{v}_{A}=\vec{v}_{C}=0\\ \omega^{-\vec{v}_{A}^{T}Q\vec{v}_{A}}\mathbbm{1},&\vec{v}_{A}=-\vec{v}_{C}\neq 0,\ \vec{v}_{B}=0\\ -\frac{\omega^{-\vec{v}_{A}^{T}Q\vec{v}_{A}}}{d^{2n}-1}X\!\!Z(\vec{v}_{B}),&\vec{v}_{A}=-\vec{v}_{C}\neq 0,\ \vec{v}_{B}\neq 0\\ 0,&\text{else}.\end{array}\right.

We now need to evaluate the LHS of Eq. 33 for each of the four cases above. In the first case, we find

\displaystyle\text{LHS}=\frac{1}{|G|}\sum_{U\in G}X\!\!Z(\vec{v}_{B})=X\!\!Z(\vec{v}_{B})

In the second case, we use Eq. 32 to simplify each summand in the LHS

	$\displaystyle U$	$\displaystyle X\!\!Z(\vec{v}_{A})U^{\dagger}UX\!\!Z(\vec{v}_{C})U^{\dagger}$
		$\displaystyle=\omega^{\eta_{A}+\eta_{C}}X\!\!Z(M\vec{v}_{A})X\!\!Z(M\vec{v}_{C})$
		$\displaystyle=\omega^{\eta_{A}+\eta_{C}+\vec{v}_{A}^{T}M^{T}QM\vec{v}_{C}}\mathbbm{1}$
		$\displaystyle=\omega^{\vec{v}_{A}^{T}\left(M^{T}QM-Q\right)\vec{v}_{A}-\vec{v}_{A}^{T}M^{T}QM\vec{v}_{A}}\mathbbm{1}$
		$\displaystyle=\omega^{-\vec{v}_{A}^{T}Q\vec{v}_{A}}\mathbbm{1}.$

Therefore, the average over the group $G$ gives $\omega^{-\vec{v}_{A}^{T}Q\vec{v}_{A}}\mathbbm{1}$ .

In the third case, we again simplify each summand using Eq. 32, but with an additional $B$ in between:

	$\displaystyle UX\!\!Z(\vec{v}_{A})U^{\dagger}X\!\!Z(\vec{v}_{B})UX\!\!Z(\vec{v}_{C})U^{\dagger}$
	$\displaystyle\ =\omega^{\eta_{A}+\eta_{C}}X\!\!Z(M\vec{v}_{A})X\!\!Z(\vec{v}_{B})X\!\!Z(M\vec{v}_{C})$
	$\displaystyle\ =\omega^{\eta_{A}+\eta_{C}+\vec{v}_{A}^{T}M^{T}Q\vec{v}_{B}-\vec{v}_{B}^{T}QM\vec{v}_{A}-\vec{v}_{A}^{T}M^{T}QM\vec{v}_{A}}X\!\!Z(\vec{v}_{B})$
	$\displaystyle\ =\omega^{\vec{v}_{A}^{T}M^{T}P\vec{v}_{B}-\vec{v}_{A}^{T}Q\vec{v}_{A}}X\!\!Z(\vec{v}_{B}).$

The average over $\vec{h}$ does not affect this sum, so we only need to consider the average over $M$ . We evaluate the average by realizing that if $d$ is prime, the Clifford group sends every non-identity Pauli string to every other non-identity Pauli string uniformly. Thus, letting $M$ run over all symplectic matrices makes $M\vec{v}_{A}$ run uniformly over all vectors $M\vec{v}_{A}\in\mathbbm{Z}_{d}^{2n}\setminus\{0\}$ . Therefore, the LHS is given by

	LHS	$\displaystyle=\frac{1}{d^{2n-1}}\sum_{\vec{v}\neq 0}\omega^{\vec{v}^{T}P\vec{v}_{B}-\vec{v}_{A}^{T}Q\vec{v}_{A}}X\!\!Z(\vec{v}_{B})$
		$\displaystyle=-\frac{\omega^{-\vec{v}_{A}^{T}Q\vec{v}_{A}}}{d^{2n-1}}X\!\!Z(\vec{v}_{B})\left[1-\sum_{\vec{v}}\omega^{\vec{v}^{T}P\vec{v}_{B}}\right]$
		$\displaystyle=-\frac{\omega^{-\vec{v}_{A}^{T}Q\vec{v}_{A}}}{d^{2n-1}}X\!\!Z(\vec{v}_{B})$

where in the final step, we used Fact 5.

In the last case, we have that each summand is of the form

	$\displaystyle UX\!\!Z(\vec{v}_{A})U^{\dagger}X\!\!Z(\vec{v}_{B})UX\!\!Z(\vec{v}_{C})U^{\dagger}$
	$\displaystyle\ =\omega^{\eta_{A}+\eta_{C}+\vec{v}_{A}^{T}M^{T}Q\vec{v}_{B}-\vec{v}_{B}^{T}QM\vec{v}_{A}-\vec{v}_{A}^{T}M^{T}QM\vec{v}_{A}}X\!\!Z(\vec{v}_{B})$
	$\displaystyle\ =\omega^{\vec{h}^{T}(\vec{v}_{A}+\vec{v}_{C})+(\cdots)}X\!\!Z(\vec{v}_{B})$

where $(\cdots)$ represents terms that are independent of $\vec{h}$ . We can again apply Fact 5 to find that the sum over $\vec{h}$ gives zero. We have thus proved $\text{LHS}=\text{RHS}$ for each of the four cases, which establishes Eq. 33.

Appendix C Leakage RB irreps

Let $G$ be a unitary group indexed by $b\in B$ ,

	$\displaystyle G$	$\displaystyle=\{U_{b,\sigma}:b\in B\ \sigma=\pm 1\}$
		$\displaystyle=\{U_{1,b}\oplus\sigma U_{2,b}:b\in B,\ \sigma=\pm 1\},$

where $G_{1}=\{U_{1,b}:b\in B\}$ and $G_{1}=\{U_{2,b}:b\in B\}$ are each unitary 1-designs on their respective subspaces. First, we prove that $|\mathbbm{1}_{1}\rangle\rangle$ and $|\mathbbm{1}_{2}\rangle\rangle$ are the only trivial irreps of the natural representation of $G$ . Next, we prove that if $G_{1}$ and $G_{2}$ are in addition unitary 2-designs and $d_{1}\neq d_{2}$ then $\mathcal{H}_{1\perp}$ is irreducible and multiplicity-free.

We start with the trivial irreps. It is clear that both $|\mathbbm{1}_{1}\rangle\rangle$ and $|\mathbbm{1}_{2}\rangle\rangle$ are trivial irreps. The trivial irrep has $\chi_{0}(U)=1$ , so Fact 4 gives

	$\displaystyle a_{0}$	$\displaystyle=\frac{1}{\|G\|}\sum_{U\in G}\chi(U)$
		$\displaystyle=\frac{1}{2\|B\|}\sum_{\begin{subarray}{c}b\in B\\ \sigma=\pm\end{subarray}}\operatorname{Tr}(U_{b,\sigma}\otimes U_{b,\sigma}^{*})$
		$\displaystyle=\frac{1}{2\|B\|}\sum_{\begin{subarray}{c}b\in B\\ \sigma=\pm\end{subarray}}\left[\begin{array}[]{l}\operatorname{Tr}(U_{1,b}\otimes U_{1,b}^{})+\sigma\operatorname{Tr}(U_{1,b}\otimes U_{2,b}^{})\\ +\sigma\operatorname{Tr}(U_{2,b}\otimes U_{1,b}^{})+\operatorname{Tr}(U_{2,b}\otimes U_{2,b}^{})\end{array}\right]$
		$\displaystyle=\frac{1}{\|B\|}\sum_{b\in B}\left[\operatorname{Tr}(U_{1,b}\otimes U_{1,b}^{})+\operatorname{Tr}(U_{2,b}\otimes U_{2,b}^{})\right]$
		$\displaystyle=\int dU_{1}\operatorname{Tr}(U_{1}\otimes U_{1}^{})+\int dU_{2}\operatorname{Tr}(U_{2}\otimes U_{2}^{})$

where in the last line we used the fact that $G_{1}$ and $G_{2}$ are unitary 1-designs. These integrals just give the number of trivial irreps of the full unitary group on $\mathcal{H}_{1}$ and $\mathcal{H}_{2}$ , respectively, which are known to be $1$ . Thus, there are only two trivial irreps of the full unitary group.

Now, we consider $\mathcal{H}_{1\perp}$ . First, we show $\mathcal{H}_{1\perp}$ is irreducible by using Fact 3. Noting $\chi_{1,\perp}(U_{b,\pm})=\left(|\operatorname{Tr}(U_{1,b})|^{2}-1\right)$ , we have

	$\displaystyle\frac{1}{\|G\|}\sum_{U\in G}\|\chi_{1\perp}(U)\|^{2}$	$\displaystyle=\frac{1}{2\|B\|}\sum_{\begin{subarray}{c}b\in B\\ \sigma=\pm\end{subarray}}\left(\|\operatorname{Tr}(U_{1,b})\|^{2}-1\right)^{2}$
		$\displaystyle=\frac{1}{\|B\|}\sum_{b\in B}\left(\|\operatorname{Tr}(U_{1,b})\|^{2}-1\right)^{2}$
		$\displaystyle=\int dU_{1}\ \left(\|\operatorname{Tr}(U_{1})\|^{2}-1\right)^{2}$
		$\displaystyle=1$

where the third equality follows from the unitary 2-design property, and the fourth follows from the fact that $\mathcal{H}_{1\perp}$ is an irrep of the natural representation of the full unitary group on $\mathcal{H}_{1}\otimes\mathcal{H}_{1}$ . Thus, we have $\mathcal{H}_{1\perp}$ irreducible.

To finish, we must prove that no other irrep of the natural representation is isomorphic to $\mathcal{H}_{1\perp}$ . Every irrep of the natural representation is a subrepresentation of $\mathcal{H}_{1}\otimes\mathcal{H}_{1}$ , $\mathcal{H}_{1}\otimes\mathcal{H}_{2}$ , $\mathcal{H}_{2}\otimes\mathcal{H}_{1}$ , or $\mathcal{H}_{2}\otimes\mathcal{H}_{2}$ , since these subspaces are all invariant under the action of $G$ . We know that the decomposition of $\mathcal{H}_{1}\otimes\mathcal{H}_{1}$ into irreps is $\mathcal{H}_{1}\otimes\mathcal{H}_{1}\simeq\mathcal{H}_{10}\otimes\mathcal{H}_{1\perp}$ , by our work above, and thus no irreps in $\mathcal{H}_{1}\otimes\mathcal{H}_{1}$ can be isomorphic to $\mathcal{H}_{1\perp}$ besides $\mathcal{H}_{1\perp}$ itself. Similarly, we know that the decomposition of $\mathcal{H}_{2}\otimes\mathcal{H}_{2}$ into irreps is $\mathcal{H}_{2}\otimes\mathcal{H}_{2}\simeq\mathcal{H}_{20}\otimes\mathcal{H}_{2\perp}$ . We can ensure $\mathcal{H}_{1\perp}\not\simeq\mathcal{H}_{2\perp}$ by requiring $d_{1}\neq d_{2}$ , as in the main text. We then have that no isomorphic representation exists in $\mathcal{H}_{2}\otimes\mathcal{H}_{2}$ . For $\mathcal{H}_{1}\otimes\mathcal{H}_{2}$ , and similarly for $\mathcal{H}_{2}\otimes\mathcal{H}_{1}$ , we note that the character of the subrepresentation $\mathcal{H}_{1}\otimes\mathcal{H}_{2}$ is given by $\chi_{12}(U_{b,\sigma})=\sigma\operatorname{Tr}(U_{1,b})\operatorname{Tr}(U_{2,b})^{*}$ , and use Fact 4:

	$\displaystyle\frac{1}{\|G\|}$	$\displaystyle\sum_{U\in G}\chi_{1\perp}^{*}(U)\chi_{12}(U)$
		$\displaystyle=\frac{1}{2\|B\|}\sum_{\begin{subarray}{c}b\in B\\ \sigma=\pm\end{subarray}}\sigma(\|\operatorname{Tr}(U_{1,b_{1}})\|^{2}-1)\operatorname{Tr}(U_{1,b})\operatorname{Tr}(U_{2,b})^{*}$
		$\displaystyle=0$

which shows that $\mathcal{H}_{1\perp}$ is an irrep with multiplicity $1$ .

Note that we could also consider a group

G^{\prime}=\{U_{b,\phi}:b\in B\}=\{U_{1,b}\oplus(e^{i\phi}U_{2,b}):b\in B\}

with an arbitrary phase between subspaces $1$ and $2$ rather than simply a $\pm 1$ phase; the proof is identical. Many experimental platforms can easily implement a random phase between two subspaces, especially if the leakage subspace is at a different energy than the computational subspace, making this group potentially easier to sample from. We can also still compute $F_{\Lambda,1}$ with $\{U_{2,a}\}$ only a unitary $1$ -design, provided $\mathcal{H}_{2}\otimes\mathcal{H}_{2}$ does not contain an irrep isomorphic to $\mathcal{H}_{1\perp}$ . Finally, in the case that $d_{1}=d_{2}$ , we can instead simply require that there exists some $b\in B$ such that $|\operatorname{Tr}(U_{1,b})|^{2}\neq|\operatorname{Tr}(U_{2,b})|^{2}$ , a much weaker condition that still suffices to ensure $\mathcal{H}_{1\perp}\not\simeq\mathcal{H}_{2\perp}$ .

	$\displaystyle\hat{\Lambda}_{G}\hat{U}$	$\displaystyle=\frac{1}{\|G\|}\sum_{U^{\prime}\in G}\hat{U}^{\prime\dagger}\hat{\Lambda}\hat{U}^{\prime}\hat{U}$
		$\displaystyle=\frac{1}{\|G\|}\sum_{U^{\prime}\in G}\hat{U}\hat{U}^{\dagger}\hat{U}^{\prime\dagger}\hat{\Lambda}\hat{U^{\prime}}\hat{U}$
		$\displaystyle=\hat{U}\frac{1}{\|G\|}\sum_{(U^{\prime}U)\in G}(\hat{U}^{\prime}\hat{U})^{\dagger}\hat{\Lambda}(\hat{U}^{\prime}\hat{U})=\hat{U}\hat{\Lambda}_{G}$

	$\displaystyle S_{i}(N)$	$\displaystyle=\frac{1}{\text{dim}(\mathcal{\overline{H}}_{\overline{i}})}\langle\langle M_{i}\|\hat{\Lambda}_{M}\hat{\Lambda}\left(\hat{\Lambda}_{G}\right)^{N}\hat{P}_{\overline{i}}\hat{\Lambda}_{P}\|\rho_{i}\rangle\rangle$
		$\displaystyle=\frac{1}{\text{dim}(\mathcal{\overline{H}}_{\overline{i}})}\langle\langle M_{i}\|\hat{\Lambda}_{M}\hat{\Lambda}\left(\bigoplus_{i^{\prime}}\hat{Q}_{i^{\prime}}\otimes\mathbbm{1}_{i^{\prime}}\right)^{N}\hat{P}_{\overline{i}}\hat{\Lambda}_{P}\|\rho_{i}\rangle\rangle$
		$\displaystyle=\frac{1}{\text{dim}(\mathcal{\overline{H}}_{\overline{i}})}\langle\langle M_{i}\|\hat{\Lambda}_{M}\hat{\Lambda}\left(\hat{Q}_{i}^{N}\otimes\mathbbm{1}_{i}\right)\hat{P}_{\overline{i}}\hat{\Lambda}_{P}\|\rho_{i}\rangle\rangle$

	$\displaystyle\frac{1}{\|G\|}\sum_{U\in G}\|\chi_{T\perp}(U)\|^{2}$	$\displaystyle=\frac{1}{3\|G_{T}\|}\sum_{\begin{subarray}{c}U_{T}\in G_{T}\\ \eta=0,1,2\end{subarray}}\|\chi_{T\perp}(U_{T})\|^{2}$
		$\displaystyle=\frac{1}{\|G_{T}\|}\sum_{G_{T}}\left(\|\operatorname{Tr}(U_{T})\|^{2}-1\right)^{2}$
		$\displaystyle=\int dU_{\alpha}\ \left(\|\operatorname{Tr}(U_{\alpha})\|^{2}-1\right)^{2}$
		$\displaystyle=1$

	$\displaystyle\frac{1}{\|G\|}\sum_{U\in G}\|\chi_{ST}(U)\|^{2}$	$\displaystyle=\frac{1}{3\|G_{T}\|}\sum_{\begin{subarray}{c}U_{T}\in G_{T}\\ \eta=0,1,2\end{subarray}}\|\operatorname{Tr}(U_{T})\|^{2}$
		$\displaystyle=\int dU_{T}\ \|\operatorname{Tr}(U_{T})\|^{2}$
		$\displaystyle=1$

	$\displaystyle\|e_{0}\rangle\rangle$	$\displaystyle=\frac{S}{\sqrt{d_{1}}(L+S)}\|\hat{\mathbbm{1}}_{1}\rangle\rangle+\frac{L}{\sqrt{d_{2}}(L+S)}\|\hat{\mathbbm{1}}_{2}\rangle\rangle$
	$\displaystyle\|\overline{e}_{0}\rangle\rangle$	$\displaystyle=\sqrt{d_{1}}\|\hat{\mathbbm{1}}_{1}\rangle\rangle+\sqrt{d_{2}}\|\hat{\mathbbm{1}}_{2}\rangle\rangle$
	$\displaystyle\|e_{1}\rangle\rangle$	$\displaystyle=\sqrt{d_{2}}\|\hat{\mathbbm{1}}_{1}\rangle\rangle-\sqrt{d_{1}}\|\hat{\mathbbm{1}}_{2}\rangle\rangle$
	$\displaystyle\|\overline{e}_{1}\rangle\rangle$	$\displaystyle=\frac{L}{\sqrt{d_{2}}(L+S)}\|\hat{\mathbbm{1}}_{1}\rangle\rangle-\frac{S}{\sqrt{d_{1}}(L+S)}\|\hat{\mathbbm{1}}_{2}\rangle\rangle$
	$\displaystyle\lambda$	$\displaystyle=1-L-S$

	$\displaystyle L$	$\displaystyle=(1-B)(1-\lambda)$		(19)
	$\displaystyle S$	$\displaystyle=B(1-\lambda)$		(20)

Character randomized benchmarking for non-multiplicity-free groups with applications to subspace, leakage, and matchgate randomized benchmarking

Abstract

I Introduction

II Mathematical Preliminaries

III The generalized character randomized benchmarking procedure

III.1 Deriving the decays

Fact 1 (Schur’s Lemma).

Fact 2 (Projection formula).

Theorem 1 (Form of GG-twirls).

Proof.

III.2 Computing the fidelity

III.3 Scaling and Feasibility

IV Application: Subspace randomized benchmarking

IV.1 Constructing the benchmarking group

IV.2 Irreps of the benchmarking group

Fact 3 (Schur normalization).

IV.3 Benchmarking GG

V Application: Leakage randomized benchmarking

V.1 Deriving LL and SS

V.2 Deriving FΛ,1F_{\Lambda,1}

V.3 Example: Two-qubit logical encodings

VI Application: Matchgate RB

VI.1 The matchgate group

Claim 1.

Claim 2.

VI.1.1 Proof of claims

Proof of Claim 1.

Proof of Claim 2.

VI.2 Irreps of the matchgate group

Claim 3.

Proof.

VI.3 Benchmarking the matchgate group

VII Conclusion and Discussions

Acknowledgements

References

Appendix A Gate-dependent errors

Theorem 2.

A.1 The generalized Fourier transform and its application to character RB

A.2 Simplifying the decay

Fact 4 (Schur orthonormality).

A.3 Computing the average fidelity

Appendix B The generalized Clifford group is a unitary 2-design

Fact 5.

B.1 Degree 1 polynomials

B.2 Degree 2 polynomials

Appendix C Leakage RB irreps

Theorem 1 (Form of $G$ -twirls).

IV.3 Benchmarking $G$

V.1 Deriving $L$ and $S$

V.2 Deriving $F_{\Lambda,1}$