Fast verification of fault tolerant magic state preparation circuits with the CNOT + T gate set

Introduction

Fault-tolerant quantum computation is now an experimental reality with much recent work demonstrating some or all of the required ingredients to perform a universal computation, such as state preparation [1], single qubit gates [1, 2] and entangling gates [3]. However, since the Eastin-Knill theorem [4] states that it is not possible to carry out a universal gate set transversally, some of these operations must consume considerable quantum resources, which leads to a performance bottleneck in fault tolerant quantum computing.

One popular approach is to choose a code in which the Clifford group can be implemented transversally, and carry out a single non-Clifford gate such as T or CCZ by injecting a magic state into the code, which is sufficient for universality [5]. When magic states [5] were introduced, it was proposed to prepare magic states fault tolerantly by means of magic state distillation, in which multiple noisy copies of magic states are prepared and then ‘distilled’ to produce a single magic state with a lower error rate. Since the introduction of magic state distillation [5], many improvements and alternative methods have been proposed [6, 7, 8], leading to lower estimates of the resources that these methods consume, in some cases giving improvements in the resources consumed of orders of magnitude. These results have greatly improved the outlook of fault tolerant quantum computing.

However, these methods mostly concern magic state preparation in the asymptotic regime, often limited to the surface code. In the near term, one is often interested in fixed codes of small size, such as small instances of the colour code [9], as these have high rates compared to surface codes of similar size, allowing for more interesting demonstrations to be performed with them on contemporary quantum computers [1, 3]. With such codes, the prospect exists of preparing the magic states directly without distillation, which is to say without the requirement of preparing several magic states at once. Preparing magic states this way naturally reduces the spatial overhead associated with magic state preparation [10, 11, 12] and, in some cases, makes the number of gates required comparable to that of syndrome extraction [13].

It is challenging to construct circuits that prepare magic states fault tolerantly, owing to the fact that within non-Clifford circuits, Pauli faults do not necessarily result in output errors that are also Pauli, which leads to a more complicated analysis than is the case for Clifford circuits [14].

One popular approach is to ‘Cliffordise’ the circuit, which is to say replacing T gates with S gates [12, 14, 15]. While this can lead to a great reduction in computational effort when simulating, it can potentially lead us to miss problematic faults and come to the wrong conclusion about a circuit’s fault tolerance, as we discuss in Appendix A.

Additionally, it is hard to reason analytically and computationally about the fault tolerance of circuits beyond distance three owing to the potentially very large number of combinations of faults in a circuit, each of which must be checked so that it does not cause an error that violates the requirements of fault tolerance.

In this work, we aim to make this task easier by noting that the CNOT + T gate set is capable of preparing a wide range of magic states fault-tolerantly, and showing that it is easier to verify the fault tolerance of circuits that are restricted to the CNOT + T gate set than those that are Clifford + T. We develop a library called Comb that verifies the fault tolerance of a magic state preparation circuit by iterating over the set of combinations of faults that appear in such circuits as efficiently as possible. We illustrate the workflow of checking the fault-tolerance of a CNOT + T circuit in Figure 1. By verifying the fault tolerance of a magic state preparation circuit quickly in software, we aim to facilitate the iterative design of such circuits, which is important as few (if any) analytic or software tools satisfying a similar function exist for any gate set that is not Clifford.

In Section II, we will define some mathematical preliminaries. In § 1.2, we review the properties of the CNOT + T gate set with a view to preparing magic states fault-tolerantly, and discuss how faults propagate through circuits composed of this gate set. In Section IV, we discuss in detail how to computationally verify the fault tolerance of CNOT + T circuits. In § 1.4, we demonstrate the performance of Comb, a library that implements the results of Section IV, on a test suite of circuits, and show that it is feasible to verify the fault tolerance of magic state preparation circuits of reasonable size, including a novel circuit that was discovered using Comb which prepares a magic state fault-tolerantly to distance five. Additionally, we simulate some of the test suite circuits with circuit-level noise and check that these simulations are in agreement with Comb.

Preliminaries

We will consider the Pauli group that contains the Pauli matrices $X,Y,Z$ together with phase,

and the Pauli group on $n$ qubits $\mathcal{P}_n$ to be the $n$-fold product

of the Pauli group. Subsequently, we will often use the term ‘Pauli’ to refer to an element of this Pauli group $\mathcal{P}_n$ on $n$ qubits, and not just elements of the single qubit group $\mathcal{P}$. We will refer to Paulis as phaseless or lacking phase if they do not possess a leading non-trivial phase when written out as products of $X$, $Y$, $Z$ and phase. For example, the Pauli $XZ$ is phaseless while the Pauli $iXY$ is not. We also define Paulis of $X$ type, or $X$ Paulis for short, to be those that are products of the identity or the Pauli $X$, together with phase,

We define Paulis of $Z$ type similarly.

Since $Y=iXZ$, we note that it is always possible to write a Pauli $P$ as a product of a Pauli of $X$ type $P^X$ and a Pauli of $Z$ type $P^Z$,

a decomposition that we will sometimes make use of. It is possible to assume that either $P^X$ or $P^Z$ is phaseless, or reverse their order in the product, without loss of generality. We will quotient out phase in defining the Clifford group $\mathcal{C}_n$,

so that in particular it is finite, and elements of it can be stored in memory on a computer. We refer the interested reader to [16] for a treatment of how to store and manipulate Paulis and Cliffords computationally.

We will also define the $S=\operatorname{diag}(1,i)$ and $T=\operatorname{diag}(1,e^{i\pi/4})$ gates to avoid confusion with definitions in other works. We note that $S=T^2$ and $T^\dagger=T^7$, meaning that $S$, $T$ and their inverses are all generated by $T$.

We will also make use of the CNOT gate, which is

adding subscripts when appropriate to make clear what the indices of the control and target qubits are. When we refer to measurements, they will be assumed to be Pauli measurements in either the $X$ or the $Z$ basis, making clear which one is intended when necessary.

The CNOT + T gate set

Applicability

The CNOT + T gate set is already widely used to prepare magic states fault-tolerantly though this fact is rarely noted explicitly. Example are the flagged fault-tolerant construction of [10], magic state cultivation [12] and distillation [17]. It also has been studied in the context of reducing the T count of a logical circuit [18], and polynomial-time simulation [19].

An illustration of fault tolerantly preparing a magic states in the CNOT + T gate set is shown schematically in Figure 2, based on the flagged fault-tolerant construction of [10]. The sketched circuit prepares a magic state $\overline{T}|+\rangle$ in an $[[n,1,d]]$ code. In this circuit, as in all CNOT + T circuits that we consider in this paper, qubits are initialised in either the $X$ or $Z$ eigenstates $|0\rangle$ and $|+\rangle$ and any measurements that are postselected on are made in the $X$ or $Z$ bases.

DCP operators

We now state and prove some results concerning the propagation of faults through magic state preparation circuits, which will be useful in verifying the fault tolerance of such circuits.

The CNOT + T gate set preserves a set of operators under conjugation. We define a Diagonal Clifford and Pauli (DCP) operator $C_{\mathrm{DCP}}$ to be the product of a diagonal Clifford $C_{\mathrm{D}}$ and a Pauli $P$,

We will assume, without loss of generality, that the Pauli $P$ is of $X$ type, as any $Z$ Pauli can be absorbed into the definition of the diagonal Clifford $C_{\mathrm{D}}$. We can write the diagonal Clifford $C_{\mathrm{D}}$ as the exponential of a linear combination of Paulis of $Z$ type $D_i$, which are diagonal, so that

Expanding the exponential would yield an expression for the DCP error as a sum of Pauli strings.

Under conjugation by CNOT or T, any DCP operator is mapped to another DCP operator. Using Equation (8) we see that conjugating $C_{\mathrm{DCP}}$ by $CNOT_{i,j}$ to obtain $C'_{\mathrm{DCP}}$ gives

As the CNOT gate preserves the X or Z type of a Pauli under conjugation, the terms $CNOT_{i,j}D_lCNOT_{i,j}$ that appear in the sum in the exponential are also Paulis of Z type. Therefore, the resulting operator $C'_{\mathrm{DCP}}$ is also DCP.

Now, conjugating $C_{\mathrm{DCP}}$ by a T gate $T_i$ to obtain $C'_{\mathrm{DCP}}$ gives that

