Compilation-informed probabilistic logical-error cancellation

TL;DR

This paper presents a quantum error mitigation scheme that reduces resource requirements for high-precision quantum computations by leveraging compilation information, circuit-agnostic methods, and linear programming, aiding near-term quantum hardware.

Contribution

It introduces a compilation-informed, circuit-agnostic quantum error mitigation method that achieves unbiased estimates with constant overhead, independent of circuit size and code distance.

Findings

Significantly reduces quantum resource needs for high-precision tasks.

Demonstrates effectiveness through numerical simulations and knot polynomial estimation.

Provides a practical approach for fault-tolerant quantum computation in the near term.

Abstract

The potential of quantum computers to outperform classical ones in practically useful tasks remains challenging in the near term due to scaling limitations and high error rates of current quantum hardware. While quantum error correction (QEC) offers a clear path towards fault tolerance, overcoming the scalability issues will take time. Early applications will likely rely on QEC combined with quantum error mitigation (QEM). We introduce a QEM scheme against both compilation errors and logical-gate noise that is circuit-, QEC code-, and compiler-agnostic. The scheme builds on quasi-probability methods and uses information about the circuit's gates' compilations to attain an unbiased estimation of noiseless expectation values incurring a constant sample-complexity overhead. Moreover, it features maximal circuit size and code distance both independent of the target precision, in contrast to strategies based on QEC alone. We formulate the mitigation procedure as a linear program, demonstrate its efficacy through numerical simulations, and illustrate it for estimating the Jones polynomials of knots. Our method significantly reduces quantum resource requirements for high-precision estimations, offering a practical route towards fault-tolerant quantum computation with precision-independent overheads for fixed circuit complexity and code distance.

Tables

Table 1. Table 1 : Three different bases for CIPEC. The noisy unitary basis elements V ~ ( L ) \widetilde{V}^{(L)} use only Clifford gates { I , H , S , S † , X , Y , Z , CNOT } \{I,H,S,S^{\dagger},X,Y,Z,\text{CNOT}\} . The span of basis B ~ 1 \widetilde{B}_{1} contains the space of completely positive and trace-preserving (CPTP) maps for 2 2 -qubits and is overcomplete. It consists of the full 2 2 -qubit Clifford group plus state preparation channels Γ ~ = { I , 𝒫 | + ⟩ , 𝒫 | + y ⟩ , 𝒫 | 0 ⟩ } ⊗ 2 \widetilde{\Gamma}=\{I,\mathcal{P}_{\ket{+}},\mathcal{P}_{\ket{+y}},\mathcal{P}_{\ket{0}}\}^{\otimes 2} that prepare the corresponding state | ψ ⟩ \ket{\psi} from any ρ \rho . The basis B ~ 2 \widetilde{B}_{2} uses the same state preparation channels Γ ~ \widetilde{\Gamma} , but with a minimal number of linearly independent compositions of Clifford gates necessary so its span contains the space of (CPTP) maps, with dim ​ ( CPTP ) = 241 \text{dim}(\text{CPTP})=241 [ 32 ] . The basis B ~ 3 \widetilde{B}_{3} was proposed in [ 22 ] (see also [ 30 ] ) and is also minimal, using Clifford gates conjugated with trace non-increasing channels π z \pi_{z} that project ρ \rho into the | 0 ⟩ \ket{0} state. Its span contains the space of CP maps, with dim ​ ( CP ) = 256 \text{dim}(\text{CP})=256 . The reported values of c ∗ c_{*} assume a local depolarizing noise model after each operation with strengths 10 − 6 10^{-6} for products of single-qubit gates and 10 − 5 10^{-5} for the others, leading to ϵ 𝚀 = 10 − 5 \epsilon_{\mathtt{Q}}=10^{-5} .

$\tilde{B}$	${\tilde{V}}^{\left(L\right)}$	$\tilde{\Gamma }$	$\left|\tilde{B}\right|$	$D$	$c_{\ast }$
${\tilde{B}}_1$	Clifford group	State prep.	$11535$	$17$	$4.47$
${\tilde{B}}_2$	Cliffords (min)	State prep.	$241$	$4$	$156.2$
${\tilde{B}}_3$	Cliffords (min)	State proj.	256	$10$	$88.0$

Equations

$${ϵ}_{\mathtt{Q}}:-\underset{j\in [|\mathcal{A}|]}{\max }\Vert {\tilde{\mathcal{A}}}_j-{\mathcal{A}}_j{\Vert }_{\diamond }$$