where we have used the fact that diagonal operators commute in order to push the $T_i$ gate through the exponential. If the Pauli $P$ does not have support on the $i$-th qubit, then it will commute with $T_i$ and so $C_{\mathrm{DCP}}=C'_{\mathrm{DCP}}$, meaning that $C'_{\mathrm{DCP}}$ is also DCP. If the Pauli $P$ does have support on the $i$-th qubit, then

which again means that $C'_{\mathrm{DCP}}$ is DCP, as $e^{-\frac{i\pi}{4}}S_i$ is a diagonal Clifford. Since CNOT and T gates take DCP operators to DCP operators under conjugation, so will any product of them.

The DCP operators of fixed size $n$ form a group, and they contain the $n$-qubit Paulis $P_n$ as a subgroup. They play the same role with respect to the CNOT + T gate set that the Pauli group does for the Clifford gate set, as both are preserved under conjugation by the respective gate set. Since Paulis are DCP, and CNOT + T preserves DCP operators, it follows that Pauli faults in a CNOT + T circuit will propagate to output errors which are DCP, which hints at the practical importance of the DCP group for fault tolerant analysis.

However, for the purposes of verifying fault tolerance, it is not the group-theoretic properties of the DCP group that are most relevant, but rather the representational properties of the group. As the DCP group lies inside the Clifford group, it can be represented up to phase in matrix form using $O(n^2)$ bits [16], and multiplication of DCP operators is roughly equivalent to multiplication of these matrices.

This is not much less efficient than the Pauli group, for which Paulis are often represented up to phase with $2n$ bits, and multiplication of Paulis is equivalent to addition mod 2, which is of complexity $O(n)$. The efficient storage and manipulation of DCP operators will also prove useful for fault tolerant analysis.

The discussion so far neglects the role of phase. If one wanted to take into account the phases that result from multiplying Cliffords, this representation would not suffice, and the analysis would be more complicated [20]. However, when analysing the fault tolerance of a circuit, the phase of an output error can be neglected altogether as it cannot change the weight of the error or whether it is detected, so we will settle for representing DCP operators in the manner discussed.

Finally, we note that DCP operators are an example of Pauli-Square-Root Clifford (PSC) operators, also introduced in the context of magic state preparation [14]. We note that DCP operators naturally arise by considering Paulis propagating within CNOT + T circuits, whereas the more general class of PSC operators results from propagating Paulis through Clifford + T circuits that satisfy some technical conditions; for example, the only controlled-$U$ gates that can be performed must have that $U^2=1$ [14]. We posit that our working within the CNOT + T gate set precludes the need to explicitly derive and work around constraints on the circuits we must consider, while retaining sufficient expressive power to prepare magic states for many codes.

Verifying fault tolerance

Definition of fault tolerance

We refer to fault-tolerance in the sense given by Gottesman [21]. We decompose a circuit that prepares some state in terms of its circuit components. Any of these circuit components, which include bare wires, measurements, state preparation, and unitaries may fail, in which case a fault may occur immediately after the circuit component on the qubits where it acts. A combination of faults will cause some corresponding error at the output.

A CNOT + T state preparation circuit is fault-tolerant to distance $d$ if it obeys the condition:

(*) any combination of $s \leq t=\left\lfloor\frac{d}{2}\right\rfloor$ faults produces an output error that is either detected by a pattern of unexpected ancilla measurements or is a sum of Pauli strings each of weight at most $s$.

The weight of the Pauli strings is taken modulo multiplication by the stabiliser group of the code that the state belongs to. This is because multiplying an error by a stabiliser gives an equivalent error which has precisely the same effect on the state. Consequently, an error that appears to be of high weight may be equivalent to one of weight less than or equal to $s$ upon multiplication by a stabiliser, and such errors are not problematic from the perspective of fault tolerance.

It is common to refer to combinations of faults, or the corresponding output errors that they give rise to, as ‘good’ if they satisfy this condition and ‘bad’ otherwise.

Only if all fault combinations are good can a circuit be fault-tolerant, and the existence of a single bad fault combination is sufficient to break fault tolerance.

The condition ($\star$) extends the typical definition of fault-tolerance for Pauli errors to DCP errors, whose weight is less straightforward to define. It ensures that good output errors are still correctable, since the individual components of the error are correctable and hence by linearity the total output error is correctable.

In Appendix B we show that this definition is equivalent to the more traditional criteria that appear in e.g. [21].

Without loss of generality, we will henceforth take all faults to be Pauli, a common assumption which has the advantage that the set of possible faults associated to each circuit component becomes finite. These Pauli faults will propagate to DCP output errors.

It is apparent from the definition of fault tolerance that the number of combinations of faults that must be considered may be very large. Assuming that there are $M$ circuit components, and that each one has at least $f$ faults associated with it, then for each $s$ one must check that at least $\binom{M}{s}f^{s} = \Theta(M^{s}f^{s})$ combinations of $s$ faults from distinct locations are good. This quantity grows very quickly with $M$ even if $f$ and $s$ are not too large.

Algorithm

We propose to verify fault tolerance simply by checking the definition, which is to say, given a circuit, we loop over all possible combinations of $s \le t$ faults, and verify that each error is either detected by a measurement or its weight is lower than or equal to $s$.

Although the number of fault combinations that must be considered is vast, there are many opportunities to save on computational time and space that we will take advantage of.

Storing output errors

First of all, we will address the issue of computing the output error that results from each combination of faults. Given $s$ faults $P_1,\ldots,P_s$, we write $U_1,\ldots,U_s$ for the unitaries that result from the part of the circuit subsequent to each of the $s$ faults until any postselections on ancilla qubits, which are all postponed to the end of the circuit. The error $E_i$ at the output immediately prior to any postselected measurements that results from considering the fault $P_i$ individually is then

The output error $E$ that results from the combination of the faults is equal to the product of the output errors $E_i$ that result from each fault. There is a subtlety that does not arise with Clifford circuits in that the time ordering of the faults becomes relevant. This is due to the fact that output errors are not necessarily Pauli and hence do not necessarily commute or anticommute with each other. Labelling the faults $P_1,\ldots,P_s$ so that when $i<j$, the circuit component corresponding to $P_i$ does not occur later than that of $P_j$, the combined output error $E$ is

Since $P_i$ is Pauli, from the results of Section III B, the output error $E_i$ is DCP. This means that it can be evaluated efficiently. Furthermore, each output error $E_i$ that results from each fault need only be computed exactly once, and then subsequently saved for future reference. This means that computing each output error combination $E$ only involves multiplying together DCP operators that have already been computed, which, as discussed in Section III B, is efficient in time and space.

We note that if we were to work in a more general gate set such as Clifford + T, then each Pauli fault $P_i$ would propagate to an output error $E_i$ that would not necessarily be Pauli or DCP, in which case storing even one output error would take space exponential in the number of qubits, let alone the cost of multiplying such operators, which illustrates one computational advantage of limiting oneself to the CNOT + T gate set.

Postselecting an output error combination

In order to check fault tolerance, we must verify whether an output error combination $E$ is good or bad. There are two parts to checking the goodness of a fault combination. We must determine whether the output error is detected, and if it is not, then we must check whether the weight of the output error is below $s$ or not.

The first of these tasks, determining whether an output error is detected, is considerably more complicated than for Clifford circuits. Given a Clifford circuit, a fault $P$ results in an output error $E$ which is also Pauli. Such errors will either flip or leave invariant any ancilla measurements, depending on whether they commute or anticommute with the measurement operator. Consequently, when one postselects on the results of some ancilla measurements, a Pauli error is either always detected or always undetected, and its effect on the remaining qubits is known.

In a CNOT + T circuit, however, the errors that result from faults are DCPs, not Paulis. As a result, they do not have a deterministic effect on measurement outcomes, as DCPs do not necessarily either commute or anticommute with the measurement operators. Postselecting on ancilla measurements in the presence of a DCP error leads to a conditional error on the data qubits that may not be the same as the error prior to postselection.

We take advantage of the structure of a DCP operator to break the problem into two cases. A DCP operator $C_{\mathrm{DCP}} = C_{\mathrm{D}}P$ can be written as the product of a diagonal Clifford $C_{\mathrm{D}}$ and a Pauli $P$ of $X$ type. Consider first any ancilla measurements that are made in the $Z$ basis. Since these measurement operators will commute with any diagonal operator, we have that a DCP error will flip the measurement deterministically depending on whether $P$ commutes or anticommutes with the measurement operator, just as is the case for Pauli errors, meaning that they can be dealt with the same way. Provided that $P$ is not detected, the conditional $X$ error $E'_X$ on the data qubits that results from postselecting on the $Z$ basis measurements will be $P$ restricted to those qubits.

We then consider postselecting on ancilla measurements that are made in the $X$ basis. The Pauli $P$ will always commute with the measurement operators since it is of $X$ type, which means that it will have no effect on the measurement outcome and can be ignored when considering postselection in the $X$ basis. As a result, we can assume for the rest of this analysis of the effect of postselection that $C_{\mathrm{DCP}} = C_{\mathrm{D}}$ is strictly diagonal Clifford, with no Pauli $X$ component.

In Appendix C, we show that, given an error $E = C_{\mathrm{D}}$ which is diagonal Clifford, the conditional error $E'$ that results from postselecting on a Pauli measurement made in the $X$ basis is either another diagonal operator or vanishes (in which case the error is always detected), and can be computed efficiently. We can apply this result repeatedly to then postselect on the remaining $X$ basis ancilla measurements, and obtain the conditional diagonal error $E'_Z$ that results on the data qubits due to the $X$ basis measurements. The resulting full conditional error $E' = E'_X E'_Z$ is the product of the two conditional errors that result from considering the $X$ and $Z$ basis measurements separately, and is also DCP.

Checking the weight of a conditional output error

Now that we have the conditional output error $E'$ on the data qubits, provided that it does not vanish we must verify whether the weight of the error $E'$ is below the number of faults $s$ or not in order to decide whether the fault combination is good or not. This is made complicated by the fact that $E'$ is not necessarily Pauli. We make use of the fact that, if two errors $A$ and $B$ are correctable, then their coherent sum $A+B$ must be correctable also [21]. We show also in Appendix C that a diagonal Clifford $C$ may be expressed as a linear combination of $2^r$ Paulis of $Z$ type that belong to the coset of some $r$-dimensional vector space $V = \langle P^Z_1,\ldots,P^Z_r\rangle$ whose generators $P^Z_1,\ldots,P^Z_r$ can be computed efficiently. This allows one to express the conditional output error $E' = C_{\mathrm{D}}P$ as a linear combination of Paulis simply by multiplying that of the diagonal Clifford $C_{\mathrm{D}}$ on the right by $P$.

We then verify that the conditional output error $E'$ is good by iterating over this sum and verifying that each element is low weight. We make use of lookup tables to compute the weight of each element in constant time, which we discuss in detail in Appendix § 1.11. Although each element can be checked in constant time, there are $2^r$ elements of the vector space, which is to say this operation takes time exponential in $r$.

Although it is unfortunate to finally encounter an overhead that is not polynomial in some parameter of the problem, the situation is not as bad as it appears.

First of all, there exists the opportunity to iterate over the elements of the vector space in parallel, which reduces the overhead in time. Secondly, the parameter $r$ is not necessarily that large. Without loss of generality, we are able to place the generators of the vector space in canonical form, as detailed in Appendix § 1.11², which limits the size of $r$. For the Golay code, a code which has a higher number of qubits and distance than the other codes considered in this paper, $r$ cannot exceed 12, as also shown in Appendix § 1.11², which is to say there is a factor of $2^{12}=4096$ in the worst case, which is not so prohibitively large as to make computation infeasible in realistic circumstances.

There is one final issue in that, as alluded to in Appendix B, in contrast to the preparation of Pauli eigenstates, magic states are usually stabilised by a set of Clifford operators that are not also Paulis. For example, the magic state $T|+\rangle$ is stabilised by the operator $e^{-\frac{i\pi}{4}}SX$. If the error combination $E' = e^{-\frac{i\pi}{4}}SX$ were to result in a circuit that purports to prepare this magic state $|T\rangle$, then the procedure we have so far outlined would incorrectly conclude that the error combination $E'$ is bad as the $X$ part of it is logical, whereas in reality the error combination $E'$ is good because as it stabilises the magic state it has no effect on it and is thus indistinguishable from no error having taken place.

We resolve this issue by considering the product of the error combination with all possible Cliffords $C_i$ that stabilise the magic state,

and we consider the error combination $E'$ to be good if $E''_i$ satisfies the procedure we have outlined for at least one choice of $i$.

In the case of preparing the usual magic state $T|+\rangle$ on a code with a single logical qubit, then the only Cliffords to consider are the identity and the operator $e^{-\frac{i\pi}{4}}SX$. However, as the number of logical qubits $k$ increases, there will be more Cliffords to consider, which may lead to an unfavourable scaling with $k$.

Demonstration

We use these theoretical insights to design and implement Comb [22], an open-source library that iterates over a set of fault combinations arising from some CNOT + T circuit that purports to prepare a magic state fault-tolerantly and verifies that each of the fault combinations is good. It is designed to be fast so as to aid in the iterative design of fault tolerant magic state preparation circuits of distance greater than three. It takes as input a set of fault locations, a set of DCP operators that correspond to the output error that results from propagating each Pauli fault through the circuit, as well as a description of the magic state that the circuit prepares, that includes lookup tables to verify that the weight of any output error is small.

Name	Code	No. of 2Q gates	No. of faults	$t_1$	No. of combinations of two faults	$t_2$	$t_f$
Steane	$[7, 1, 3]$	38	788	0.006 s	-	-	13 s
$d = 5$, unoptimised	$[17, 1, 5]$	298	5040	0.747 s	12 662 112	8m 3s	11m 17s
$d = 5$, optimised	$[17, 1, 5]$	196	3408	0.486 s	5 781 600	3m 42s	4m 55s
$d = 3$	$[17, 1, 5]$	161	2712	0.260 s	(100 292)	4.05 s	32.9 s
Non-fault tolerant	$[17, 1, 5]$	24	(18)	0.003 s	-	-	12 s

Comb is intended solely to perform the iteration over combinations of faults very quickly, and verify that they are either detectable or low weight. Consequently, in order to use it for a particular circuit, one must choose an appropriate set of fault locations and faults and propagate them through the circuit to obtain the DCP errors, since Comb itself cannot propagate faults through a circuit, and then supply them to Comb. One must also calculate the lookup tables using the algorithm outlined in Appendix § 1.11.

We demonstrate the performance of Comb in Table I when it comes to verifying the fault tolerance of several magic state preparation circuits. All the circuits tested have the general structure shown in Figure 2, based on the construction of [10, 11]. For the $[17, 1, 5]$ code we consider several versions of the construction. One that is non-fault-tolerant and includes no logical stabiliser measurements or error detection cycles. A second that has circuit fault distance $d = 3$ meaning it can handle a single fault but not two, this is achieved by appending one cycle of the logical stabiliser measurement and error detection. Finally, we consider two circuits that achieve the same circuit fault distance as the distance of the code $d = 5$ by including two cycles of logical stabiliser measurements and the error detection. A Dell Latitude 7430 laptop with a Intel(R) Core(TM) i7-1265U CPU and 16GB of RAM was used. However, as Comb is single-threaded at the time of writing the number of cores is not expected to have a bearing on the fraction of the running time spent in Comb, though the faults themselves were propagated in parallel to calculate the output errors.

For each circuit, the CNOT and T gates were assumed to be noisy two- and single-qubit gates respectively that are immediately followed by a fault location that may spawn faults equal to any Pauli on the support of the qubits involved in the gates. Comb was used to verify the fault tolerance criterion that each possible combination of $s$ faults (for appropriate $s \le\left\lfloor\frac{d}{2} \right\rfloor$ chosen according to the distance of the code), leads to either a detectable error or an error of weight less than or equal to $s$ at the output.

One important use case we envision for Comb is to help design fault-tolerant magic state preparation circuits with lower resource requirements, we give a simple illustrative example of that using the $[17, 1, 5]$ colour code. In Table I we show the performance of Comb on two $d = 5$ variants of the magic state preparation circuit for the $[17, 1, 5]$ code. The smaller of these circuits was obtained from the first using a simple greedy iterative optimization process with Comb certifying at each step that the reduced circuit remained fault-tolerant. The final circuit, along with additional details of the optimisation, are given in Appendix E.