$$\begin{array}{rl}& \underset{{\mathbf{b}}_i\in {\mathbb{R}}^{|\tilde{B}|}}{\min }{\gamma }_i:-1+\text{norm}{{\mathbf{b}}_i}_1\\ \text{s.t.}& \{{\mathcal{U}}_i={\tilde{\mathcal{V}}}_{c,i}^{(L_i)}+\sum _{j\in [|\tilde{B}|]}{b}_{i,j}{\tilde{\mathcal{B}}}_j;\end{array}$$

$${ϵ}_{\text{char}}$$

$$\le \underset{j\in [|A|]}{\max }\Vert {\mathcal{E}}_j-I{\Vert }_{\diamond }+\underset{j\in [|A|]}{\max }\Vert I-{\mathcal{E}}_j^{\prime }{\Vert }_{\diamond }$$

$$={ϵ}_{\mathtt{Q}}+\underset{j\in [|A|]}{\max }\Vert I-{\mathcal{E}}_j^{\prime }{\Vert }_{\diamond }.$$

$$\Vert \mathcal{C}{\Vert }_{\tilde{B}}\le c_{\ast }\Vert \mathcal{C}{\Vert }_{\diamond },c_{\ast }:-\underset{\text{norm}{\mathcal{C}}_{\diamond }=1}{\max }\text{norm}{\mathcal{C}}_{\tilde{B}}.$$

$$c_{\text{low}}\text{norm}{\mathcal{C}}_{\diamond }\le \text{norm}{\mathcal{C}}_{\tilde{B}}\le c_{\text{up}}\text{norm}{\mathcal{C}}_{\diamond }.$$

$$c_{\text{low}}\le \text{norm}{\mathcal{C}}_{\tilde{B}}\le c_{\text{up}}.$$

$$\begin{array}{r}{c}_{\text{low}}:-\underset{\text{norm}{\mathcal{C}}_{\diamond }=1}{\min }\text{norm}{\mathcal{C}}_{\tilde{B}}\\ c_{\text{up}}:-\underset{\text{norm}{\mathcal{C}}_{\diamond }=1}{\max }\text{norm}{\mathcal{C}}_{\tilde{B}}.\end{array}$$

$$\Vert \mathcal{U}-{\tilde{\mathcal{V}}}_c^{(L)}{\Vert }_{\diamond }$$

$$\le \Vert \mathcal{U}-{\mathcal{V}}_c^{(L)}{\Vert }_{\diamond }+\Vert {\mathcal{V}}_c^{(L)}-{\tilde{\mathcal{V}}}_c^{(L)}{\Vert }_{\diamond }$$

$$\le {ϵ}_c+L{ϵ}_{\mathtt{Q}},$$

$$S_{\text{CIPEC}}:-{\gamma }^2\Vert O{\Vert }^2\log (2/\delta )/(2{\epsilon }^2)$$

$$L{ϵ}_{\mathtt{Q}}\le \left(1-\frac{1}{\xi }-\frac{1}{\eta }\right)\epsilon .$$

$$|\text{tr}(O\mathcal{U}(\rho ))-\text{tr}(O{\tilde{\mathcal{V}}}_c^{(L)}(\rho ))|\le \sum _{i\in [G]}(L_i{ϵ}_{\mathtt{Q}}+{ϵ}_{c,i}).$$

$$|\overline{y}-\text{tr}(O\mathcal{U}(\rho ))|\le |\overline{y}-\text{tr}(O{\mathcal{V}}_c^{(L)}(\rho ))|+\sum _{i\in [G]}(L_i{ϵ}_{\mathtt{Q}}+{ϵ}_{c,i}),$$

$$S_{\text{QEC}}:-{\gamma }^2\Vert O{\Vert }^2\log (2/\delta )/(2{\epsilon }^2)$$

$$|{\overline{O}}^{\prime }-⟨O⟩|\le \epsilon +L\Vert O\Vert {ϵ}_{\text{char}}.$$

$$\text{norm}{\prod _{i\in [G]}{\tilde{\mathcal{B}}}_{{\alpha }_i}^i-\prod _{i\in [G]}{\tilde{\mathcal{B}}}_{{\alpha }_i}^{\prime i}}_{\diamond }$$

$$\le \sum _{i\in [G]}{L}_i\underset{j\in [|\mathcal{A}|]}{\max }\text{norm}{{\mathcal{E}}_j^{\prime }-{\mathcal{E}}_j}_{\diamond }$$

$$=L{ϵ}_{\text{char}},$$

$$\text{norm}{\sum _{\alpha }p(\alpha )\left(\prod _{i\in [G]}{\tilde{\mathcal{B}}}_{{\alpha }_i}^i-\prod _{i\in [G]}{\tilde{\mathcal{B}}}_{{\alpha }_i}^{\prime i}\right)}_{\diamond }\le L{ϵ}_{\text{char}}.$$

$$|{\overline{O}}^{\prime }-⟨O⟩|\le \epsilon +L\Vert O\Vert {ϵ}_{\text{char}}.$$

$$L=c_1{\log }^{c_2}(1/{ϵ}_c),$$

$$J_{\text{P}(\Sigma )}\left(e^{i\frac{2\pi }{5}}\right)=(-e^{-i\frac{3\pi }{5}}{)}^{3{\omega }_{\Sigma }}{\varphi }^{\frac{n_s}{2}-1}\text{expectationvalue}s|U_{\Sigma }|s,$$

Full text

Introduction.

Quantum computers have the potential to exponentially outperform classical computers in certain tasks [1]. These quantum advantages presuppose deep and fault-tolerant circuits on a large number of qubits: current estimates are hundreds of logical qubits at $10^{-6}$ error rates for scientific applications, and thousands at $10^{-12}$ rates for industrial applications [2]. However, real-world devices are inherently noisy, with current error rates well above those. In principle, these error rates can be reduced via quantum error-correction (QEC) techniques by using multiple physical qubits to encode each logical one [3]. Moreover, this comes also with an increase in circuit depth, as the circuit has to be compiled into the universal gate set of the code.

One of the cornerstones of quantum computing is that QEC and compilation overheads are modest: both width [4] and depth [5] are increased only by polylog factors in the number of qubits, physical noise rate and compilation error, respectively (for certain quantum LDPC codes also logarithmic in the number of logical qubits [6]). Yet, in practice, such overheads can cause years of delay for implementations. These considerations led to the development of various tools to reduce the resource requirements of quantum algorithms, including quantum error mitigation (QEM) techniques [7] and early fault-tolerant (EFT) hybrid algorithms. In QEM, the goal is to reduce the impact of noise on expectation values without actually correcting the errors by running noisy quantum circuits multiple times and classically post-processing the outcomes. In turn, EFT algorithms aim at de-quantizing (e.g., via randomization [8, 9, 10, 11]) certain components of a quantum algorithm and reserving quantum hardware strictly to crucial sub-routines. Both approaches incur an overhead in the statistical sample complexity due to the need for (possibly high-precision) estimates of random variables on quantum hardware.

Given the limitations of near-term quantum hardware, early applications will likely rely on simple QEC codes with a moderate qubit overhead followed by QEM at the level of the encoded qubits [12, 13]. This approach was explored in [14, 15], where encoded Clifford gates were protected using QEC while $T$-gate noise was mitigated using the probabilistic error-cancellation (PEC) method [16, 17]. The resulting sample overhead scales exponentially with the number of $T$ gates. Similarly, [18] used PEC with the inversion method to mitigate decoding and compilation errors, which incurs also an overhead in gate complexity. The same setting has also been explored for probabilistically synthesizing continuous logical quantum operations [19, 20, 21], which again requires a significant sample overhead.

In this Letter, we propose compilation-informed probabilistic error cancellation (CIPEC), a logical-error mitigation scheme that simultaneously removes biases from compilation errors $\epsilon_{c}$ and logical-gate errors $\epsilon_{\mathtt{Q}}$ in expectation-value estimates. The scheme requires the characterization of noisy logical gates up to an accuracy proportional to the target one, but is agnostic to circuit, basis, compiler, and QEC code. Our scheme exploits a quasi-probability decomposition similar to the compensation method of PEC [22], where unitary gates are represented as affine combinations of gate sequences drawn from the noisy native set, but at the level of logical qubits relative to a QEC code instead of physical ones. For a fixed code distance and a circuit composed of $G$ two-qubit gates, by requiring $\epsilon_{c}$ to be $\mathcal{O}(1/G)$ (independent of $\varepsilon$), CIPEC estimates expectation values to arbitrary precision $\varepsilon$ incurring a constant sample-complexity overhead. In other words, by combining standard QEC and compilation methods with quasi-probabilities, it is possible to make fault-tolerance overheads depend only on the circuit size, not the target precision. This resource reduction can be instrumental to unlock useful applications on EFT hardware, specially for high-precision estimations. Moreover, we present linear programs to find optimal quasi-probability decompositions of unitary channels in our framework efficiently and illustrate the method with numerical experiments for Jones polynomial estimations [23].

In addition, we show that, for realistic target precisions, CIPEC allows one to solve problem instances orders of magnitude larger than those achievable using only QEC (see Fig. 2). We also provide explicit estimates of scaling factors for compilation errors in diamond norm (see App. E) and worst-case negativity of random 2-qubit gates with respect to three different basis of noisy channels (see Tab. 1).

Notation.

Given a unitary operator $U$, we denote by $\mathcal{U}$ the corresponding unitary channel, acting on a state $\rho$ as $\mathcal{U}(\rho):=U\rho\,U^{\dagger}$, and by $\widetilde{\mathcal{U}}$ a noisy channel approximating $\mathcal{U}$. More generally, $\widetilde{\Lambda}$ denotes the noisy realization of any quantum channel $\Lambda$. Errors in quantum channels are quantified by the diamond norm $\|\Lambda\|_{\diamond}$, and the composition of $L$ channels $\Lambda_{L}\circ\cdots\circ\Lambda_{1}$ is denoted as $\prod_{i=1}^{L}\Lambda_{i}$. We also use $\|\boldsymbol{a}\|_{1}=\sum_{\alpha}|a_{\alpha}|$ for the $\ell_{1}$-norm of a vector $\boldsymbol{a}$, $[K]:=\{1,\ldots,K\}$ for the list of positive integers up to $K$, and $|A|$ for the cardinality of a set $A$.

Setup.

Our goal is to solve the following problem:

Problem 1** (Expectation-value estimation).**

Given a $n$-qubit state $\rho$, an observable $O$, and the description of a unitary $U$ in terms of $G$ two-qubit unitary gates $U_{i}$, estimate $\langle O\rangle\coloneq\tr(O\,U\rho\,U^{\dagger})$ up to precision $\varepsilon>0$ with probability at least $1-\delta$.

For that, we assume access to a device with limited QEC capabilities that can perform imperfect operations at the logical level, as follows.

Definition 2** (Noisy logical quantum device).**

A noisy logical quantum device $\mathtt{Q}$ relative to a chosen QEC code consists of:

$n_{p}$* physical qubits encoding $n_{\ell}\leq n_{p}$ logical qubits;* 2. 2.

a set $\widetilde{A}\coloneq\widetilde{V}\cup\widetilde{\Gamma}$ of noisy logical channels whose elements approximate ideal operations from a set $A\coloneq V\cup\Gamma$ of completely positive (CP) maps that contains a universal unitary gate set $V$ and non-unitary operations $\Gamma$. We refer to these as implementable operations and denote by

[TABLE]

their worst-case diamond-norm error.

The non-unitary operations $\Gamma$ may include, e.g., state preparation channels, state projectors, and projective measurements. As we shall see in the following, in the presence of noise these are needed to span the space of two-qubit unitaries. Without loss of generality, when considering the noisy version $\widetilde{\Gamma}$ of $\Gamma$ we assume that state preparation and measurement (SPAM) errors are negligible as they can always be incorporated by adding the relevant noisy measurement and noisy state-preparation channels to the set of implementable operations. The cost of logical gate characterization is addressed in App. A.

Implicitly, Def. 2 assumes stationary gate-dependent errors, meaning that each $\widetilde{\mathcal{A}}_{j}$ is characterized by the ideal operation $\mathcal{A}_{j}$ alone, regardless of its position in a circuit or the gates preceding it. This restriction is not strictly necessary and the construction applies to more general noise models as long as they are characterized in diamond norm. The implementable logical operations include effective logical error channels $\mathcal{E}_{j}$, i.e. $\widetilde{\mathcal{A}}_{j}=\mathcal{E}_{j}\circ\mathcal{A}_{j}$. Although Def. 2 assumes these errors are perfectly characterized, in App. D we discuss the effect of imperfect characterization. The $\mathcal{E}_{j}$ can appear as the result of physical gates and a few cycles of error correction that project the physical state back to the logical space. Moreover, we assume $\mathcal{E}_{j}$ to act locally in at most a small neighbourhood of the qubits acted on by $\mathcal{A}_{j}$. For the unitary gates in $V$, common physical error models manifest as stochastic Pauli noise in the logical level [24].

The set of universal gates $V$ is determined by the particular QEC code. Here, we focus on $V\coloneq\{I,H,S,S^{\dagger},X,Y,Z,T\}^{\otimes 2}\cup\{\text{CNOT}\}$, which is a typical gate set for stabilizer codes [25]. The QEC code also determines the relation between the error suppression and the number of physical and logical qubits. A popular example is that of distance-$d$ surface codes, for which (for a physical error rate $p$ below the code’s threshold $p_{\text{th}}$) the logical error rate is suppressed exponentially in $d$ as $\epsilon_{\mathtt{Q}}\propto\,(p/p_{\text{th}})^{(d+1)/2}$ by employing $2d^{2}-1$ physical qubits [2] per logical one. Similar relations hold for other codes (e.g., quantum LDPC codes [26]). Notice that our requirements refer to diamond norm errors, while the gate error $p$ is frequently reported as the average gate infidelity, which can be efficiently estimated via randomized benchmarking. However, the two are related as $(1+1/m)\,p\leq\epsilon_{\diamond}\leq\sqrt{(m+1)\,m\,p}$ [27], where $m$ is the dimension of the space acted on by the error channel, here assumed to be composed of two qubits only. Moreover, the lower bound is saturated for stochastic Pauli noise [27], which is often assumed for obtaining the error threshold of QEC codes [28, 26] 111Importantly, Pauli noise is known to be a good effective logical error model [24]. Moreover, the expression $\epsilon_{\mathtt{Q}}\propto\,(p/p_{\text{th}})^{(d+1)/2}$ gives the error rate per EC cycle to store a qubit. In turn, the error per logical gate can be approximated as the number of cycles the gate takes to be implemented ($2d$ cycles for the CNOT and Hadamard gates, for example) times the error per cycle [43, 18]. The proportionality factor depends on the device and code specifications, including the physical error model..

Importantly, in a scenario where the number of available physical qubits is limited, logical errors cannot be arbitrarily suppressed. In solving Problem 1, this implies the target precision $\varepsilon$ may not be achievable through QEC alone even if statistical fluctuations were absent, in which case combining error-correction with error mitigation methods can prove useful. We resort to probabilistic error cancellation (PEC), a technique proposed to mitigate errors in the context of NISQ devices [16, 30]. Previous approaches to PEC employ the so called error inversion method [18, 31], where the inverse noise channels $\mathcal{E}_{j}^{-1}$ are decomposed as linear combinations of noisy implementable channels (assuming the noise is invertible) and applied after each noisy gate in the circuit, which may cause a significant increase in circuit depth. Here, instead, we consider the compensation method for PEC [22], in which case the ideal gates $U_{i}$ in the circuit are decomposed in a similar way. For that to be possible, one must ensure that the implementable operations span the space of all $2$-qubit unitary channels (as opposed to the inversion method, which in general requires the larger space of all completely positive trace-preserving maps). Notice that, for arbitrary noise models $\mathcal{E}_{j}$, composing the noisy unitary operations $\widetilde{V}_{j}$ alone may not suffice, and we need to include non-unitary operations. This is formalized in the following definition.

Definition 3** (CIPEC-capable device).**

Let $\mathtt{Q}$ be a noisy logical quantum device as in Def. 2, ${V}^{(D)}\subset\left\{\prod_{i=1}^{D}{\mathcal{V}}_{\alpha_{i}}\,|\,\alpha_{1},\ldots,\alpha_{D}\in[|V|]\right\}$ be a set of unitary channels obtained by composing $D>0$ operations from the universal set $V$, and $\widetilde{V}^{(D)}$ be the implementable version of $V^{(D)}$. We say that $\mathtt{Q}$ is capable of running CIPEC if there exists $\widetilde{V}^{(D)}$ such that the span of $\widetilde{B}\coloneq\widetilde{V}^{(D)}\cup\,\widetilde{\Gamma}$ contains the subspace of $2$-qubit unitary channels on the logical space. We refer to $\widetilde{B}$ as the basis of implementable operations and $D$ as the maximum sequence length (or depth) of $\widetilde{B}$.

The basis $\widetilde{B}$ in Def. 3 can be overcomplete. In App. B, we discuss the effects of the basis size $|\widetilde{{B}}|$ and maximal sequence length $D$ on practical implementation. The following definition introduces the negativity, which accounts for the sample overhead of error mitigation, and its upper bound $c_{*}$, that will play an important role in the analysis of CIPEC.

Definition 4** (Worst-case negativity).**

Given a basis $\widetilde{B}$ of implementable operations as in Def. 3, the negativity of a channel $\mathcal{C}\in\text{span}(\widetilde{B})$ in this basis is $\|\mathcal{C}\|_{\widetilde{B}}:=\min_{b}\left\{\sum_{j}|b_{j}|\text{ s.t. }\mathcal{C}=\sum_{j\in|\widetilde{B}|}b_{j}\,\widetilde{\mathcal{B}}_{j}\right\}$. We denote by $c_{*}$ the worst-case negativity over all the unitary channels, $c_{*}\coloneq\max_{\norm{\mathcal{C}}_{\diamond}=1}\norm{\mathcal{C}}_{\widetilde{B}}$.

In App. B we prove that $c_{*}$ connects the negativity and the diamond norm, i.e., $\|\mathcal{C}\|_{\widetilde{B}}\leq c_{*}\|\mathcal{C}\|_{\diamond}$.

A minimal channel basis for two-qubit CP maps was proposed in [22] using Clifford gates conjugated with trace non-increasing channels $\pi_{z}(\rho)\coloneq\left(\frac{1+Z}{2}\right)\rho\left(\frac{1+Z}{2}\right)$ that project any state $\rho$ into the $\outerproduct{0}{0}$ state, but these make it inconvenient for PEC due to the need for postselection. Here, we introduce two different bases whose span contain the space of CPTP maps and investigate their influence on the statistical overhead of PEC as measured by the constant $c_{*}$ of Def. 4. For each noisy basis $\widetilde{B}$, we empirically estimate $c_{*}$ by taking the worst-case negativity over an ensemble of $10^{4}$ Haar random unitaries. In particular, we explicitly construct bases from state preparation channels $\mathcal{P}_{\ket{\psi}}(\rho)\coloneq\outerproduct{\psi}{\psi}$, which prepare the state $\ket{\psi}$ from any given $\rho$, and sequences of Clifford operations only, with no need for $T$ gates, so that all its unitary elements are classically simulable. We first construct a minimal basis of maximal sequence length $D=4$ using sequences of implementable gates and state preparations. From all possible sequences of implementable operations, we use a greedy search to choose $241$ linearly independent elements to form a minimal basis spanning the space of 2-qubit channels. We also investigate the advantage of using an overcomplete basis given by the full $2$-qubit Clifford group and state preparation channels, which reduces the negativity of the decomposition but uses longer sequences with $D=17$. The main features are summarized in Tab. 1. Since the runtime and performance guarantees of CIPEC depend on the chosen basis, an open question is whether one can go beyond the above heuristics and construct an optimal basis, i.e., one with the smallest value of $c_{*}$ while having the minimum number of elements.

Main results.

Here we introduce a new algorithm, which we refer to as Compilation-Informed Probabilistic Error Cancellation (CIPEC), to mitigate errors in the estimation of observables that remain after logical operations with only partial quantum error-correction. The idea is to statistically simulate the compiled circuit by using a hybrid classical/quantum procedure based on randomly sampling $j\in|\widetilde{B}|$ according to its importance in the quasi-probability distribution in Eq. (LABEL:eq:lp) for each gate. The algorithm is described in Alg. 5 and resembles PEC with the compensation method [18], where the compensation term is given by the implementable version of the compilation of each gate $U_{i}$ in the circuit into the universal gate set $V$ realized by the device. Using this term to augment the basis $\widetilde{B}$ realized by a CIPEC-capable device Q (see Def. 3), the negativity of the corresponding quasi-probability decomposition can be controlled. The framework is conceptually different from the standard setup of [18], which was designed for noisy-intermediate scale quantum (NISQ) devices and uses a direct noisy implementation of the desired gates $U_{i}$ on hardware (e.g., parametrized rotation gates) themselves as the compensation term. In our case, compilation and circuit noise are mitigated together with a single procedure. Moreover, while in NISQ devices error mitigation can be restricted to entangling gates, which are typically the noisiest operations [33], this is no longer valid in a fault-tolerant setting. For instance, $T$ gates require gate teleportation and magic state distillation [34], which may result in effective noise levels comparable to those of $2$-qubit gates. Therefore, our method mitigates errors from all gates in the circuit.

Algorithm 5** (CIPEC).**

*.*Input: $\{U_{i}\}_{i\in[G]}$, $O,\rho,\varepsilon,\delta$ as in Problem 1; a CIPEC-capable device $\mathtt{Q}$ with basis $\widetilde{B}\coloneq\big\{\widetilde{\mathcal{B}}_{j}\big\}_{j\in[|\widetilde{B}|]}$, worst-case negativity $c_{*}$, and able to measure in the eigenbasis of $O$; and $\omega_{1}>1$. Output: $\overline{O}$ s.t. $\big|\overline{O}-\langle O\rangle\big|\leq\varepsilon$ with prob. $1-\delta$.

For each $i\in[G]$:

1. (a)

(Classical) compile $U_{i}$ up to diamond-norm error $\epsilon_{c,i}=\log(\omega_{1})/(2c_{*}G)$ using a sequence of $L_{i}$ noiseless gates from the universal set $V$, and use the resulting unitary $\mathcal{V}^{(L_{i})}_{c,i}$ to define the augmented noisy basis $\widetilde{B}^{i}\coloneq\widetilde{B}\,\cup\big\{\widetilde{\mathcal{V}}^{(L_{i})}_{c,i}\big\}$; 2. (b)

(Classical) decompose $\mathcal{U}_{i}$ in the basis $\widetilde{B}^{i}$ by solving the following optimization problem

[TABLE] 2. 2.

Set $\gamma\coloneq\prod_{i\in[G]}\gamma_{i}$ (total negativity) and $M\coloneq\frac{1}{2}\gamma^{2}\log(2/\delta)\,\|O\|^{2}\varepsilon^{-2}$ (number of samples); 3. 3.

For $s\in[M]$:

1. (a)

(Classical) sample indices $j_{i,s}\in\big[|\widetilde{B}^{i}|\big]$ for each $i\in[G]$ according to the probability distribution $p_{i}(j)\coloneq\gamma_{i}^{-1}|b_{i,j}|$ defined by the solution to Eq. (LABEL:eq:lp) (here we include also $b_{i,|\widetilde{B}|+1}\coloneq 1$). Then compute the total sign $\sigma_{s}\coloneq\prod_{i\in[G]}\mathrm{sgn}(b_{i,j_{i,s}})$; 2. (b)

(Quantum) run the noisy circuit $\widetilde{\mathcal{C}}_{s}\coloneq\prod_{i\in[G]}\widetilde{\mathcal{B}}^{i}_{j_{i,s}}$ with $\rho$ as the initial state and measure $O$ on it. Using the measurement outcome $o_{s}\in\text{spectrum}(O)$, record a sample of the random variable $x_{s}\coloneq\gamma\,\sigma_{s}\,o_{s}$;

Return:* $\overline{O}\coloneq\frac{1}{M}\sum^{M}_{s=1}x_{s}$ (empirical mean of $x_{s}$).*

The optimization problem in Eq. (LABEL:eq:lp) has $v=|\widetilde{B}|$ variables and $c=256$ constraints (for $2$-qubit gates) corresponding to the matrix entries of $U_{i}$ and can be solved using standard convex optimization solvers [35] in time $\text{poly}\big(v,c\big)=\text{poly}\big(|\widetilde{B}|\big)$. We used cvxpy [36, 37] as the modeling interface to the mosek [38] solver. We chose to keep $b_{i,|\widetilde{B}|+1}\coloneq 1$ out of the optimization variables since this allows simpler proofs when bounding the total negativity in Lemma 9, and a clear choice of compilation error in step 1. (a). In our numerics, a solution in the largest basis $B_{1}$ of $\sim 11k$ elements took about $5$ minutes on a standard laptop, while for bases $B_{2}$ and $B_{3}$ the time was negligible. For the compilation step, the method supports any two-qubit gate synthesis algorithm. Here we illustrate it for a strategy as follows: we first decompose $U_{i}$ into CNOT, $\sqrt{X}$, and $R_{z}$ gates, and then use gridsynth [39] to synthesize each $R_{z}$ into Clifford + $T$ gates, from which we extract the sequence lengths $L_{i}$ appearing in Theorem 6. This is described in detail in App. E where, as a side result which may be of independent interest, we also numerically estimate the average sequence length $L$ of this compilation strategy for Haar random $2$-qubit unitaries. We observe a Solovay-Kitaev-like polylog scaling $L=c_{1}\,\log^{c_{2}}\big(1/\epsilon_{c}\big)$ with the diamond-norm error $\epsilon_{c}$, where $c_{1}\approx 210$ and $c_{2}\approx 0.75$. Improved compilation strategies such as probabilistic synthesis [21] are likely to yield better constants.

Remarkably, for a target precision $\varepsilon$, Alg. 5 solves Problem 1 with $\varepsilon$-independent logical circuit sizes and physical-qubit overhead. This is made possible by using a noisy compilation of each gate $U_{i}$ with $\varepsilon$-independent compilation error in Eq. (LABEL:eq:lp), which propagates to an $\varepsilon$-independent negativity $\gamma$ as shown by Lemma 9 in App. C. This is the core of our main theorem, proven in App. C:

Theorem 6** (CIPEC).**

Let $\{U_{i}\}_{i\in[G]},\varepsilon,\delta,\rho,O$ as in Problem 1, $\mathtt{Q}$ be a CIPEC-capable device with worst-case error $\epsilon_{\mathtt{Q}}$ and worst-case negativity $c_{*}$, and $\omega_{1},\omega_{2}>1$ be constants. Denote by $L_{i}$ the length of a logical-gate sequence that compiles $U_{i}$ up to diamond-norm error $\epsilon_{c,i}\leq\log(\omega_{1})/(2c_{*}G)$ into the universal set $V$ realized by $\mathtt{Q}$, and let $L\coloneq\sum_{i\in[G]}L_{i}$ be the total number of logical gates. Then, if $L\leq L_{\text{max}}\coloneq\log(\omega_{2})/(2c_{*}\epsilon_{\mathtt{Q}})$, Alg. 5 solves Problem 1 on $\mathtt{Q}$ for any target precision $1/\varepsilon$ with constant sample overhead $\gamma^{2}\leq\omega_{1}\omega_{2}$, i.e., using at most $\frac{1}{2}\gamma^{2}\,\norm{O}^{2}\varepsilon^{-2}\log(2/\delta)$ samples.

This result unlocks the possibility of solving instances of Problem 1 that cannot be solved with QEC alone or with standard PEC under the same assumptions given a fixed budget of quantum resources. More precisely, with QEC alone, after compilation up to precision $\epsilon_{c}=\mathcal{O}(\varepsilon)$, the circuit has $L_{\text{QEC}}(\varepsilon)$ logical gates, each of which having noise strength at most $\epsilon_{\mathtt{Q}}$. Intuitively, one expects any estimation built out of this noisy circuit to succeed only as long as $L_{\text{QEC}}(\varepsilon)\,\epsilon_{\mathtt{Q}}=\mathcal{O}(\varepsilon)$ (see Lem. 10 in App. C for a formal argument). To satisfy the latter, given a fixed physical error rate $p$, one should either restrict to logical-circuits of maximal size $L_{\text{QEC}}(\varepsilon)=\mathcal{O}(\varepsilon/\epsilon_{\mathtt{Q}})$ or scale the logical error rate down with $L_{\text{QEC}}(\varepsilon)$ as $\epsilon_{\mathtt{Q}}=\mathcal{O}(\varepsilon/L_{\text{QEC}}(\varepsilon))$ by increasing the code distance, and hence the physical-qubit overhead. For an error-mitigated strategy based on standard PEC, the negativity scales as $\gamma\leq e^{L_{\text{PEC}}(\varepsilon)\,\epsilon_{\mathtt{Q}}}$, which implies that Problem 1 is only solvable with constant sample overhead for $\varepsilon$ above a threshold value (see Lemma 12 in App. C). In contrast, by virtue of Thm. 6, CIPEC is limited by the $\varepsilon$-independent condition $L\,\epsilon_{\mathtt{Q}}\leq\log(\omega_{2})/(2c_{*})$, which is advantageous in high-precision regimes. Importantly, if this condition is not satisfied for a given $\omega_{2}$, CIPEC still solves Problem 1, but with a sample overhead now growing as $\gamma^{2}\leq e^{\frac{1}{2}\log(\omega_{1})+c_{*}\,L\,\epsilon_{\mathtt{Q}}}$ (see App. C). This is remarkable because there exist (see Corollary 11 in App. C) high-precision instances of Problem 1 that cannot be solved using QEC only or PEC even in the limit of infinite statistical samples. Moreover, apart from physical-qubit overhead, note that CIPEC enables also a logical-circuit depth reduction, with $L=\mathcal{O}(G\,\text{polylog}(G))$ being $\varepsilon$-independent while QEC and PEC require $\mathcal{O}(G\,\text{polylog}(G/\varepsilon))$ depth. Finally, although Thm. 6 assumes perfect noise characterization (cf. Eq (1)), in Thm. 13 in App. D we prove that CIPEC is stable against characterization errors, tolerating deviations in the characterization of the logical noise channels of up to diamond norm $\varepsilon/(2L\norm{O})$. We stress that CIPEC is agnostic to the choice of compiler, basis, and QEC code – any improvement in each of these is automatically inherited by CIPEC.

Application: Jones polynomial estimation.

As an illustration, we apply our framework to the problem of estimating the Jones polynomial $J_{K}(q)$ of a knot (or link) $K$ [40]. This knot invariant is a polynomial over the complex variable $q$ with integer coefficients determined solely by the knot topology. We focus on the estimation of $J_{K}(q)$ at the special point $q=e^{2\pi i/5}$ to relative precision $\varepsilon$, which is known to be BQP-complete under certain conditions on $\varepsilon$ and $K$ (see App. F for details). As shown in [23], this quantity is proportional to the matrix element $\expectationvalue{s|U_{\Sigma}|s}$, where $\ket{s}\coloneq\ket{0101\cdots 010}$ is a $(n_{s}+1)$-qubit computational basis state and $U_{\Sigma}$ is a unitary representation of the braid word $\Sigma$ on an even number $n_{s}$ of strands describing the knot $K$ with a plat closure [41]. We use the control-free Hadamard test quantum algorithm proposed in [23] (see App. F) to statistically estimate $\expectationvalue{s|U_{\Sigma}|s}$ in Eq. (11), which falls within the scope of Problem 1 with $U\rho\,U^{\dagger}=\outerproduct{s}{s}$, $G=\mathcal{O}(|\Sigma|,n)$ two-qubit gates, and $O=\frac{1}{2}(U_{\Sigma}^{\phantom{\dagger}}+U_{\Sigma}^{{\dagger}})$ or $O=\frac{1}{2i}(U_{\Sigma}^{\phantom{\dagger}}-U_{\Sigma}^{{\dagger}})$ for the real and imaginary parts of $\expectationvalue{s|U_{\Sigma}|s}$, respectively. Thus, we can solve it for a generic knot using CIPEC as long as the total circuit size $L=\mathcal{O}\big((|\Sigma|+n)\log^{c_{2}}(|\Sigma|+n)\big)\leq L_{\text{max}}=1/(2c_{*}\epsilon_{\mathtt{Q}})$.

In Fig. 3 we show the results for the trefoil knot assuming a noisy logical device realizing Clifford + $T$ operations with local depolarizing noise after each gate. The corresponding Hadamard test circuit can be rearranged in terms of $G=9$ ideal two-qubit gates $U_{i}$ (see App. F), which were then compiled into the native Clifford + $T$ gateset using the gridsynth algorithm [39]. The depolarizing strengths were chosen as $10^{-6}$ for single-qubit gates and state preparation channels, and $10^{-5}$ for CNOT and $T$ gates, which imply $\epsilon_{\texttt{Q}}=10^{-5}$. The target relative precision (or, equivalently, the absolute precision for $\expectationvalue{s|U_{\Sigma}|s}$ in Eq. (11)) was set to $\varepsilon=10^{-2}$ and we used the minimal basis $\widetilde{B}_{2}$. We compare CIPEC (with $\omega_{1}=\omega_{2}=e$) with the QEC-only strategy of Lemma 5, SI-II (with $\xi=\eta=3$). The compilation errors are $\epsilon_{c,\text{CIPEC}}\approx 3.5\times 10^{-4}$ and $\epsilon_{c,\text{QEC}}\approx 3.7\times 10^{-4}$, giving a circuit size of $L=3902$ for QEC and $L=3962$ for CIPEC (in the worst-case). The former violates the feasibility condition $L\leq L_{\text{max}}$, while the latter does not. The total number of samples was $S_{\text{CIPEC}}\approx 5.9\times 10^{5}$ (corresponding to $\gamma^{2}\approx 2.46$) and for QEC we used $S_{\text{QEC}}\approx 1.07\times 10^{6}$ as in Lemma 5. We see that CIPEC is able to deliver an estimate within the desired precision while QEC alone is not, even when allowing QEC to collect extra samples. This is in agreement with the $\varepsilon$-independence of the CIPEC overheads and the feasibility conditions illustrated in Fig. 2 for $\epsilon_{\mathtt{Q}}=10^{-5}$.

Closing remarks.

The ability to estimate expectation values using quantum resources that are independent of the target precision can unlock important applications where high precision is required. Apart from the estimation of Jones polynomials considered here, other relevant use cases may be found, e.g., in chemistry and materials science. For instance, for ground-state energy estimation of molecules one typically requires chemical accuracy, which implies that $1/\varepsilon$ grows with the number of spin orbitals. For QEC-only strategies, this in turn implies that both the size of the compiled logical circuit and the code-distance (hence also the physical-qubit count) must explicitly grow with the number of orbitals as well. In contrast, CIPEC does the job with both quantities insensitive to $\varepsilon$, incurring a moderate sample overhead $\gamma^{2}$. Ultimately, of course, any statistical estimation strategy (including CIPEC and QEC alone) will fail for a task demanding too-high precision, given that the total sample complexity (i.e., number of runs) unavoidably grows with $1/\varepsilon^{2}$. However, for early fault-tolerant hardware, sample complexity is a far more abundant resource than logical-circuit complexity or code-distance. This means that there is a precision regime achievable only by CIPEC (see Lem. 10 and Corollary 11, in App. C). With that in mind, our findings offer a practical route towards fault-tolerant quantum computation with precision-independent overheads.

Acknowledgements.

We thank Ryuji Takagi, Suguru Endo, Fernando G. S. L. Brandão, Samson Wang, Ingo Roth, and Ariel Bendersky for discussions. D.S.F. acknowledges financial support from the Novo Nordisk Foundation (Grant No. NNF20OC0059939 Quantum for Life)

Appendix A Logical gate characterization.

To characterize the implementable operations in $\widetilde{\mathcal{A}}$, we can use standard gate set tomography (GST) methods [44, 45]. $\widetilde{A}$ includes noisy gates $\widetilde{V}$ from a universal set and the non-unitary channels $\widetilde{\Gamma}$. The former are well described by ideal unitary gates followed by stochastic Pauli noise, while the latter can be more general. For concreteness, here we focus on a device that can implement CNOT and tensor-product gates $\{I,H,S,S^{\dagger},T\}^{\otimes 2}$ on any pair of logical qubits connected by the fault-tolerant architecture in question (Pauli gates ($X,Y,Z$) can be implemented without noise by changing the Pauli frame [46]). The gate-set size is then $|\widetilde{A}|=26\tau$, where $\tau$ is a logical-architecture dependent factor that, in the foreseen fault-tolerant architectures, is expected to scale as $\tau=\mathcal{O}(n)$. For example, in superconducting-qubit platforms [47], the scaling follows from local circuit connectivity. In turn, for platforms with all-to-all connectivity due to movable qubits, such as Rydberg atoms [48], 2-qubit gates will be executed at few entangling zones while single-qubit ones will be executed at each logical qubit’s position, giving again $\tau=\mathcal{O}(n)$. Standard GST requires $\mathcal{O}(1/\epsilon_{\text{char}}^{2})$ uses to estimate a two-qubit channel up to diamond norm error $\epsilon_{\text{char}}$ and, as shown in Thm. 7 in the SI, CIPEC requires $\epsilon_{\text{char}}\leq\varepsilon/(2L\|O\|)$ for each gate. Putting everything together, this implies a total of $\mathcal{O}\left(\tau L^{2}\|O\|^{2}/\varepsilon^{2}\right)$ runs, which is prohibitive for practical applications. Fortunately, this scaling can be improved using long-sequence GST [49]. Empirically, it has been verified that, using a set of $\mathcal{O}(\log(L_{\text{GST}}))$ sequences of gates with maximum depth $\mathcal{O}(L_{\text{GST}})$ and collecting $S$ measurement samples for each, the diamond norm error in the estimated gates scales as $\epsilon_{\text{char}}=\mathcal{O}(1/(L_{\text{GST}}\sqrt{S}))$. Therefore, the time required by long-sequence GST to characterize a given 2-qubit gate set is $T_{\text{char}}=\mathcal{O}(\log(L_{\text{GST}})\,S\,L_{\text{GST}})$ per logical qubit pair. The larger the sequence size, the more precise the estimation, as long as the hardware allows it (i.e., $L_{\text{GST}}=\mathcal{O}(1/\epsilon_{\mathtt{Q}})$). Therefore, considering sequences of maximum depth $L_{\text{GST}}=L$ and $S=\mathcal{O}(\tau/\epsilon^{2})$ samples yields a time $T_{\text{char}}=\mathcal{O}(\tau L\log(L)/\varepsilon^{2})$ to characterize all the implementable operations up to the precision required to guarantee stability of CIPEC. Moreover, in case the platform allows for parallel execution of 2-qubit gates, this scaling can be further reduced to $T_{\text{char}}=\mathcal{O}({(\tau/n)}L\log(L)/\varepsilon^{2})$. For instance, for 1D connectivity, the universal set can be split into $n-1$ gate sets, which can be characterized in only two separate experiments. In turn, for the 2D square lattice, gates on $2(n-\sqrt{n})$ pairs of qubits have to be characterized, which can be done in $4$ experiments by first scanning row connections and then column ones. To assess the feasibility of $T_{\text{char}}$, we compare it against CIPEC’s observable estimation time, namely $T_{\text{est}}=\mathcal{O}\big((L/n)\gamma^{2}\|O\|^{2}/\varepsilon^{2}\big)$, since the largest circuit sampled has $L$ gates and up to $\mathcal{O}(n)$ of them can be run in parallel. We thus see that $T_{\text{char}}/T_{\text{est}}=\mathcal{O}\big(\tau\log(L)/(\gamma^{2}\|O\|^{2})\big)$, which scales linearly with the number of qubits ($\tau=\mathcal{O}(n)$) and mildly with $L$. Therefore, characterization can be done with a modest runtime overhead compared to CIPEC’s own runtime.

Moreover, if the error-correction budget is sufficiently large, CIPEC can completely bypass the need for logical noise characterization. This is because (here $I$ is the two-qubit identity channel)

[TABLE]

Therefore, if resources allow a device with noise $\epsilon_{\mathtt{Q}}\leq\varepsilon/(2L\|O\|)$, the stability condition $\epsilon_{\text{char}}\leq\varepsilon/(2L\|O\|)$ is trivially ensured by characterizing the noise channels as identity channels – in other words, the quasi-probability decomposition in Eq. (5) can be done entirely with respect to noise-free gates and basis elements. Moreover, in this regime, the overhead in error correction is strictly smaller than that of standard QEC, and our circuits are shallower than using regular compilation (since we only need to compile to $\varepsilon$-independent precision $\epsilon_{c}=\mathcal{O}(1/G)$). Note that, in this case, the mitigation is of compilation errors, similar to the proposal in [22] based on logical-level PEC with the error inversion method.

Appendix B Building a basis of 2-qubit channels from implementable operations

Here we show how to build a basis $\widetilde{B}$ of channels in Def. 3 using the set of implementable operations of a noisy logical quantum device. First, consider a two-qubit system $AB$, where $A$ labels the first qubit and $B$ the second. The linear span of the set of all CP maps in $AB$ in its Choi representation is the real vector space $\text{Herm}(AB)$ of Hermitian operators acting on the composite system, whose dimension is $4^{4}=256$. The linear span of the set of all CPTP maps, on the other hand, is given by the subspace of all $X\in\text{Herm}(AB)$ with $\tr_{B}X=\openone_{A}$, which has dimension $4^{4}-4^{2}+1=241$ (see, e.g., [50]). Therefore, given a set of two-qubit (CPTP or CP) noisy channels, one can test whether it is a spanning set for the span of CPTP or CP by computing the rank of the Gram matrix built from its elements under the Hilbert-Schmidt inner product, defined by $\tr(XY)$ for all $X,Y\in\text{Herm}(AB)$.

Given a CIPEC-capable noisy quantum device $\mathtt{Q}$ as in Def. 3, it is always possible to build a set $\widetilde{B}$ of noisy channels of the form $\widetilde{V}^{(D)}\cup\,\widetilde{\Gamma}$ that passes the Gram matrix test, where we recall that $\widetilde{V}^{(D)}$ denotes a set of sequences of noisy universal gates from $\mathtt{Q}$ of length up to $D$. However, these are far from unique, and some sets of channels might be preferable given specific conditions. To assess how suited one set of channels will be for a particular application, we need to identify relevant figures of merit. From the description of the CIPEC algorithm (Alg. 5), it is natural to worry about: $i)$ the amount of classical pre-processing required to build and store the classical description of all distinct circuits that need to be estimated on the device; $ii)$ the length (or depth) of each one of those circuits and; $iii)$ the total number of samples required for the estimation task. These concerns lead us to identify three figures of merit for a particular set $\widetilde{B}$, respectively: the first is the number of elements $|\widetilde{B}|$, which affects the number of distinct circuits that need to be implemented on a particular task; second, the maximal length $D$ of the elements of $\widetilde{B}$ affects the total length and depth of each distinct circuit; and lastly, the maximal negativity $c_{*}$ of a unitary channel decomposed with respect to $\widetilde{B}$ indirectly affects the total number of samples required for the estimation. Although it is not hard to see why the first two figures of merit address their respective concerns, the last one requires further explanation.

Let $\widetilde{B}$ be a spanning set for the linear span of CPTP (or CP) channels, and let $c_{*}$ be as given in Def. 4. The following lemma shows that $c_{*}$ is the maximal amplification factor required for $\|\mathcal{C}\|_{\diamond}$ of an arbitrary channel $\mathcal{C}$ to bound negativity $\|\mathcal{C}\|_{\widetilde{B}}$, which is related to the sample complexity overhead of error mitigation.

Lemma 7** (Negativity and the diamond norm).**

Given a basis $\widetilde{B}$ for the vector space of two-qubit maps, its associated negativity satisfies

[TABLE]

Proof.

Given two norms $\norm{\cdot}_{\widetilde{B}}$ and $\norm{\cdot}_{\diamond}$ on a finite-dimensional vector space, they must be connected within constant factors of one another. Specifically, there exist two real numbers $0<c_{\text{low}}\leq c_{\text{up}}\coloneq c_{*}$ such that, for all $\mathcal{C}$ in this vector space, we have

[TABLE]

This inequality is trivially true for $\mathcal{C}=0$. For $\mathcal{C}\neq 0$, we can restrict to the case where $\norm{\mathcal{C}}_{\diamond}=1$, obtaining

[TABLE]

By the extreme value theorem, a continuous function (in this case $\norm{\cdot}_{\widetilde{B}}$) on a compact set (the closed and bounded set of unitary channels defined by $\norm{\mathcal{C}}_{\diamond}=1$) must achieve a maximum and minimum value on this set, therefore