As an independent check of Comb’s correctness we perform noisy circuit simulations of all the $[17, 1, 5]$ circuits tested with Comb and compute their logical error rates, plotted as a function of physical error rate $p$ in Figure 3. In each case we observe scaling of the logical error rates with the physical error rate that is consistent with the findings of Comb. We additionally plot the post-selection acceptance rates for each circuit. These simulations were carried out using SymFT [23]. Additional details of the simulations are given in Appendix F. It is interesting to note that in these simulations we see better logical error rate for the unoptimised distance 5 circuit than the optimised circuit. We attribute this to the simplicity of the error model used in the simulations, which does not include memory errors such as from coherent dephasing or qubit leakage. In experiments we would expect to see advantages from shorter depth resource state preparation circuits that are not captured in these simulations, including increased robustness to memory errors as well as decreased idling times for other logical qubits in the computation. The shorter circuit depths combined with improved acceptance rates has the potential to improve the rate of producing magic states while decreasing the number of ancilla qubits needed to manufacture them.

Discussions and outlook

In this paper, we make the case to explicitly restrict magic state preparation circuits to the CNOT + T gate set, and present Comb, a library that implements our algorithm insights and is designed to verify the fault tolerance of such circuits.

One future direction is to use Comb to design new fault tolerant circuits, which is an area of active research. The distance-five magic state preparation circuit that is presented in Figure 6 achieves a reduction of over 30% in the number of CNOT gates over the circuit it was derived from, and the fact that it was obtained simply by greedily deleting ancillas from another fault tolerant state preparation circuit gives us hope that there exist many fault tolerant magic state preparation circuits with much lower resource consumption than currently known distillation- or cultivation-based circuits, that may be verified using Comb.

!Fault tolerant preparation circuit for $T\lvert+\rangle$ in the [17, 1, 5] colour code bbox=0.10,0.12,0.90,0.95

Work is underway to analyse pytket [24] circuits with Comb, which will greatly reduce the barrier to entry of using Comb for QEC researchers.

We anticipate that an autonomous agent, perhaps powered by AI, may be able to propose new magic state preparation circuits and then use Comb to check them. This suggests the possibility of closed-loop generation of magic state preparation circuits, which may lead to novel circuits that a human designer would struggle to propose in reasonable time.

There have been a plethora of recent results [25, 26, 27, 28] concerning fault tolerant rewrites for the ZX calculus, which allow one to change a fault-tolerant circuit using local rewrites that do not destroy its fault tolerance properties. Such techniques, relying as they do on local operations, are inherently scalable and promise to realise a practical way to reduce the resources that a fault tolerant quantum computation consumes at larger distances. To our knowledge, such results have been done exclusively on Clifford circuits to date. Our results concerning the practicality of the CNOT + T gate set for fault tolerant magic state preparation motivate the investigation of fault tolerant rewrites for ZX diagrams that correspond to CNOT + T circuits.

There exist many potential performance improvements to Comb that have not yet been implemented, such as use of parallelism, with both multiple cores and SIMD. We expect such improvements to have a great impact on performance, as verifying fault tolerance is inherently a parallel problem, as different combinations of faults can be divided into groups and checked independently in parallel.

Such performance improvements will not only aid in the exploration of distance-three and distance-five magic state preparation circuits, but will plausibly, combined with the use of HPC, make possible the exploration of distance-seven magic state preparation circuits, which is an exciting prospect.

Acknowledgments

The authors are employed by Quantinuum. We would like to thank Asmae Benhamou, Serban Cercelescu, Ben Criger and Pablo Andrés-Martínez for productive discussions, Alec Edgington for invaluable help with C++, Marek Kochanowicz for useful insights regarding program optimisation, and Mateusz Kupper for assistance with the noisy simulations.

References

Counterexample to checking fault-tolerance by replacing $T$ gates with $S$ and $Z$

A trick that is sometimes used when designing and testing magic state preparation circuits is to replace the $T$ gates with $S$ or $Z$ gates [12, 15]. For some circuits (though not generically) this will replace the preparation of $T|+\rangle$ with $S|+\rangle$ or $Z|+\rangle$, without changing the circuit structure. Checking whether these CNOT + $S$ or CNOT + $Z$ circuits are fault-tolerant is much easier because they are Clifford circuits, though they are not preparing the desired magic state.

In this appendix we consider circuits $C(Q)$ with the property that there are a set of $Q$ type gates in the circuit such that setting $Q=Z$, $S$ or $T$ prepares the logical state $Q|+\rangle$ and we ask what, if any, links we can draw between the fault-tolerance of $C(S)$, $C(Z)$ and $C(T)$. An example of such a circuit is shown in Figure 4(a), for the $[[7,1,3]]$ Steane code.

$C(S)$ and $C(Z)$ FT $\not\Rightarrow C(T)$ FT

We might hope that by finding a way to make $C(S)$ and $C(Z)$ fault-tolerant we will discover circuits that will also be fault-tolerant for $C(T)$.

Here we present a circuit $C(Q)$ where $C(S)$ and $C(Z)$ are fault-tolerant preparations of $S\lvert+\rangle$ and $Z\lvert+\rangle$ while $C(T)$ is a non-fault-tolerant preparation of $T\lvert+\rangle$. This demonstrates that $C(S)$ and $C(Z)$ being fault-tolerant does not in general guarantee that $C(T)$ is also. This agrees with [12], which found in small statevector simulations that the logical error rate for preparing $T\lvert+\rangle$ in a magic state cultivation scheme was about twice as high as $S\lvert+\rangle$ or $Z\lvert+\rangle$.

The circuit we consider is shown in Figure 4(a), it is a non-fault-tolerant state preparation circuit for the $[[7,1,3]]$ Steane code. Setting the highlighted gates to $Q=Z$, $S$ or $T$ prepares the logical state $Q\lvert+\rangle$. It is composed of two parts: a non-fault-tolerant encoding circuit, followed by a flagged fault-tolerant measurement of the logical stabiliser [10]. The logical stabiliser measurement catches logical errors. Nevertheless, it is non-fault-tolerant because there are single faults that propagate to errors on the output with weight 2 or higher.

We can identify minimal size sets of stabilisers that make the circuits fault-tolerant if we measure them at the end of the circuits. The $Z\lvert+\rangle$ circuit can be made fault-tolerant by measuring a single stabilier element:

$IXXIIXX$.

For the $S\lvert+\rangle$ circuit there are no single stabilisers that catch all bad errors, but there are pairs of stabilisers that will make the circuit fault-tolerant:

$YYIYIIY, YIYYIYI$

$YYIYIIY, IYYIIYY$

$ZZIZIIZ, IYYIIYY$

$YIYYIYI, IYYIIYY$

$ZIZZIZI, IYYIIYY$

$IYYIIYY, YYIIYYI$

$IYYIIYY, ZZIIZZI$

$IYYIIYY, YIYIYIY$

$IYYIIYY, ZIZIZIZ$

$YYIIYYI, YIYIYIY$.

Finally, there are triples of stabilisers that can be appended to the $T\lvert+\rangle$ circuit to make it fault-tolerant:

$ZZIZIIZ, IIIZZZZ, IXXIIXX$

$ZZIZIIZ, IYYIIYY, XXIIXXI$

$ZZIZIIZ, ZZIIZZI, IXXIIXX$

$IIIZZZZ, ZZIIZZI, IXXIIXX$.

These stabiliser elements can be measured fault-tolerantly at the end of the circuits using flag gadgets [29].

Appending three fault-tolerant measurement gadgets to the end of the circuits to measure the $IXXIIXX$, $YYIYIIY$ & $YIYYIYI$ stabilisers yields a circuit $C(Q)$ such that $C(S)$ and $C(Z)$ are fault-tolerant but $C(T)$ is not.

Although there are triples of stabilisers that we can measure to make the circuit in Figure 4(a) a fault-tolerant preparation of $T\lvert+\rangle$ this choice made by considering the $S$ and $Z$ variants is not sufficient.

An intuitive understanding of this failure can be gained by considering the Pauli propagation of a problematic fault, Figure 4(b). A $YI$ error after the CNOT gate highlighted (bold, red) in Figure 4(a) propagates through the subsequent CNOT gates to become $YZZIIII$ before the first layer of $T$ gates. The $T$ gate on qubit 1 acts as a linear combination of an $S$ gate and a $Z$ gate, which we draw as two paths leading to $XZZIIII$ and $YZZIIII$. To see this, consider the conjugation of a Pauli $P$ that anticommutes with $Z$ through a $T$ gate