[TABLE]

∎

Appendix C Main proofs

Here we prove our main statements, namely Theorem 6 and a similar statement (Lemma 10 below) for a naive estimator that leverages partial quantum error-correction only without quantum error mitigation. For that, we need to introduce two auxiliary lemmas. The following lemma gives the diamond norm error of two-qubit gate synthesis after replacing ideal gates from the universal set $V$ by noisy gates from the set of implementable operations $\widetilde{V}$.

Lemma 8** (Two-qubit Unitary Synthesis using noisy channels).**

Let $\mathcal{U}$ be a two-qubit unitary gate, $\mathcal{V}_{c}^{(L)}$ be a synthesis of $\mathcal{U}$ with $L$ ideal gates from the set $V$, and $\epsilon_{c}$ denotes the compilation error $\big\|\mathcal{U}-\mathcal{V}_{c}^{(L)}\big\|_{\diamond}\leq\epsilon_{c}$. Let $\widetilde{\mathcal{V}}_{c}^{(L)}$ be the version of $\mathcal{V}_{c}^{(L)}$ with $L$ noisy operations from the implementable set $\widetilde{V}$. Then $\big\|\mathcal{U}-\widetilde{\mathcal{V}}_{c}^{(L)}\big\|_{\diamond}\leq\epsilon_{c}+L\,\epsilon_{\mathtt{Q}}$.

Proof.

[TABLE]

where $\big\|\mathcal{U}-\mathcal{V}_{c}^{(L)}\big\|_{\diamond}\leq\epsilon_{c}$ by assumption, and $\big\|\mathcal{V}_{c}^{(L)}-\widetilde{\mathcal{V}}_{c}^{(L)}\big\|_{\diamond}\leq L\,\epsilon_{\mathtt{Q}}$ follows by repeated use of the inequality $\|\mathcal{C}_{1}\mathcal{C}_{0}-\mathcal{D}_{1}\mathcal{D}_{0}\|_{\diamond}\leq\|\mathcal{C}_{0}-\mathcal{D}_{0}\|_{\diamond}+\|\mathcal{C}_{1}-\mathcal{D}_{1}\|_{\diamond}$. ∎

This result can then be used to prove the following upper bound on the negativity of CIPEC, which, in turn, is the main ingredient for proving Theorem 6.

Lemma 9** (CIPEC negativity upper bound).**

Let $\epsilon_{\mathtt{Q}}$ as in Eq. (1), $\mathcal{U}_{i}$ a two-qubit unitary channel, and $\mathcal{U}_{i}=\widetilde{\mathcal{V}}_{c,i}^{(L_{i})}+\sum_{j\in[|\widetilde{B}|]}b_{i,j}\widetilde{\mathcal{B}}_{j}$ be its compilation-informed quasi-probability decomposition Eq. (LABEL:eq:lp). Then the total negativity of $\mathcal{U}=\prod_{i\in[G]}\mathcal{U}_{i}$ is $\gamma\leq e^{\lambda}$ with $\lambda:=\sum_{i\in[G]}c_{*}\big[\epsilon_{c,i}+L_{i}\,\epsilon_{\mathtt{Q}}\big]$.