For our purposes here we can neglect the phases on these terms. Through the CNOTs of the logical stabiliser measurement the second term gains a $Z$ on the first ancilla to become $YZZIIII Z$, this term is a product of a $Z$ logical and a single qubit $X$ error and so the logical stabiliser measurement has been triggered. This term splits again at the second layer of $T^\dagger$ gates but this does not remove the support on the ancilla and so both are detected when we measure the ancillas. The $XZZIIII$ term, in contrast, does not trigger the logical stabiliser measurement and splits at the $T^\dagger$ layer into two terms $YZZIIII$ and $XZZIIII$. The first of these anti-commutes with our choice of final stabiliser measurements and is detected, but the second commutes with the stabilisers and so is undetected at the output of the circuit, becoming a high weight bad error. By inspecting what paths have given rise to detected and undetected errors it is clear what has happened: the $(S,S^\dagger)$ branch is detected, we have guaranteed this when considering $C(S)$. The problematic error has arisen from the $(S,Z)$ branch, which is absent when considering $C(S)$ and $C(Z)$.

Does $C(S)$ and $C(Z)$ non-FT $\Rightarrow C(T)$ non-FT?

If $C(S)$ and $C(Z)$ being non-fault-tolerant guaranteed that $C(T)$ also cannot be fault-tolerant this would still be a useful property as it would offer a fast rejection method for candidate circuits. Using the fact that $C(S)$ and $C(Z)$ are Clifford circuits, their fault-tolerance can be checked with standard Pauli propagation and weight minimisation methods. We could quickly test $C(S)$ and $C(Z)$, and only proceed to the more expensive testing of $C(T)$ if they both pass.

We leave this as an open problem. We are not able to provide a counter example, but we cannot prove that it must always be true.

Fault tolerance condition

Here we show that the condition given for the fault-tolerance of CNOT + T state preparation circuits in section IV A is equivalent to those commonly given in the literature [10, 21]. We specifically compare to the definition given by Gottesman in [21].

For completeness, the condition we use is: a CNOT + T state preparation circuit is fault-tolerant to distance $d$ if it obeys:

$(\star)$ any combination of $s \leq t = \left\lfloor\frac{d}{2} \right\rfloor$ faults produces an output error that is either detected by a pattern of unexpected ancilla measurements or is a sum of Pauli strings each of weight at most $s$.

The weight of a Pauli string must be taken modulo the stabilisers. This is achieved with a weight lookup table, as discussed in appendix § 1.11.

The fault-tolerance conditions for state preparation circuits presented in [21] make use of $r$-filters and ideal decoders. This machinery is more general than we needed so we avoided it in our definition, opting instead for a simpler (but equivalent) condition based on Pauli weights.

An $r$-filter for a QEC code projects onto the subspace of errors spanned by states $E|\psi\rangle$ where $|\psi\rangle$ is any code-word and $E$ is a Pauli error of weight at most $r$. We use the terminology passing an $r$-filter to mean that the state is unchanged by applying an $r$-filter.

An ideal decoder is a quantum channel that takes an encoded state, noiselessly extracts the syndrome of any errors before correcting them and returning the decoded logical state. If the errors on the encoded state are un-correctable the ideal decoder will return a state with a logical error.

Fault-tolerance for a state-preparation circuit is defined in [21] by the two properties:

1. Preparation Error Propagation: any combination of $s \leq t$ faults passes an $s$-filter,

2. Preparation Correctness: the noisy output state will be decoded by an ideal decoder to the logical state that the noiseless circuit would prepare, provided $s \leq t$ faults occur,

where $t = \left\lfloor\frac{d}{2} \right\rfloor$ is the number of correctable errors.

We first show that if $(\star)$ is satisfied then (1) is satisfied. Let us write the output DCP error as

where the $P_i$ are Pauli strings. The weight $E$ is taken modulo the stabilisers. To be explicit the weight is minimised with respect to the full stabiliser group of the code $S$ enlarged with the logical Clifford $C$ that the noiseless output state is the $+1$ eigenstate of, e.g. for the $T|+\rangle$ logical state $C = \overline{T}\overline{X}\overline{T}^{-1} = \overline{S}\overline{X}$. Consider the minimisation under $S$ first, each $P_i$ can be written

where $Q_i$ is the Pauli string with minimum weight and $S_i$ is an element of the stabiliser group. The action of $E on the noiseless encoded output state $|\overline{\Psi}\rangle$ is

where we have used the fact that $|\overline{\Psi}\rangle$ is stabilised by $S_i$. The noiseless state $|\overline{\Psi}\rangle$ is also stabilised by $C$ so we can consider the alternative DCP error