Proof.

For each $i\in[G]$ the minimum negativity satisfies $\gamma_{i}=1+\norm{\mathcal{U}_{i}-\widetilde{\mathcal{V}}_{c,i}^{(L_{i})}}_{\widetilde{B}}\leq 1+c_{*}\norm{\mathcal{U}_{i}-\widetilde{\mathcal{V}}_{c,i}^{(L_{i})}}_{\diamond}$, according to Lemma 7. Using Lemma 8 we get $\gamma_{i}\leq 1+\lambda_{i}$ with $\lambda_{i}=c_{*}[\epsilon_{c,i}+L_{i}\,\epsilon_{\mathtt{Q}}]$, and hence $\gamma=\prod_{i\in[G]}\gamma_{i}\leq\prod_{i\in[G]}(1+\lambda_{i})\leq e^{\lambda}$ with $\lambda\coloneq\sum_{i\in[G]}\lambda_{i}$. ∎

With the results above, we can finally prove our main statement, Theorem 6, which we restate here for convenience

Theorem 6.

Let $\{U_{i}\}_{i\in[G]},\varepsilon,\delta,\rho,O$ as in Problem 1, $\mathtt{Q}$ be a CIPEC-capable device with worst-case error $\epsilon_{\mathtt{Q}}$ and worst-case negativity $c_{*}$, and $\omega_{1},\omega_{2}>1$ be constants. Denote by $L_{i}$ the length of a logical-gate sequence that compiles $U_{i}$ up to diamond-norm error $\epsilon_{c,i}\leq\log(\omega_{1})/(2c_{*}G)$ into the universal set $V$ realized by $\mathtt{Q}$, and let $L\coloneq\sum_{i\in[G]}L_{i}$ be the total length of the logical-circuit number of logical gates. Then, if $L\leq L_{\text{max}}\coloneq\log(\omega_{2})/(2c_{*}\epsilon_{\mathtt{Q}})$, Alg. 5 solves Problem 1 on $\mathtt{Q}$ for any target precision $1/\varepsilon$ with constant sample overhead $\gamma^{2}\leq\omega_{1}\omega_{2}$, i.e., using at most $\frac{1}{2}\gamma^{2}\,\norm{O}^{2}\varepsilon^{-2}\log(2/\delta)$ samples.

Proof.

Since, by construction, the quasiprobability decomposition in Eq. (LABEL:eq:lp) is an unbiased estimator for each $\mathcal{U}_{i}$, the total estimation error of the expectation value in Alg. 5 comes only from statistics and can be arbitrarily decreased by collecting more samples. In particular, a target error $\varepsilon$ requires a statistical accuracy $\epsilon_{s}=\varepsilon/\gamma$, where $\gamma$ is the total negativity of the quasiprobability decomposition of $\mathcal{U}$. By Hoeffding inequality, this can be achieved (with probability $1-\delta$) by collecting

[TABLE]

samples, where $\|O\|$ is the spectral norm of $O$, which gives the range of the random variable under our assumption of measuring in the eigenbasis of $O$. In Lemma 9, we showed that $\gamma\leq e^{\lambda}$ with $\lambda:=\sum_{i\in[G]}c_{*}\big[\epsilon_{c,i}+L_{i}\,\epsilon_{\mathtt{Q}}\big]$. This can be made constant, more precisely $\gamma\leq\sqrt{\omega_{1}\omega_{2}}=\mathcal{O}(1)$, by ensuring $\lambda\leq\log(\sqrt{\omega_{1}\omega_{2}})$ as follows: $i)$ compile each $\mathcal{U}_{i}$ to precision $\epsilon_{c,i}\leq\log(\sqrt{\omega_{1}})/(c_{*}G)$ such that $\sum_{i\in[G]}c_{*}\,\epsilon_{c,i}\leq\log(\sqrt{\omega_{1}})$; this fixes the total sequence length $L\coloneq\sum_{i\in[G]}L_{i}$; $ii)$ the condition $\lambda\leq\log(\sqrt{\omega_{1}\omega_{2}})$ then becomes $c_{*}\,L\,\epsilon_{\mathtt{Q}}\leq\log(\sqrt{\omega_{2}})$, which is true by assumption. ∎

Lemma 10** (Quantum resource requirements for solving Problem 1 using QEC only).**

Let $\eta>0$ and $\xi>1$ be given constants. In the same setup of Theorem 6, compile each $\mathcal{U}_{i}$ into a sequence $\mathcal{V}^{(L_{i})}_{c,i}$, with compilation error $\epsilon_{c,i}\leq\varepsilon/(\eta G)$, such that $\mathcal{V}^{(L)}_{c}=\prod_{i\in[G]}\mathcal{V}^{(L_{i})}_{c,i}$ is a compilation of $\mathcal{U}$ of total sequence length $L\coloneq\sum_{i\in[G]}L_{i}$. Then, statistically estimate $\operatorname{tr}\big(O\,\widetilde{\mathcal{V}}^{(L)}_{c}(\rho)\big)$, where $\widetilde{\mathcal{V}}^{(L)}_{c}$ is the noisy version of $\mathcal{V}^{(L)}_{c}$, using $\xi^{2}\|O\|^{2}\log(2/\delta)/(2\varepsilon^{2})$ samples. This algorithm solves Problem 1 if and only if

[TABLE]

Proof of Lemma 10.

Repeated use of Lemma 8 followed by application of the inequality $\|\mathcal{C}_{1}\mathcal{C}_{0}-\mathcal{D}_{1}\mathcal{D}_{0}\|_{\diamond}\leq\|\mathcal{C}_{0}-\mathcal{D}_{0}\|_{\diamond}+\|\mathcal{C}_{1}-\mathcal{D}_{1}\|_{\diamond}$ allow us to conclude that

[TABLE]

Therefore, if $y$ is the random variable representing output of measuring $O$ on the noisy state $\widetilde{\mathcal{V}}^{(L)}_{c}(\rho)$, and $\overline{y}$ is its empirical mean, the total estimation error can written as

[TABLE]

or in other words, $\overline{y}$ is an unbiased estimator for $\tr(O\mathcal{V}^{(L)}_{c}(\rho))$, but a biased estimator for $\tr(O\,\mathcal{U}(\rho))$ with bias given by $\sum_{i\in[G]}(L_{i}\epsilon_{\mathtt{Q}}+\epsilon_{c,i})$. Now let $\epsilon_{s}$ be an upper bound for the estimation error of $\tr(O\mathcal{V}^{(L)}_{c}(\rho))$. Then, in order to have $|\overline{y}-\tr(O\,\mathcal{U}(\rho))|\leq\varepsilon$ we need the following conditions: $i)$ for every $i\in[G]$, set the compilation error to $\epsilon_{c,i}\leq\varepsilon/(\eta\,G)$ – this fixes the total sequence length $L\coloneq\sum_{i\in[G]}L_{i}$; $ii)$ collect enough statistics so as to make the statistical error $\epsilon_{s}=\varepsilon/\xi$ – using Hoeffding’s inequality, this can be achieved with success probability $1-\delta$ using

[TABLE]

samples; $iii)$ then, imposing $\varepsilon\,(1/\eta+1/\xi)+L\,\epsilon_{\mathtt{Q}}\leq\varepsilon$ we obtain the remaining condition given by the lemma, which concludes the proof. ∎

Corollary 11**.**

Under the conditions of Lemma 10 and given a fixed $\eta>0$, if $L>(1-\frac{1}{\eta})\frac{\varepsilon}{\epsilon_{\mathtt{Q}}}$ then correctness of the QEC-only strategy in solving Problem 1 cannot be guaranteed even in the limit of infinitely many samples. Also, given a fixed precision and a fixed amount of samples, there must be a maximal $\eta$ (and therefore a maximal $L$) for which Problem 1 can be guaranteed to be solvable in this device.

Proof.

The first conclusion follows directly from Lemma 10 by taking the limit $\xi\rightarrow\infty$, while the second follows from Lemma 10 by the fact that $L$ grows as $\eta$ grows regardless of the compilation algorithm, which implies that the left-hand side of the feasibility condition grows without bounds, while the right-hand side is eventually constant. ∎

The following Lemma shows that another strategy for solving Problem 1 based on standard PEC after compiling the circuit into ideal gates from the universal set $V$ cannot succeed for arbitrary $\varepsilon$.

Lemma 12** (Solving Problem 1 with standard PEC).**

In the same setup of Theorem 6, compile each $U_{i}$ into gates from the universal set $V$ up to precision $\epsilon_{c,i}\leq\varepsilon/(2G)$, and denote by $U_{c}\coloneqq\prod_{i\in[G]}\prod_{\ell\in[L_{i}]}V_{i_{\ell}}$ the full compiled circuit, with $L_{i}=c_{1}\log^{c_{2}}(2G/\varepsilon)$. Then statistically estimate $\Tr(OU_{c}^{\dagger}\rho U_{c})$ up to precision $\epsilon_{s}=\varepsilon/2$ using standard PEC in the basis $\widetilde{B}$, i.e., decompose each unitary channel $\mathcal{V}_{i_{\ell}}$ as $\mathcal{V}_{i_{\ell}}=\sum_{j\in[|\widetilde{B}|]}a_{j}\widetilde{B}_{j}$ and sample basis elements from the distribution $p_{j}\coloneqq|a_{j}|/\|a\|_{1}$. This strategy solves Problem 1 only if $\varepsilon\geq 2G\,e^{-(c_{1}G\epsilon_{\mathtt{Q}})^{-1/c_{2}}}$.

Proof.

The compilation error $\epsilon_{c}=\varepsilon/(2G)$ for each $U_{i}$ ensures that $\|\mathcal{U}-\mathcal{U}_{c}\|_{\diamond}\leq\varepsilon/2$. The compiled circuit has a total size $L=\sum_{i\in[G]}L_{i}=G\,c_{1}\log^{c_{2}}(2G/\varepsilon)$. By applying standard PEC independently to each ideal $\mathcal{V}_{i_{\ell}}$ as $\mathcal{V}_{i_{\ell}}=\sum_{j\in[|\widetilde{B}|]}a_{j}\widetilde{B}_{j}$, the negativity of each term should be trivially controlled assuming $\widetilde{\mathcal{V}}_{i_{\ell}}$ (the noisy version of the gate $\mathcal{V}_{i_{\ell}}$ we are trying to expand) is included in the basis, with $\|\mathcal{V}_{i_{\ell}}-\widetilde{\mathcal{V}}_{i_{\ell}}\|_{\diamond}\leq\epsilon_{\mathtt{Q}}$ (see Eq. (1)). This is a reasonable assumption, given that the noisy gate itself is one of the simplest implementable operations. As a result, the total negativity is $\gamma\leq(1+\epsilon_{\mathtt{Q}})^{L}\leq e^{L\epsilon_{\mathtt{Q}}}$, which can be made constant, $\gamma\leq e$, as long as $L\,\epsilon_{\mathtt{Q}}\leq 1$. Since $L=G\,c_{1}\log^{c_{2}}(2G/\varepsilon)$ is fixed, in terms of $\varepsilon$ this means that the strategy only succeeds for $\varepsilon\geq 2G\,e^{-(c_{1}G\epsilon_{\mathtt{Q}})^{-1/c_{2}}}$. ∎

As an illustration, for a device with $\epsilon_{\mathtt{Q}}=10^{-5}$ and a circuit with $G=10^{2}$ two-qubit gates before compilation, using our estimated $c_{1}=210$ and $c_{2}=0.75$ from Sec. IV of the SI, Lemma 12 gives $\varepsilon\geq 0.07$, while CIPEC succeeds for any $\varepsilon$.

Appendix D Stability analysis

Here we analyze the stability of CIPEC under noise characterization errors. Recall that the implementable operations in Eq. (1) have the form $\widetilde{\mathcal{A}}_{j}=\mathcal{E}_{j}\circ\mathcal{A}_{j}$ where $\mathcal{E}_{j}$ are the noise channels affecting the device. It is important to notice that we will never be able to perfectly characterize each $\mathcal{E}_{j}$ – the best one can hope for are estimates ${\mathcal{E}}^{\prime}_{j}$ obtained via gate-set tomography [44, 51]. Therefore, we must ensure that our method is stable to characterization errors. This is proven in the following theorem. In particular, solving Problem 1 with a target precision $\varepsilon$ using CIPEC requires characterizing the error channels to diamond-norm precision $\varepsilon/(2L\norm{O})$ and running CIPEC with a target error $\varepsilon/2$.

Theorem 13** (Stability of CIPEC).**

Let $\langle O\rangle$ be the exact expectation value, $\overline{O}^{\,\prime},\overline{O}$ be the outputs of Alg. 5 using, respectively, the characterized noise channels $\widetilde{\mathcal{A}}^{\prime}_{j}=\mathcal{E}^{\prime}_{j}\circ\mathcal{A}_{j}$ and the true noisy channels $\widetilde{\mathcal{A}}_{j}=\mathcal{E}_{j}\circ\mathcal{A}_{j}$, and $\epsilon_{\text{char}}\coloneqq\max_{j\in[|\mathcal{A}|]}\norm{\mathcal{E}^{\prime}_{j}-{\mathcal{E}}_{j}}_{\diamond}$ be the maximum characterization error. Then

[TABLE]

Proof.