where the sum will generically contain a different number of Paulis. Minimising the weight of each $P'_j$ under the stabilisers produces an alternative set of minimum weight Paulis $\{Q'_j\}$. Condition (⋆) requires that the weight of all the Paulis $\{Q_i\}$, or all the $\{Q'_j\}$, is at most $s$ when $s \leq t$ faults occur. If this condition holds then the output state Eqn. (B3) will pass an $s$-filter, since $|\overline{\Psi}\rangle$ is a codeword, hence condition (1) is satisfied.

Next we show that if (⋆) is satisfied then (2) is satisfied. If the weight of all the $Q_i$ is at most $s \leq t$ then all the $Q_i$ are correctable errors [21]. In this case an ideal decoder is sure to decode the output state Eqn. (B3) to $|\overline{\Psi}\rangle$, the logical state of the noiseless circuit, and so (2) is satisfied.

Diagonal coherent error analysis

A great advantage of DCP errors is that they are able to compactly describe coherent errors which, when expressed as a sum of Paulis, may be very large.

In fault tolerant analysis, one is interested in whether errors are detected by ancilla measurements, and if not, what the residual conditional error on the remaining qubits is. It is known how to do this efficiently for errors that are Pauli. For a general coherent error, one is compelled to expand the coherent error as a sum of Paulis and check each term in the sum in turn. If the sum is very large, then even initially writing down the sum (or storing it memory on a computer) for a single coherent error will be infeasible, let alone considering a very large number of errors.

In this section, we will show that given a coherent error prior to any measurements that is DCP, it is possible to efficiently calculate and describe the conditional error that results from postselecting on the results of ancilla measurements, without having to expand any sums explicitly. This result is of great importance in ensuring that DCP errors can be verified to be good or bad in a reasonable timeframe.

Operator state duality

As it is possible to deal with the Pauli X part of a DCP error in much the same way as a Pauli error, we assume for the rest of this section, without loss of generality that the error $E$ that results immediately prior to any ancilla measurements is strictly diagonal Clifford.

As discussed, we are interested in dealing with the error $E$ without having to express it directly as a sum of Paulis. In order to do so, we will make use of a operator-state duality. We define the error state $|\psi_E\rangle$ to be the state

The error state can be defined not only for the error $E$ that is diagonal Clifford, but for all diagonal operators. The error state contains all information about the error $E$; it is possible to reconstruct the operator $E$ from the state $|\psi_E\rangle$ by observing that the coefficients of the state $|\psi_E\rangle$ in the computational basis form the diagonal of the operator $E$ in that basis.

It will turn out that many operations of interest that one would like to carry out on the coherent error $E$, such as determining the conditional coherent error $E'$ that results from postselecting on an ancilla measurement, are best expressed in terms of operations on the corresponding error state $|\psi_E\rangle$.

Before turning to them, we will first note that since the error $E$ is diagonal Clifford, the error state $|\psi_E\rangle$ is a Pauli stabiliser state, as it is defined as the application of a Clifford $E$ to the state $|+\rangle^{\otimes n}$, which is a Pauli stabiliser state. Furthermore, a choice of stabiliser generators for the error state $|\psi_E\rangle$ can be determined easily from a description of the diagonal Clifford $E$, as they are simply the images of the Pauli operators $X_i$ under conjugation by the error $E$.

Postselection

In order to detect such diagonal coherent errors, one or more ancillas may be measured in the $X$ basis, and the state discarded if one of the measurement results, or a parity of them, is non-trivial. We denote $P_X$ to be the Pauli of $X$ type corresponding to one such measurement.

Writing the error $E$ as a sum of diagonal Paulis,

the error $E'$ that results will be of the form

which is to say that the errors that anticommute with the measurement operator are detected, while those that do not are undetected and survive to form the coherent error conditional on postselection.

Eliminating terms from $E$ in this drastic manner means that $E'$ may not necessarily be Clifford; indeed it may vanish altogether. Despite this, it is certainly diagonal, which means that it possesses an error state $|\psi_{E'}\rangle$, and there exists a simple relation between the error states $|\psi_E\rangle$ and $|\psi_{E'}\rangle$.

Applying Equation (C2) and Equation (C3) to the state $|+\rangle^{\otimes n}$ gives that

and

The terms $P_i|+\rangle^{\otimes n}$ that appear in Equation (C5) are those terms appearing in Equation (C4) for which $P_i$ commutes with $P_X$. These terms are precisely the positive eigenstates of $P_X$, which means that they are preserved under the projector onto the positive eigenspace of $P_X$, while the unwanted terms are eliminated by this projector.

Concretely, we get that

As previously discussed, the error state $|\psi_E\rangle$ is a Pauli stabiliser state. The above equation gives that the state $|\psi_{E'}\rangle$ is equal to that which results from measuring the Pauli operator $P_X$ to be $+1$ on the state $|\psi_E\rangle$. But as the state $|\psi_E\rangle$ is a Pauli stabiliser state, a choice of stabiliser generators of the state $|\psi_{E'}\rangle$ can be evaluated efficiently by updating those of the state $|\psi_E\rangle$[21]. We note also that the resulting state may vanish, in which case the conditional error $E'$ vanishes also. This case corresponds to the situation in which the error $E$ is completely detected by the measurement.

Although the error $E'$ is not necessarily Clifford, its error state $|\psi_{E'}\rangle$ is a Pauli stabiliser state just as the state $|\psi_E\rangle$ is. This means that we can apply the above technique again to $|\psi_{E'}\rangle$ to evaluate the error state that results from postselecting on any number of ancilla measurements, not just one, simply by measuring the corresponding measurement operators in turn on the error state.

In this way, one can describe in terms of its error state the conditional coherent error $E_C$ that results from all postselections made in the $X$ basis. Such a state will necessarily have much fewer terms in its expansion than the original error $E$.

Canonicalisation

The conditional coherent error $E_C$, when expressed as a linear combination of diagonal Paulis, may contain distinct terms which are equal up to a stabiliser. It is important to eliminate such terms in the sum by stabiliser multiplication and cancellation in order to accurately evaluate whether the error is good (low weight) or not.

To illustrate this more clearly, we give a hypothetical example. A circuit that purports to prepare the magic state $T|+\rangle$ in the Steane code might have a coherent error, after postselection, of the form

Both of these terms are logical operators, and so they cannot be equivalent to a weight one error up to multiplication by any $Z$ stabiliser, since the code is distance three. Consequently, they are both undetectable bad errors and a circuit that produced either of these errors in the presence of a single fault would have to be rejected as not fault tolerant.

However, they are both equal after multiplication by the stabiliser generator $IIIZZZZZ$. Consequently, the error $E_C$ vanishes, which necessarily means that it is good, as a vanishing error is one that was already detected by postselecting on an ancilla measurement.

Consequently, when considering coherent errors, one must take care to ensure that errors that are equivalent up to multiplication by a stabiliser have already been cancelled out in the sum, otherwise checking each term in the sum in turn would lead to a spurious conclusion that some good coherent errors are bad.

Surprisingly, it is possible to perform such a cancellation by measuring the error state $|\psi_{E_C}\rangle$ without having to expand the sum at all.

We note first that, following Appendix § 1.11, every Pauli operator has a unique canonical form that results from multiplying by a choice of stabiliser generators. Writing the operator

we note that if every operator $P_i$ is replaced by its canonical form, then all cancellations in the sum will take place automatically, as operators that are equivalent up to a stabiliser will be mapped to the same canonical form.

Every term in the sum $E_C$ can be canonicalised by multiplying on the left by $E_C$ by an appropriate choice of operators. In order to help illustrate this, we will give a worked example. In order to canonicalise operators with respect to the Steane code, operators which have support on the first qubit must be multiplied by the stabiliser generator $S_Z = ZIZIZIZ$, and those that do not are left unchanged. Doing the same for subsequent stabiliser generators of the Steane code ensures that all operators have vanishing support on the first, second and fourth qubits, and are so left in canonical form.

Considering for now only the effect of canonicalising on the first qubit, we must separate the error $E_C$ into a part $E_C^+$ which has trivial support on the first qubit, and a part $E_C^-$ that doesn’t. It will prove simpler to discuss these operators in terms of their error states. The error states according to these operators can be obtained straightforwardly by projecting into the positive and negative eigenspaces of $X_1$ respectively,

The terms in $E_C^+$ have vanishing support on the first qubit and do not need to be changed. However, all the terms in $E_C^-$ have support on the first qubit, and must all be multiplied by $S_Z$ in order to canonicalise them. This can be accomplished by applying the operator $S_Z$ to the error state $|\psi_{E_C^-}\rangle$. Once this is done, both halves of the sum can be recombined. The error $E'_C$ that results is

Substituting (27) and (28) gives that

This expression can be simplified by observing that $S_Z$ anticommutes with $X_1$, as $S_Z$ has support on the first qubit. One gets, after pushing $S_Z$ past $X_1$ and factorising, that

This corresponds to measuring two Pauli operators to be $+1$ on the error state $|\psi_{E_C}\rangle$, and generalises to canonicalising on any qubit in any CSS code, not just the Steane code.

As the state $|\psi_{E_C}\rangle$ is a Pauli stabiliser state, just as before the state $|\psi_{E'_C}\rangle$ that results from measuring Paulis on it is also a Pauli stabiliser state, and a description of it can be obtained efficiently.

Error states to cosets

In this way, one can efficiently determine the final coherent error $E_F$ that results from postselecting on all relevant ancilla measurements and eliminating terms that are equal up to a stabiliser multiplication. The next step is to expand the error $E_F$ as a linear combination of Paulis and ensure that each of them is low weight.

The previous results have shown how to efficiently obtain a description of the error $E_F$ in terms of its error state $|\psi_{E_F}\rangle$. We now turn to the problem of obtaining an expansion of $E_F$ in terms of Paulis from its error state. Although the number of terms in the sum will turn out to be exponential in some parameter, they all belong to the coset of a vector space which can be described and calculated efficiently.

Let us first relate the expansion of $E_F$ to the error state. Writing the error

as a linear combination of distinct Paulis $P_i$ of $Z$ type, and applying both sides to the state $|+\rangle^{\otimes n}$, we get that

Consider the states $P_i|+\rangle^{\otimes n}$. The effect of $P_i$ on the state $|+\rangle^{\otimes n}$ is to flip $|+\rangle$ to $|-\rangle$ on qubits on which $P_i$ has support, and leave it invariant otherwise. Consequently, one can read off $P_i$ from the patterns of $+$ and $-$ on such states. This means that all the information about the sum - in particular the Pauli $P_i$ that make up the sum - can be obtained by expressing the error state $|\psi_{E_F}\rangle$ in the basis of X stabiliser states $B=\{|+\rangle,|-\rangle\}^{\otimes n}$ and seeing which basis terms form the combination.

Since the error state $|\psi_{E_F}\rangle$ is a Pauli stabiliser state, it is possible to compactly characterise the set of diagonal Paulis $P_i$.

We first fix a choice of stabiliser generators $S_1,\ldots,S_n$ for the state $|\psi_{E_F}\rangle$. In order to make contact with the basis $B$, we will fix some X basis state $|C\rangle\in B$ and write

where we have related the error state $|\psi_{E_F}\rangle$ to an X computational basis state $|C\rangle$ by applying the projector onto the error state $|\psi_{E_F}\rangle$ to the state $|C\rangle$. As the X computational basis states are a basis, there must always exist at least one choice of state $|C\rangle$ such that the constant of proportionality is not zero i.e. at least one state $|C\rangle$ must have non-zero overlap with the error state $|\psi_{E_F}\rangle$. We will turn later to the question of how to choose the state $|C\rangle$, and assume for now that it has been chosen appropriately by some means.

Let us now consider the effect of projecting the state $|C\rangle$ onto the positive eigenspace of the stabiliser generators $S_i$. In order to simplify the analysis, we will consider two special cases. If $S_i$ is of Z type and not the identity, then one gets that

where the state $S_i|C\rangle$ is a X computational basis state that is distinct from $|C\rangle$.

If $S_i$ is of X type, then $S_i|C\rangle$ will be equal to $|C\rangle$ up to a sign. Correspondingly, acting on $|C\rangle$ with the projector $\frac{I+S_i}{2}$ will either give zero or the state $|C\rangle$, depending on the sign.

Of course, the stabiliser generators $S_i$ are not known to be of either X or Z type, and are unlikely to be. However, armed with these observations, we can cast the stabiliser generators $S_i$ into a different set $S'_i$ via Gaussian elimination, that will allow us to apply these observations.

Concretely, we perform Gaussian elimination on the Z parts of the stabiliser generators $S_i$. A new set of stabiliser generators $S'_i$ stabilising the same state is produced, such that for some known $r$, for all $i \le r$ the Z parts of the stabiliser generators $S'_i$ are linearly independent, and for all $i > r$, the Z parts of the stabiliser generators $S'_i$ vanish.

Substituting the new set of stabiliser generators $S'_i$ into Equation (C16), we get

First, we apply the projectors corresponding to $i > r$ to the state $|C\rangle$. The stabiliser generators $S'_i$ are all of X type by construction, and so as discussed they either eliminate the state $|C\rangle$ or preserve it. It is always possible to choose the state $|C\rangle$ such that none of the projections eliminate the state by writing $|C\rangle= P_C|+\rangle^{\otimes n}$, where $P_C$ is some Pauli of Z type, and choosing $P_C$ to anticommute or commute with $S'_i$ as appropriate to ensure that the projections do not eliminate the state $|C\rangle$. This can always be done efficiently by solving a linear system of equations over $\mathbb{F}_2$.

Once this has been done, the projections act trivially on the state $|C\rangle$ and we are left with

The operators $S_1'^{a_1}\cdots S_r'^{a_r}$ that appear in the sum are not necessarily of either X or Z type, so when they act on the state $|C\rangle$ they can both flip the sign and take $|C\rangle$ to a different X computational basis state $|C_{a_1,\ldots,a_r}\rangle$, depending on the Z parts of the operator $S_1'^{a_1}\cdots S_r'^{a_r}$. However, as the Z parts of the operators $S'_i$ are all linearly independent by construction, the states $|C_{a_1,\ldots,a_r}\rangle$ are all unique, meaning that there can be no cancellations in the resulting sum. Furthermore, the states $|C_{a_1,\ldots,a_r}\rangle$ can be obtained by acting on the state $|C\rangle= P_C|+\rangle^{\otimes n}$ by all the $2^r$ ways of forming products of the Z parts of the stabiliser generators $S'_i$.

Consequently, one gets that the coherent error $E_F$ can be expressed as a linear combination of $2^r$ Paulis that form the coset of a vector space

where the operator $S_i^{\prime Z}$ is the Z part of the stabiliser generator $S'_i$.

We conclude by remarking that although the number of errors in the coherent sum is $2^r$, it is at least possible to calculate and describe the set efficiently, and that the elements of the coset may be iterated over in constant memory via a Gray code for example.

Weight lookup tables

We make liberal use of lookup tables to fetch the weight of an operator up to multiplications by a stabiliser from some code. In this Appendix, we discuss novel theoretical and technical work related to lookup tables that is key to the performance of Comb and may be of general interest to those working with stabiliser codes.

For simplicity of discussion, we will limit ourselves to CSS codes as we do not consider non-CSS codes elsewhere in this work, though much of the discussion does generalise to non-CSS codes.

Lookup tables

A lookup table is a simple data structure that stores the result of some computation, so that it need not be computed more than once. Given a computable but expensive function $f$ that acts on a finite set $R$ of size $N$, we produce a lookup table for $f$ by first labelling the members $r_i \in R : 0 \le i < N$ and then computing and storing the values $f(r_i)$ in an array of size $N$. Subsequent evaluations of $f$ may be made by looking up the appropriate entry in the array, which is an operation that takes constant time.

This dramatic saving in time naturally comes at the cost of a great increase in space. If the set $R$ is too large, it may be too impractical to store it in memory. Furthermore, we must evaluate $f$ on every element of $R$ in order to form the lookup table. This will not be worthwhile unless the number of evaluations of $f$ is likely to exceed $N$. Since by assumption $f$ is expensive to evaluate on individual elements, then evaluating $f$ on the entire set $R$ is likely to be very expensive indeed. However, if $R$ is known, then the lookup table only needs to be computed once, and depending on the details of evaluating $f$, forming the lookup table may not be prohibitively expensive.

We argue that lookup tables are appropriate to compute the weight of operators modulo a stabiliser code of reasonable size. Firstly, determining the weight of an operator $P$ modulo a stabiliser code is a computationally costly problem. Given a code with $s$ stabiliser generators of a code, one way to compute the weight of $P$ is to form all possible products of $P$ with a stabiliser and find the minimum weight over the products. However, as there are $2^s$ different stabilisers, this is not very efficient, and prohibitive for codes of even moderate size. We note also that, for the purposes of checking fault tolerance, the number of output error combinations that must be considered is extremely large, and accordingly we expect to check the weight of a very large number of operators. Any potential saving in the cost of evaluating the weight of an operator must therefore be considered.

We present two results that greatly facilitate the use of tables to look up the weight of an operator; firstly, a trick to compress the size of a lookup table, and secondly an algorithm to compute lookup tables in reasonable time.

Compression trick

Given a $[n,k,d]$ CSS code, we define the CSS weight lookup tables $T_i^X$ and $T_i^Z$,

where $P_X$ and $P_Z$ are Paulis of $X$ and $Z$ type respectively, $h$ is the Hamming weight, $C_X$ and $C_Z$ are the spaces spanned by the $X$ and $Z$ stabilisers of the code $C$ respectively. One can verify whether the weight of a Pauli $P_X$ of $X$ type is of weight lower than or equal to $s$ by looking up the appropriate entry in the table $T_i^X(P_X)$, and similarly for a Pauli $P_Z$ of $Z$ type.

In practice, we will limit ourselves to the case $i\le\left\lfloor\frac{d}{2}\right\rfloor$, since this is sufficient for verifying fault tolerance.

We define the weight lookup tables this way for two reasons: first of all, when we check whether a given output error combination is good, we are not interested so much in the weight of an output error combination as much as verifying that its weight is below some threshold $s$. Consequently, the weight lookup tables will contain all the information - one bit per Pauli - that is required and no more, which is optimal from the perspective of storing information.

The Paulis $P_X$ and $P_Z$ can be put into correspondence with bit strings of length $n$, which mean that each of these lookup tables are of size $2^n$. Given a code of reasonable size such as the Golay code, with $n=23$ and $d=7$, each of the six lookup tables will occupy 1MiB of space. These can easily fit into memory, which illustrates the practicality of lookup tables for this problem.

However, one can compress the lookup tables by a considerable degree by casting the stabiliser generators of a code into a particular form.

To use the Golay code as an example, its stabiliser generators of $Z$ type can be written,

ZIZIIZIIZZZZZIIIIIIIIII
ZZZZIZZIZIIIIZIIIIIIIII
IZZZZIZZIZIIIIZIIIIIIII
IIZZZZIZZIZIIIIZIIIIIII
IIIZZZZIZZIZIIIIZIIIIII
ZIZIZIZZZIIZIIIIIZIIIII
ZZZZIIIZIIZZIIIIIIZIIII
ZZIZZZIIIZZIIIIIIIIZIII
IZZIZZZIIIZZIIIIIIIIZII
ZIIZIIZZZZZIIIIIIIII...

where the tableau restricted to the final eleven qubits is the identity. We note that this can be done for any stabiliser code up to a qubit permutation by performing Gaussian elimination.

Given an operator $P_Z$ on 23 qubits, one can multiply $P_Z$ by these stabiliser generators as appropriate to eliminate any $Z$ operators on the last eleven qubits, thus obtaining another operator $P'_Z$ which acts non-trivially on at most the first twelve qubits. We will refer to $P'_Z$ as the canonical form of $P_Z$.

Since multiplying an operator by a stabiliser preserves its weight modulo the stabiliser code, instead of using $P_Z$ to index into the lookup table we can instead use $P'_Z$. This means that the lookup tables now each occupy $2^{12}$ bits, or 512 bytes, which is a massive saving, and makes the use of lookup tables practical for larger codes.

This saving in space is in fact likely to deliver a saving in time. A lookup table of 512 bytes can fit into the L1 cache of many CPUs, which, assuming that the rest of the program does not occupy too much space, is likely to make memory accesses considerably faster.

Computing a weight lookup table

The compression trick discussed in the last section gives rise to an inductive algorithm for computing the lookup tables $T_i^Z(P_Z)$. Computing the tables $T_i^X(P_X)$ will follow by symmetry.

Considering the first base case $i=0$, we observe that $T_0^Z(P'_Z)$ is 1 if and only if $P'_Z$ is the identity. This is because only stabilisers may have weight zero modulo the code, and all stabilisers have the identity as their canonical form.

We will proceed to generate the tables $T_i^Z(P'_Z)$ inductively. Given tables $T_0^Z,\ldots,T_i^Z$, we will show how to generate table $T_{i+1}^Z$.

Iterate over all elements $P$ such that $T_i^Z(P)=1$ i.e. those that have weight less than or equal to $i$ modulo the code. For each $1\le j\le n$, we form the product $PZ_j$. By construction, $PZ_j$ must have weight less than or equal to $i+1$. We then check to see if $T_i^Z(PZ_j)=0$, noting that this involves computing the canonical form of $PZ_j$. If so, then we know that $PZ_j$ is of weight $i+1$, so we add it to $T_{i+1}^Z$ and continue. Every operator of weight $i+1$ must be of this form for some $P$ and $j$, so we are guaranteed to reach all operators of weight $i+1$ this way.

Consequently, once we have finished the iteration, $T_{i+1}^Z(P)$ will contain all operators of weight exactly $i+1$ modulo the code. We then set $T_{i+1}^Z(P)=1$ for all $P$ in $T_i^Z$ to ensure that $T_i^Z$ contains all operators of weight at most $i+1$, which completes the induction.

Although this algorithm is not polynomial in time or space, it need only be evaluated once per code, and it does not take excessive time to evaluate it for many codes of reasonable size such as the aforementioned Golay code. It is to be hoped that it is practical to compute and use lookup tables for larger codes, though this will not scale indefinitely.

Greedy circuit optimisation of the distance five circuit for $[\![17,1,5]\!]$

In Section V we use Comb to test the fault-tolerance of a set of magic state preparation circuits for the $[\![17,1,5]\!]$ colour code, including two variants that are fault-tolerant to the distance of the code $d=5$. The second $d=5$ circuit is obtained from the first by a circuit optimisation procedure. This optimisation is made possible by the high speed that Comb is able to test candidate circuits. In this appendix we describe that optimisation process in detail. The code definitions we adopt for $[\![17,1,5]\!]$ are shown in Figure 5.

!Code definitions for the [[17, 1, 5]] colour code bbox=0.52,0.07,0.93,0.25

As a starting point we begin with a circuit following the structure given in [10]. This is a circuit composed of (a) non-fault-tolerant preparation of $T|+\rangle$, (b) flagged measurement of the logical stabiliser $S_X$, (c) full quantum error detection cycle. To make the circuit fault-tolerant (b) and (c) are repeated giving a final sequence (a)-(b)-(c)-(b)-(c). The quantum error detection cycle we employ measures each stabiliser generator of the code once, using a separate flag gadget for each. If any stabiliser measurement, or its flag, has a non-trivial outcome the circuit would be discarded. We first verify with Comb that this circuit is fault-tolerant, as expected.

The number of two-qubit gates in the circuit is reduced using a greedy approach that removes ancilla qubits. These ancillas correspond to the measurement of a stabiliser eigenvalue or a flag on that measurement. An ancilla is selected and removed from the circuit along with any gates that touch it. After each removal we use Comb to test if the remaining circuit is still fault-tolerant. If it is, the circuit is accepted, else it is rejected and that ancilla remains in the circuit. In either case we proceed to test subsequent ancillas until we have tested all of them.

Applying this procedure, we are able to reduce the number of two-qubit gates in the fault-tolerant $T|+\rangle$ magic state preparation for the $[\![17,1,5]\!]$ colour code from 298 to 193. The final circuit we find is shown in Figure 6.

Independent verification of the $[\![17,1,5]\!]$ circuits by simulation

To independently test Comb’s findings we simulate the magic state preparation circuits for the $[\![17,1,5]\!]$ code using SymFT [23] and compute the scaling of their logical error rates with physical error rate. In this appendix we give full details of these simulations. The code definitions we adopt for $[\![17,1,5]\!]$ are shown in Figure 5.

The simulation scheme we use is sketched in Figure 7(b). It begins by simulating the state preparation circuit with circuit level noise. The circuit noise makes all operations faulty with the same physical error rate $p$: that is we insert bit-flip qubit initialisation and measurement errors, as well as single qubit and two qubit depolarising channels on gates, all with rate $p$. Following this we apply a sequence of noiseless measurements and rotations that allows us to compute the logical fidelity.

!Simulation scheme outline and numerical results for the optimised $d = 5$ circuit for the $[[17, 1, 5]]$ code bbox=0.10,0.06,0.92,0.40

The logical fidelity is the overlap between the output state of the circuit and the ideal state $T|+\rangle$. Characterising the logical fidelity by reconstructing the logical Bloch vector $(\langle\bar{X}\rangle,\langle\bar{Y}\rangle,\langle\bar{Z}\rangle)$ and computing its projection onto the $T|+\rangle$ state would require a very large number of samples: $\sim10^{14}$ for logical error rates $\sim10^{-7}$. To avoid this sampling cost we measure directly in the $\{T|+\rangle,T|-\rangle\}$ basis, Figure 7(c).

Measuring in the $\{T|+\rangle,T|-\rangle\}$ basis is achieved by adding the following noiseless operations to the end of the simulation, Figure 7(b): (1) extracting the full syndrome $\vec{S}$, (2) applying a logical $T^\dagger$ rotation and (3) measuring in either the logical $\bar{X}$ or $\bar{Y}$ basis as determined by $\vec{S}$. We decode the syndrome $\vec{S}$ with a code-capacity decoder (described below), this yields a physical Pauli error $E$ that we believe to be equivalent, up to stabilisers, to the error that is present on the state. If we encounter a syndrome that is not in the code-capacity decoder we reject the sample as an additional post-selection. If $E$ commutes with the logical $\bar{Z}$ then the logical $T^\dagger$ acts correctly and $T|+\rangle$ is rotated back to $|+\rangle$ and so we measure in the $\bar{X}$ basis. However, if $E$ anti-commutes with $\bar{Z}$ then $T^\dagger$ is flipped to $T$ and $T|+\rangle$ is rotated instead to $|+i\rangle$, in this case we measure in the $\bar{Y}$ basis. To avoid the need to dynamically choose the final measurement basis at runtime we run the simulation twice: once with the measurement in the $\bar{X}$ basis and again in the $\bar{Y}$ basis. When we post-

process the data, we use the syndrome of each sample to determine whether its measurement corresponds to the $\{T\lvert+\rangle,T\lvert-\rangle\}$ basis or the rotated $\{T\lvert+i\rangle,T\lvert-i\rangle\}$ basis and aggregate the statistics accordingly.

The code capacity decoder we use for the final noiseless decoding step is constructed as a lookup table over all the correctable Pauli errors of the code. The set of all correctable Pauli errors is taken to be all Paulis with CSS weight $\le2$, meaning that if we split them $P = P^{X}P^{Z}$ as is Equation (4) we require the weight (the number of non-identity elements) of both $P^{X}$ and $P^{Z}$ to be $\le2$. To give a clarifying example the CSS weight of $X_{1}X_{4}X_{5}$ would be 3 while $X_{1}Z_{4}X_{5}$ would only be 2. We iterate over all correctable errors and calculate their syndrome. If the syndrome is not already in the lookup table we add it, else we do nothing.

The logical error rate is the number of samples where we measure $T\lvert-\rangle$ when measuring in the $\{T\lvert+\rangle,T\lvert-\rangle\}$ basis. This is plotted versus physical error rate $p$ for the optimised $d=5$ circuit in Figure 7(d), and for the full set of $\{17,1,5\}$ circuits in Figure 3. Error bars shown on the logical error rates are 95% Wilson confidence intervals. We plot the measurements in the two complementary basis of the optimised $d=5$ circuit in Figure 7(e), with error bars showing the 95% Wald confidence intervals. Finally the acceptance rates of post-selection are plotted for the optimised $d=5$ circuit in Figure 7(a), and for all circuits in Figure 3.

We simulate the circuits at error rates $p = 10^{-3}$, $10^{-2.666}$, $10^{-2.333}$, $10^{-2}$ with the following numbers of samples: for the non-FT and $d=3$ circuits $10^{7}$, $10^{6}$, $10^{6}$, $10^{6}$; for the $d=5$ circuits $10^{8}$, $5 \times10^{7}$, $10^{7}$, $10^{7}$. The most expensive simulation, the $d=5$ unoptimised circuit at $p = 10^{-3}$, took 63 minutes on an Apple Macbook Pro with Apple M4 processor.