Here for simplicity we denote by $\alpha_{i}=\mathbf{b}_{i}\cup\{1\}$ the extended vector of coefficients in Eq. (2) that includes also the unit coefficient of the compilation term. Let $p(\alpha)=\prod_{i=1}^{G}p_{i}(\alpha_{i})$ be the distribution defined by the coefficients in Alg. 5 with respect to the perfectly characterized $\mathcal{E}_{j}$, and $\prod_{i\in[G]}\widetilde{\mathcal{B}}^{i}_{\alpha_{i}}$ be the corresponding sampled circuit, with $\widetilde{\mathcal{B}}^{i}_{\alpha_{i}}$ denoting the element of the augmented basis $\widetilde{B}^{i}$ with coefficient $\alpha_{i}$ in the decomposition of $U_{i}$. Similarly, denote by $\prod_{i\in[G]}\widetilde{\mathcal{B}}^{\prime\,i}_{\alpha_{i}}$ the corresponding circuit obtained by replacing each $\mathcal{E}_{j}$ by the imperfectly characterized noise channel $\mathcal{E}^{\prime}_{j}$. It follows from repeated use of the diamond norm inequality $\|\mathcal{C}_{1}\mathcal{C}_{0}-\mathcal{D}_{1}\mathcal{D}_{0}\|_{\diamond}\leq\|\mathcal{C}_{0}-\mathcal{D}_{0}\|_{\diamond}+\|\mathcal{C}_{1}-\mathcal{D}_{1}\|_{\diamond}$ that

[TABLE]

where the second inequality uses the fact that each $\widetilde{\mathcal{B}}^{i}_{\alpha_{i}}$ contains a sequence of either $L_{i}$ or $D<L_{i}$ implementable operations $\widetilde{\mathcal{A}}_{j}=\mathcal{E}_{j}\circ\mathcal{A}_{j}$ and $\norm{\mathcal{A}_{j}}_{\diamond}=1$. As a result,

[TABLE]

Given that each statistical run outputs a random variable with absolute value at most $\|O\|$, it follows that $\big|\overline{O}^{\,\prime}-\overline{O}\big|\leq L\|O\|\epsilon_{\text{char}}$. Since $\left|\overline{O}-\langle O\rangle\right|\leq\varepsilon$, the claim follows immediately by triangle inequality,

[TABLE]

In particular, $\big|\overline{O}^{\,\prime}-\langle O\rangle\big|\leq 2\varepsilon$ for $\epsilon_{\text{char}}\leq\frac{\varepsilon}{L\norm{O}}$. ∎

Thus, we see that the protocol is stable with respect to imperfect characterization. Furthermore, if the (logical) noise in the device does not change over time, we only need to perform this characterization once to implement all circuits having a number of gates $G$ and target precision $\varepsilon$ of the same order.

Appendix E Estimates for $c_{1}$ and $c_{2}$

Here we study the sequence length $L$ required for the Clifford+T compilation of Haar random two-qubit unitaries to a prescribed diamond norm precision $\epsilon_{c}$. Assuming a Solovay-Kitaev-like power-log scaling for the average case sequence length,

[TABLE]

the goal is to obtain numerical estimates for the constants $c_{1}$ and $c_{2}$. To the best of our knowledge, there is no straightforward algorithm to synthesize arbitrary two-qubit unitaries using the Clifford+T gateset with an end-to-end scaling analysis in diamond norm error. For that reason, we adopt a heuristic strategy described below.

We draw random unitaries $U_{\text{Haar}}$ from the Haar distribution and, for each of them, proceed as follows: $i)$ decompose it in terms of $15\,R_{z}$, $10\,\sqrt{X}$, and $3\,\mathrm{CNOT}$ gates as illustrated in Fig. 4 using the qiskit transpiler [52]; $ii)$ sample a target error upper bound $\epsilon^{\prime}_{c}$ from the uniform distribution in the interval $\big[\log\log 10^{-2},\log\log 10^{-9}\big]$; $iii)$ synthesize each $\mathrm{R}_{z}$ into Clifford+T gates to precision $\epsilon^{\prime}_{c}/15$ using the gridsynth algorithm [39] and each $\sqrt{\mathrm{X}}$ by its exact Clifford decomposition $\sqrt{\mathrm{X}}=e^{i\pi/4}\mathrm{S}^{\dagger}\mathrm{H}\mathrm{S}^{\dagger}$; $iv)$ reduce the resulting word using trivial Clifford+T identities (e.g., $\mathrm{S}^{\dagger}\mathrm{S}=\mathrm{S}\mathrm{S}^{\dagger}=\openone$) and record the word length $L$ as well as the resulting diamond norm compilation error $\epsilon_{c}\leq\epsilon^{\prime}_{c}$ returned by gridsynth. Since gridsynth is non-deterministic, we repeat steps $iii)$ and $iv)$ above $50$ times and keep only the pair $(L,\epsilon_{c})$ corresponding to the smallest $L$ (i.e., best-case synthesis of $U_{\text{Haar}}$). Fig. 5 shows a scatter plot of $L$ vs. $\log(1/\epsilon_{c})$ (in log-log scale) for $10^{4}$ Haar random unitaries, from which the constants $c_{1}$ and $c_{2}$ (for the average case) can be estimated by fitting a linear model, from which we obtain $c_{1}=210.36$ and $c_{2}=0.75$.

Appendix F Quantum circuits for Jones polynomial estimation

Here we briefly review some basic properties of knots and the control-free Hadamard test circuits introduced in [23] for estimating the Jones polynomial [40] of a knot on a quantum computer.

Recall that any knot can be represented as the closure of some braid $\Sigma$ [40, 53] on a number $n_{s}$ of strands. Any such braid can be expressed as a product of braid generators $\sigma_{i}^{\pm 1}$ acting on the $(i,i+1)$ pair of strands, where $i\in[n_{s}-1]$ and the number $|\Sigma|$ of terms in the product gives the total number of braid crossings. These braids can be closed in multiple ways, the most used ones being the so-called Markov and plat closures, $\text{M}(\Sigma)$ and $\text{P}(\Sigma)$ respectively [54, 53]. The former exists for any number of strands and is constructed by connecting the endpoints of each strand in $\Sigma$ through simple non-crossing arcs; the latter, on the other hand, requires a braid $\Sigma$ with an even number of strands and is obtained by connecting pairs of endpoints sequentially (see Fig. 1 in [23] for an illustration). For instance, the trefoil knot illustrated in Fig. 3 can be represented either as $\text{M}(\Sigma)$ for the $2$-strand braid $\Sigma=\sigma_{1}^{3}$ or as $\text{P}(\Sigma)$ for the $4$-strand braid $\Sigma=\sigma_{1}\,\sigma_{2}^{3}\,\sigma_{1}^{-1}$.

The Jones polynomial of a knot (or link) $K$, denoted $J_{K}(q)$, is a polynomial over a complex variable $q$ with integer coefficients determined solely by the topology of $K$. We focus on the estimation of $J_{K}(q)$ at the special point $q=e^{2\pi i/5}$ to relative precision $\varepsilon$. For $K=\text{M}(\Sigma)$ and $K=\text{P}(\Sigma)$ this problem is known to be DQC1-complete [55] and BQP-complete [56], respectively, for $\varepsilon=1/\mathcal{O}(\text{poly}(n_{s}))$ and $|\Sigma|=\mathcal{O}(\text{poly}(n_{s}))$. For simplicity, here we focus on the BQP-complete case. As shown in [23], the Jones polynomial can be expressed as

[TABLE]

where $\phi\coloneq\frac{1}{2}(1+\sqrt{5})$ is the golden ratio, $\omega_{\Sigma}\coloneq\sum_{\sigma\in\Sigma}\text{sign}(\sigma)$ is called the writhe of the braid $\Sigma$, $\ket{s}\coloneq\ket{0101\cdots 010}$ is $(n_{s}+1)$-qubit a computational basis state, and $U_{\Sigma}$ is the Fibonacci unitary representation of the braid $\Sigma$ [41] (see Fig. 4 of [23] for a specific realization in terms of single- and two-qubit gates). We use the quantum algorithm proposed in [23], which is based on a Monte Carlo estimation of the matrix element $\expectationvalue{s|U_{\Sigma}|s}$ in Eq. (11) using a control-free Hadamard test. The associated circuit $U$ (see Fig. 5 of Ref. [23]) contains $3|\Sigma|$ parametrized two-qubit gates, $\mathcal{O}(n)$ controlled-NOT gates, and $\mathcal{O}(|\Sigma|)$ single-qubit rotations. This falls within the scope of Problem 1 with $U\rho\,U^{\dagger}=\outerproduct{s}{s}$, $G=\mathcal{O}(|\Sigma|,n)$ two-qubit gates, and $O=\frac{1}{2}(U_{\Sigma}^{\phantom{\dagger}}+U_{\Sigma}^{{\dagger}})$ or $O=\frac{1}{2i}(U_{\Sigma}^{\phantom{\dagger}}-U_{\Sigma}^{{\dagger}})$ for the real and imaginary parts of $\expectationvalue{s|U_{\Sigma}|s}$, respectively. After estimating $\expectationvalue{s|U_{\Sigma}|s}$, the Jones polynomial estimate is obtained by post-processing the result by the factor in Eq. (11). Notice that, since each of these gates needs to be compiled into a sequence of $L_{i}=\mathcal{O}\left(\log^{c_{2}}(G)\right)$ universal gates from the set $V$, the total circuit size is $L=\mathcal{O}\!\left(G\log^{c_{2}}G\right)=\mathcal{O}\big((|\Sigma|+n)\log^{c_{2}}(|\Sigma|+n)\big)$.

In the example of the trefoil knot reported in Fig. 3, we group single- and two-qubit gates acting on the same target qubits into a new two-qubit gate $U_{i}$, so that the resulting circuit $U$ is composed of only $G=9$ such $U_{i}$ acting on distinct qubit pairs.