Three-dimensional color code thresholds via statistical-mechanical mapping
Aleksander Kubica, Michael E. Beverland, Fernando Brandao, John Preskill, Krysta M. Svore
TThree-dimensional color code thresholds via statistical-mechanical mapping
Aleksander Kubica, Michael E. Beverland,
1, 2
Fernando Brand˜ao,
2, 1
John Preskill, and Krysta M. Svore Institute for Quantum Information & Matter, California Institute of Technology, Pasadena, CA 91125, USA Station Q Quantum Architectures and Computation Group, Microsoft Research, Redmond, WA 98052, USA (Dated: August 25, 2017)Three-dimensional (3D) color codes have advantages for fault-tolerant quantum computing, suchas protected quantum gates with relatively low overhead and robustness against imperfect mea-surement of error syndromes. Here we investigate the storage threshold error rates for bit-flip andphase-flip noise in the 3D color code on the body-centererd cubic lattice, assuming perfect syndromemeasurements. In particular, by exploiting a connection between error correction and statisticalmechanics, we estimate the threshold for 1D string-like and 2D sheet-like logical operators to be p (1)3DCC (cid:39) .
9% and p (2)3DCC (cid:39) . PACS numbers: 03.67.Pp, 03.67.Lx, 11.15.Ha, 75.40.Mg, 75.50.Lk
Some approaches to building scalable quantum com-puters are more practical than others due to their morefavorable noise and resource requirements. The two-dimensional (2D) surface code approach [1–3] has verydesirable features: (1) geometrically local syndrome mea-surements, (2) a high accuracy threshold and (3) fault-tolerant Clifford gates with low overhead. Unfortunately,the surface code is not known to admit a (4) fault-tolerantnon-Clifford gate with low overhead. The formidablequbit overhead cost of state distillation [4, 5] for the nec-essary non-Clifford gate motivates the quest for alterna-tives to the surface code with all features (1)–(4).Such alternatives may be sought in the general class oftopological codes [1, 2, 6–8], which includes the surfacecode as a special case. By definition, topological codesrequire only geometrically local syndrome measurementsand tend to have high accuracy thresholds. Topologicalcodes often admit some fault-tolerant transversal gates(implemented by the tensor product of single-qubit uni-taries), which have low overhead cost. However, no quan-tum error-detecting code (whether topological or not),has a universal transversal encoded gate set [9, 10].Here we focus on the 3D topological color codes [11, 12]closely related to the 3D toric code [13], which comein two types. The stabilizer type has 1D string-like Z and 2D sheet-like X logical operators, and a logicalnon-Clifford gate T = diag(1 , e iπ/ ) is transversal. Inthe subsystem type, there are 1D string-like X and Z dressed logical operators, and all logical Clifford gatesare transversal. Moreover, in the subsystem color codeit is possible to reliably detect measurement errors in asingle time step [14, 15]. By fault-tolerantly switching be-tween the stabilizer and subsystem color codes [12, 16],one can combine the desirable features (1), (3) and (4).In this work, we address feature (2) for the 3D colorcodes by finding thresholds p (1)3DCC (cid:39) .
9% and p (2)3DCC (cid:39) .
6% for phase-flip Z and bit-flip X noise, respectively.Our results assume optimal decoders for independent X and Z noise with perfect measurements, and thereby givefundamental error-correction bounds against which effi- cient, but suboptimal decoders (such as that studied in[15]) can be compared. These thresholds are comparableto the analogous thresholds for the cubic lattice 3D toriccode: p (1)3DTC (cid:39) .
3% and p (2)3DTC (cid:39) .
5% [17–19], butcompare unfavorably to p (cid:39) .
9% for the squarelattice 2D toric code [20].Our approach extends techniques known for othercodes [3, 8, 21–25] in order to relate the 3D color codethresholds to phase transitions in two new 3D statistical-mechanical models: the 4- and 6-body random couplingIsing models (RCIM). We use large-scale parallel tem- s e1 s e2 s e5 s e6 s e4 s e2 s v1 s v4 s v2 s v3 FIG. 1. The disorder-temperature ( p, T )-phase diagrams ofthe 4-body (top) and 6-body (bottom) 3D random couplingIsing models. Both models are defined on the 3D body-centered cubic lattice built of tetrahedra. The 4- and 6-bodymodels have spins on vertices and edges, respectively. Theerror correction threshold p c can be found as the intersectionof the Nishimori line (blue line) with the anticipated phaseboundary (red dotted line). a r X i v : . [ qu a n t - ph ] A ug pering Monte Carlo simulations [26] and analyze specificheat, sublattice magnetization and Wilson loop operatorsto map the relevant parts of the disorder-temperature( p, T )-phase diagram; see Fig. 1. The 6-body RCIM isan example of a lattice gauge theory with a local (gauge) Z × Z symmetry, which makes this model both inter-esting and challenging to study.
3D stabilizer color code.—
Let L be a three-dimensional lattice built of tetrahedra such that its ver-tices are 4-colorable, i.e., vertices connected by an edgeare of different colors. An example of such a lattice isthe body-centered cubic (bcc) lattice obtained from twointerleaved cubic lattices; see Fig. 2(b). We denote by∆ i ( L ) the set of all i -simplices of L . Then, 0-simplicesof L are vertices, 1-simplices are edges, etc. We placeone qubit at every tetrahedron t ∈ ∆ ( L ). For every ver-tex v ∈ ∆ ( L ) and edge e ∈ ∆ ( L ) we define operators S X ( v ) and S Z ( e ) to be the product of either Pauli X or Z operators on qubits identified with tetrahedra in theneighborhood of the vertex v or edge e , namely S X ( v ) = (cid:89) t ∈ ∆ ( L ) t ⊃ v X ( t ) , S Z ( e ) = (cid:89) t ∈ ∆ ( L ) t ⊃ e Z ( t ) . (1)The 3D stabilizer [27] color code is defined by specifyingits stabilizer group [28] S = (cid:104) S X ( v ) , S Z ( e ) | v ∈ ∆ ( L ) , e ∈ ∆ ( L ) (cid:105) . (2)Using the colorability condition one can show that S isan Abelian subgroup of the Pauli group not containing − I . The code space is the +1 eigenspace of all elementsof S and the lowest-weight logical X and Z operators ofthe 3D color code are 2D sheet-like and 1D string-likeobjects; see Fig. 2(a). In general, the color code can bedefined in d ≥ d + 1)-colorable simplical d -complex [16]. (b) (c) (a) FIG. 2. (a) The 3D stabilizer color code has both 1D string-like (red) and 2D sheet-like (blue) logical operators. (b) Thebcc lattice can be constructed starting from two interleavedcubic lattices (red and blue) and filling in tetrahedra (green).Not all tetrahedra are depicted. (c) The neighborhood ofany vertex in the bcc lattice looks the same — every vertexbelongs to 24 edges, 36 triangular faces and 24 tetrahedra.The bcc lattice is 4-colorable, i.e., every vertex is colored inred, green, blue or yellow, and no two neighboring verticesare of the same color.
Error correction in CSS codes.—
Since the color codeis a CSS code [29], we choose to separately correct X -and Z -type errors, which simplifies the discussion. We also assume perfect measurements. For concreteness, wefocus on X -error correction; Z -errors can be analyzedanalogously [30].The set of all Z -type stabilizers which return − Z -type syndrome. Notethat any nontrivial Z -syndrome signals the presence ofsome X -errors in the system. Correction of X -errors ina CSS code can be succinctly described by introducing achain complex [31, 32] C ∂ −→ C ∂ −→ C X -stabilizers qubits Z -stabilizers (3)where C , C and C are vector spaces over Z with bases B = X -stabilizer generators, B = physical qubits and B = Z -stabilizer generators, respectively. The linearmaps ∂ and ∂ , called boundary operators, are chosen insuch a way that the support of any X -stabilizer ω ∈ C isgiven by ∂ ω , and the Z -syndrome corresponding to any X -error (cid:15) ∈ C can be found as ∂ (cid:15) . Note that ∂ ◦ ∂ = 0,since any X -stabilizer has trivial Z -syndrome. One canthink of the boundary operators as parity-check matrices H TX and H Z of the CSS code. In the case of the 3D colorcode, C , C , C are generated by vertices, tetrahedra,and edges respectively, i.e., B = ∆ ( L ), B = ∆ ( L )and B = ∆ ( L ). The boundary operators are defined tobe ∂ v = (cid:80) ∆ ( L ) (cid:51) t ⊃ v t and ∂ t = (cid:80) ∆ ( L ) (cid:51) e ⊂ t e for any v ∈ ∆ ( L ) and t ∈ ∆ ( L ).Let (cid:15), ϕ ∈ C be two X -errors with the same Z -syndrome, ∂ (cid:15) = ∂ ϕ . We say that (cid:15) and ϕ are equiva-lent iff they differ by some X -stabilizer ω ∈ C , namely (cid:15) + ϕ = ∂ ω . To correct errors, we need a decoder —an algorithm which takes the Z -syndrome σ ∈ C as aninput and returns a Z -correction ϕ which will restore all X -stabilizers to have +1 outcomes, i.e., ∂ ϕ = σ . Thedecoder succeeds iff the actual error (cid:15) and the correction ϕ are equivalent. An optimal decoder finds a represen-tative ϕ of the most probable equivalence class of errors ϕ = { ϕ + ∂ ω |∀ ω ∈ C } . Statistical-mechanical models.—
In this section, weprovide a brief derivation of the connection between op-timal error-correction thresholds and phase transitions[3, 8, 21–25]. In particular, we derive two new statistical-mechanical models relevant for the 3D color code.We assume bit-flip noise, i.e., every qubit is indepen-dently affected by Pauli X error with probability p . Theprobability of an X -error (cid:15) ∈ C affecting the system ispr( (cid:15) ) = (cid:89) j ∈B p [ (cid:15) ] j (1 − p ) − [ (cid:15) ] j ∝ (cid:18) p − p (cid:19) (cid:80) j ∈B [ (cid:15) ] j , (4)where [ (cid:15) ] j ∈ Z denotes the j coefficient of (cid:15) in the B basis, (cid:15) = (cid:80) j ∈B [ (cid:15) ] j j .For a general CSS code family with the chain complexin Eq. (3), the X -error correction threshold is the largest p c such that for all p < p c the probability of successfuldecoding goes to 1 in the limit of infinite system sizepr(succ) = (cid:88) (cid:15) ∈ C pr( (cid:15) )pr(succ | (cid:15) ) → . (5)With the optimal decoder, the conditional probabilitypr(succ | (cid:15) ) equals 1 if (cid:15) belongs to the most probable errorequivalence class consistent with the syndrome ∂ (cid:15) , and0 otherwise. The probability of equivalence class (cid:15) ispr( (cid:15) ) = (cid:88) ω ∈ C pr( (cid:15) + ∂ ω ) ∝ (cid:88) ω ∈ C e − β ( p ) (cid:80) j ∈B [ (cid:15) + ∂ ω ] j , (6)where we use Eq. (4) and introduce β ( p ) = −
12 log p − p . (7)To rewrite Eq. (6), we use [ ∂ ω ] j ≡ (cid:80) i ∈B ∧ ∂ i (cid:51) j [ ω ] i mod 2 and 1 − (cid:15) + ∂ ω ] j = ( − [ (cid:15) ] j ( − [ ∂ ω ] j =( − [ (cid:15) ] j (cid:81) i ∈B ∧ ∂ i (cid:51) j ( − [ ω ] i . By introducing new (clas-sical spin) variables s i = ( − [ ω ] i for all i ∈ B , we canreplace the sum over ω ∈ C in Eq. (6) by a sum overdifferent configurations { s i = ± } , yieldingpr( (cid:15) ) ∝ (cid:88) { s i = ± } e − β ( p ) H (cid:15) ( { s i } ) , (8)where we introduce the Hamiltonian H (cid:15) ( { s i } ) = − (cid:88) j ∈B ( − [ (cid:15) ] j (cid:89) i ∈B [ ∂ i ] j =1 s i . (9)We define the random coupling Ising model (RCIM)to be a classical spin s i = ± − [ (cid:15) ] j described by H (cid:15) ( { s i } ) inEq. (9). The RCIM has two independent parameters:disorder strength p (i.e., the probability of negative cou-plings) and inverse temperature β . The partition func-tion of the RCIM with disorder (cid:15) at temperature β − isgiven by Z (cid:15) ( β ) = (cid:88) { s i = ± } e − βH (cid:15) ( { s i } ) . (10)Note that for the proportionality pr( (cid:15) ) ∝ Z (cid:15) ( β ) inEq. (10) to hold one requires β = β ( p ).For the 3D color code, Eq. (9) leads to the followingtwo new statistical-mechanical models H X(cid:15) ( { s v } ) = − (cid:88) t ∈ ∆ ( L ) ( − [ (cid:15) ] t s a s b s c s d s e s f , (11) H Z(cid:15) ( { s e } ) = − (cid:88) t ∈ ∆ ( L ) ( − [ (cid:15) ] t s a s b s c s d , (12)relevant to correction of X - and Z -errors, respectively.Note that H X(cid:15) ( { s v } ) (respectively H Z(cid:15) ( { s e } )) contains 4-body (6-body) terms, which are products of vertex (edge)spins of every tetrahedron. We observe that for p = 0,i.e., the case with no disorder, these two models are self-dual in the sense that the low-temperature expansion ofeach model matches the high-temperature expansion ofthe other [33]; see the Supplemental Material. The Hamiltonian in Eq. (9) determines a thermal en-semble of excitations in the statistical mechanical model.For H X(cid:15) ( { s v } ) the excitations are 2D domain walls resid-ing on a set of tetrahedra ϕ = (cid:15) + ∂ ω ∈ C , where thesewalls terminate at the edges contained in ∂ ϕ = ∂ (cid:15) ∈ C .In the color code, this ensemble of domain walls corre-sponds to the ensemble of possible X -errors which gen-erate the same error syndrome as (cid:15) , and the Boltzmannweight of a wall configuration coincides with the proba-bility of the corresponding X -error configuration ϕ . Like-wise, for H Z(cid:15) ( { s e } ) the excitations are 1D strings termi-nating at vertices in ∂ (cid:15) , corresponding to Z -errors whichgenerate the same error syndrome as (cid:15) .To determine the storage threshold for the 3D colorcode, we investigate the disorder-temperature ( p, T )-phase diagram of the RCIM in Eq. (9). In the orderedphase, large fluctuations of domain walls (or strings) aresuppressed [3], and the free energy cost∆ λ ( (cid:15) ) = − log Z (cid:15) + λ ( β ) + log Z (cid:15) ( β ) (13)of introducing any non-trivial domain wall λ ∈ ker ∂ \ im ∂ to the system at inverse temperature β with dis-order (cid:15) should diverge in the limit of infinite system sizewhen averaged over all disorder configurations (cid:104) ∆ λ (cid:105) = (cid:88) (cid:15) ∈ C pr( (cid:15) )∆ λ ( (cid:15) ) → ∞ . (14)Correspondingly, in the color code, the error ϕ producesa syndrome ∂ ϕ which points to a unique equivalenceclass ϕ , so that the syndrome can be decoded success-fully with high probability. Indeed, we show in the Sup-plemental Material, pr(succ) → p implies (cid:104) ∆ λ (cid:105) → ∞ for the RCIM at inverse temperature β ( p ) and disorder strength p . Thus, by finding the criti-cal point along the line defined by Eq. (7) (the Nishimoriline [34]) we obtain the threshold value p c . Phase diagram.—
We describe how to map out the( p, T )-phase diagrams of the two RCIMs, H X(cid:15) ( { s v } ) and H Z(cid:15) ( { s e } ). The discontinuity in energy density acrossa first order phase transition allows for straightforwardidentification of the phase boundary in the regime of lowdisorder. However, more reliable order parameters arerequired to probe a (higher-order) phase transition closeto the critical point on the Nishimori line. Moreover,an appropriate order parameter takes symmetries of themodel into account. Note that flipping a subset of spins { s i } i ∈ I , i.e., s i (cid:55)→ − s i for i ∈ I , is a symmetry if it leavesthe Hamiltonian describing the model invariant.The 4-body RCIM in Eq. (11) has a global Z × Z × Z symmetry. An example of a symmetry operation is a si-multaneous flip of vertex spins on all red and blue ver-tices, since it leaves every term of H X(cid:15) ( { s v } ) unchanged.Due to this symmetry, the total magnetization is not agood order parameter; however the sublattice magneti-zation of spins of a single color is. Instead of using thesublattice magnetization directly, more precise estima-tions are obtained by considering the finite-size scaling -1 L=6L=8L=10L=12L=16 -2 -1 L=6L=8L=10L=12 (d) (e) (f)(a) (b) (c)
FIG. 3. Results for the 3D 4-body (a)-(c) and 6-body (d)-(f) RCIM. By finding the peak positions of specific heat c L fordifferent system sizes L and exploiting finite-size scaling we estimate for p = 0 the critical temperature of a phase transition in(a) and (d) to be T c = 8 . T c = 0 . p = 0 .
27 we identify T c = 2 . ξ L /L for different system sizes L . (c) For p = 0 .
28 there is no indication of a phase transition. In(e) and (f) we check if the Wilson loop operator W ( γ ) satisfies the perimeter law by plotting − log (cid:104) W ( γ ) (cid:105) /P ( γ ) as a functionof perimeter P ( γ ) of the square loop γ for different temperatures T . (e) For p = 0 .
018 we see a change of scaling as the systemundergoes a phase transition at T = 0 . p = 0 .
021 there is no indication of a phase transition. of the spin-spin correlation function [35]. Near the phasetransition, for fixed disorder strength p and temperatures T close to the critical temperature T c ( p ), the correlationlength ξ L is expected to scale as ξ L ( p, T ) /L ∼ f ( L /ν ( T − T c ( p ))) , (15)where L is the linear system size, f is a scaling functionand ν is the correlation length critical exponent [36]. Wecan estimate T c ( p ) by plotting ξ L ( p, T ) /L as a functionof temperature T for different system sizes L and findingtheir crossing point; see Fig. 3(a)(b). If no crossing isobserved, then we conclude that there is no phase tran-sition.The 6-body RCIM in Eq. (12) describes a lattice gaugetheory with a local Z × Z symmetry. An example of asymmetry operation is a flip of edge spins on edges froma single yellow vertex to all neighboring red and bluevertices; see Fig. 2(c). Due to Elitzur’s theorem [37], thegauge symmetry rules out existence of any local orderparameter. We define a Wilson loop operator [38, 39] W ( γ ) = (cid:89) e ∈ γ s e , (16)to be a product of edge spins along a loop γ ⊂ ∆ ( L ).For W ( γ ) to be gauge-invariant the loop γ can only becomposed of edges connecting vertices of two (out of fourpossible) colors. The phase transition can be identified by analyzing scaling of the thermal expectation value of W ( γ ) averaged over different disorder configurations (cid:104) W ( γ ) (cid:105) = (cid:88) (cid:15) ⊂ ∆ ( L ) pr( (cid:15) ) (cid:88) { s e } W ( γ ) e − βH Z(cid:15) ( { s e } ) Z (cid:15) ( β ) . (17)Namely, in the limit of large square loops [19, 21, 40], − log (cid:104) W ( γ ) (cid:105) scales linearly with the loop’s perimeter P ( γ ) in the ordered (Higgs) phase, whereas in the dis-ordered (confinement) phase it scales linearly with theminimum area A ( γ ) enclosed by γ ; see Fig. 3(d)-(f).We find the ( p, T )-phase diagrams of the 4- and 6-bodyRCIMs by performing Monte Carlo simulations with par-allel tempering technique [26]; see Fig. 1. We test equi-libration of the system by a logarithmic binning of thedata. Since we can simulate only finite-size systems, acareful analysis of finite-size effects is necessary. Param-eters of numerical simulations and additional details areprovided in the Supplemental Material. Discussion.—
We have found 3D stabilizer color codethresholds for phase-flip Z and bit-flip X noise mod-els with optimal decoding and perfect measurements: p (1)3DCC (cid:39) .
9% and p (2)3DCC (cid:39) . X -stabilizersdetecting Z -errors are the same for the 3D stabilizer andsubsystem color codes. Since the subsystem code is sym-metric under the exchange of X - and Z -generators, itsphase- and bit-flip thresholds are the same and equal to p (1)3DCC of the stabilizer color code on the same lattice fam-ily. The 3D color code threshold [41] with the (efficient)clustering decoder p (1)clust (cid:39) .
46% [15] is about a fourthof p (1)3DCC , justifying a search for efficient color-code de-coders with performance as close to optimal as for 2Dsurface and color codes.We hope that our work initiates and motivates a care-ful study of the 3D random coupling Ising models andtheir phase diagrams. We conjecture the existence of aspin-glass phase [42] in the 6-body RCIM, correspondingto a regime with intermediate noise strength in whichmemory lifetime with non-optimal decoders is polyno-mial rather than exponential in the system size.A future extension of this work might incorporate mea-surement errors in the analysis. This would require the study of 4D random coupling models and thus use morecomputational resources. If successful, this research pro-gram could provide a deeper understanding of single-shoterror correction [14, 15] from the standpoint of statisticalmechanics.We thank R. Andrist, H. Bomb´ın, N. Delfosse, L.Pryadko, B. Yoshida and I. Zintchenko for helpful discus-sions. AK would like to thank the QuArC group for theirhospitality during a summer internship. We acknowledgefunding provided by the Institute for Quantum Informa-tion and Matter, an NSF Physics Frontiers Center (NFSGrant PHY-1125565) with support of the Gordon andBetty Moore Foundation (GBMF-12500028). Appendix A: Duality of models for zero disorder
We already mentioned that the 4- and 6-body RCIM described by Eqs. (11) and (12) are dual for p = 0, i.e., the casewith no disorder . Here we say that two models are dual if the low-temperature expansion of the partition functionof one model matches the high-temperature expansion of the partition function of the other and vice versa [33]. Weobserve that for any CSS code, the two statistical-mechanical models relevant for correction of X - and Z -errors arealways dual for p = 0. In particular, if there is only one phase transition in the first model at temperature T Xc , thenthere is a unique phase transition in the dual model at temperature T Zc = − T Xc log tanh 1 T Xc . (A1)This serves as a consistency check for our results. Indeed, for zero disorder p = 0 the critical temperatures T Xc =8 . T Zc = 0 . Appendix B: Proof of implication
Here we show that successful decoding implies diverging average energy cost of introducing any non-trivial domainwall. We used this fact in the derivation of statistical-mechanical models to relate the threshold p c of optimal errorcorrection to the critical point p N on the Nishimori line. Note that this implication allows us to only infer that p c ≤ p N . However, we expect that successful decoding be possible throughout the ordered phase and thus these twovalues should coincide. Lemma 1.
Consider a CSS code described by the chain complex in Eq. (3). Let H = ker ∂ / im ∂ be the firsthomology group of finite cardinality, | H | < ∞ . If the probability of successful optimal X -error correction goes to 1 inthe limit of infinite system size pr(succ) = (cid:88) (cid:15) ∈ C pr( (cid:15) )pr(succ | (cid:15) ) → , (B1) then the average free energy cost of introducing any non-trivial domain wall λ ∈ ker ∂ \ im ∂ diverges (cid:104) ∆ λ (cid:105) = (cid:88) (cid:15) ∈ C pr( (cid:15) )∆ λ ( (cid:15) ) → ∞ . (B2) Proof.
Let (cid:15) = { (cid:15) + ∂ ω | ω ∈ C } denote the equivalence class of errors for (cid:15) ∈ C and E = { (cid:15), . . . } be the set of allequivalence classes. We define a representative of the most probable equivalence class of errors consistent with thesyndrome σ ∈ C to be ρ ( σ ) = arg max ρ ∈ C ∂ ρ = σ pr( ρ ) . (B3)The conditional probability of successful decoding using the optimal (maximum likelihood) decoder is given bypr(succ | (cid:15) ) = (cid:40) (cid:15) ∈ ρ ( ∂ (cid:15) ) , . (B4)Thus, we have pr(succ) = (cid:88) (cid:15) ∈ C pr( (cid:15) )pr(succ | (cid:15) ) = (cid:88) σ ∈ im ∂ pr( ρ ( σ )) . (B5)By rewriting the sum over all equivalence classes of errors (cid:15) ∈ E as the sum over all possible syndromes σ ∈ im ∂ anddifferent representatives λ (cid:48) ∈ H of the homology group we arrive at1 = (cid:88) (cid:15) ∈E pr( (cid:15) ) = (cid:88) σ ∈ im ∂ (cid:88) λ (cid:48) ∈ H pr( ρ ( σ ) + λ (cid:48) ) = pr(succ) + (cid:88) σ ∈ im ∂ (cid:88) (cid:54) = λ (cid:48) ∈ H pr( ρ ( σ ) + λ (cid:48) ) . (B6)We want to show two inequalitiespr(succ) ≥ (cid:88) (cid:15) ∈E pr( (cid:15) ) pr( (cid:15) ) (cid:80) λ (cid:48) ∈ H pr( (cid:15) + λ (cid:48) ) ≥ − . (B7)In order to show the first inequality (B7) note thatpr(succ) = (cid:88) σ ∈ im ∂ (cid:88) λ (cid:48)(cid:48) ∈ H pr( ρ ( σ ) + λ (cid:48)(cid:48) ) pr( ρ ( σ )) (cid:80) λ (cid:48) ∈ H pr( ρ ( σ ) + λ (cid:48) ) (B8) ≥ (cid:88) σ ∈ im ∂ (cid:88) λ (cid:48)(cid:48) ∈ H pr( ρ ( σ ) + λ (cid:48)(cid:48) ) pr( ρ ( σ ) + λ (cid:48)(cid:48) ) (cid:80) λ (cid:48) ∈ H pr( ρ ( σ ) + λ (cid:48) ) = (cid:88) (cid:15) ∈E pr( (cid:15) ) pr( (cid:15) ) (cid:80) λ (cid:48) ∈ H pr( (cid:15) + λ (cid:48) ) , (B9)where we use pr( ρ ( σ )) ≥ pr( ρ ( σ ) + λ (cid:48)(cid:48) ) for all σ ∈ im ∂ and λ (cid:48)(cid:48) ∈ H . The second inequality (B7) follows frompr(succ) = (cid:88) σ ∈ im ∂ pr( ρ ( σ )) = (cid:88) σ ∈ im ∂ pr( ρ ( σ )) pr( ρ ( σ )) (cid:80) λ (cid:48) ∈ H pr( ρ ( σ ) + λ (cid:48) ) + (cid:88) σ ∈ im ∂ pr( ρ ( σ )) (cid:80) (cid:54) = λ (cid:48) ∈ H pr( ρ ( σ ) + λ (cid:48) ) (cid:80) λ (cid:48) ∈ H pr( ρ ( σ ) + λ (cid:48) ) (B10) ≤ (cid:88) σ ∈ im ∂ (cid:88) λ (cid:48)(cid:48) ∈ H pr( ρ ( σ ) + λ (cid:48)(cid:48) ) pr( ρ ( σ ) + λ (cid:48)(cid:48) ) (cid:80) λ (cid:48) ∈ H pr( ρ ( σ ) + λ (cid:48) ) + (cid:88) σ ∈ im ∂ (cid:88) (cid:54) = λ (cid:48) ∈ H pr( ρ ( σ ) + λ (cid:48) ) (B11)= (cid:88) (cid:15) ∈E pr( (cid:15) ) pr( (cid:15) ) (cid:80) λ (cid:48) ∈ H pr( (cid:15) + λ (cid:48) ) + (1 − pr(succ)) . (B12)If pr(succ) →
1, then from inequalities (B7) we infer that (cid:88) (cid:15) ∈E pr( (cid:15) ) pr( (cid:15) ) (cid:80) λ (cid:48) ∈ H pr( (cid:15) + λ (cid:48) ) → , (B13)and thus for λ ∈ ker ∂ \ im ∂ we have (cid:88) (cid:15) ∈E pr( (cid:15) ) pr( (cid:15) + λ ) (cid:80) λ (cid:48) ∈ H pr( (cid:15) + λ (cid:48) ) → . (B14)In the last step we used the following inequalities0 ≤ (cid:88) (cid:15) ∈E pr( (cid:15) ) pr( (cid:15) + λ ) (cid:80) λ (cid:48) ∈ H pr( (cid:15) + λ (cid:48) ) ≤ (cid:88) (cid:15) ∈E pr( (cid:15) ) (cid:80) (cid:54) = λ (cid:48) ∈ H pr( (cid:15) + λ (cid:48) ) (cid:80) λ (cid:48) ∈ H pr( (cid:15) + λ (cid:48) ) = 1 − (cid:88) (cid:15) ∈E pr( (cid:15) ) pr( (cid:15) ) (cid:80) λ (cid:48) ∈ H pr( (cid:15) + λ (cid:48) ) . (B15)We rewrite (cid:104) ∆ λ (cid:105) in the following way (cid:104) ∆ λ (cid:105) = (cid:88) (cid:15) ∈ C pr( (cid:15) )∆ λ ( (cid:15) ) = (cid:88) (cid:15) ∈E pr( (cid:15) ) log pr(¯ (cid:15) )pr( (cid:15) + λ ) (B16)= (cid:88) (cid:15) ∈E pr( (cid:15) ) log pr(¯ (cid:15) ) (cid:80) λ (cid:48) ∈ H pr( (cid:15) + λ (cid:48) ) − (cid:88) (cid:15) ∈E pr( (cid:15) ) log pr( (cid:15) + λ ) (cid:80) λ (cid:48) ∈ H pr( (cid:15) + λ (cid:48) ) . (B17)Using the inequality log x ≥ − x to lower-bound the first term and Jensen inequality for the second term we obtain (cid:104) ∆ λ (cid:105) ≥ (1 − | H | ) − log (cid:88) (cid:15) ∈E pr( (cid:15) ) pr( (cid:15) + λ ) (cid:80) λ (cid:48) ∈ H pr( (cid:15) + λ (cid:48) ) → ∞ . (B18) Appendix C: Finding phase transitions
In order to map the disorder-temperature phase diagrams of the 4- and 6-body RCIM in Fig. 1 we need to reliablyidentify phase transitions. Here we describe in detail how we achieve that by analyzing specific heat, the spin-spincorrelation function and the Wilson loop operator. We aslo present additional results for the 4- and 6-body RCIM inFig. 4. -3 (d) (e) (f)(a) (b) (c) FIG. 4. Additional details about the 4-body (a)-(c) and 6-body (d)-(f) RCIM. The discontinuity in energy per spin
E/N in(a) and (d) suggests first-order phase transitions for both models for p = 0. (b) For p = 0 the normalized correlation length ξ L /L does not seem to be described well by the scaling ansatz in Eq. (15) possibly due to a transition being first-order. (c)For the disorder value p = 0 .
276 close to the critical point on the Nishimori line p N = p (2)3DCC detecting a phase transition andestimating its critical temperature becomes difficult. (e) We check if the Wilson loop operator W ( γ ) satisfies the perimeter lawby plotting − log (cid:104) W ( γ ) (cid:105) /P ( γ ) as a function of perimeter P ( γ ) of the square loop γ for different disorder values p and fixedtemperature T = 0 .
42. We see a change of scaling as the system undergoes a phase transition at p = 0 . − log (cid:104) W ( γ l ) (cid:105) /l ∼ al + b + c log l to the data for p = 0 .
018 in Fig. 3(c) and plot the fit coefficient a as a function of temperature T . We identify the critical temperature T c = 0 . a = 0.
1. Specific heat
For a second-order phase transition, the specific heat c ( T ) as a function of temperature T is expected to have adiscontinuity near a phase transition at temperature T c in the limit of infinite system size L → ∞ . However, for asystem of finite linear size L , the peak of the specific heat c L ( T ) appears at temperature T c ( L ) = arg max T c L ( T )shifted from that in the infinite system by an amount (cid:12)(cid:12)(cid:12)(cid:12) T L − T c T c (cid:12)(cid:12)(cid:12)(cid:12) ∝ L − /ν , (C1)where ν is the correlation length critical exponent [36]. A similar scaling behavior has been established for first-orderphase transitions [43–46]. Thus, we find the critical temperature T c by fitting a function T c ( L ) ∼ aL − b + T c (C2)to the position of the specific heat peaks for different system sizes and evaluating T c ( L = ∞ ).
2. Correlation function
One might not be able to identify a phase transition of higher order by looking at the specific heat. Rather, oneneeds to analyze the behavior of e.g. the order parameter correlation length ξ . In particular, for the system of finitesize L and with fixed disorder strength p we define the two-point finite-size correlation length ξ L as a function oftemperature T ξ L ( T ) = 12 sin( k / (cid:115) (cid:104) χ ( (cid:126) (cid:105)(cid:104) χ ( (cid:126)k ) (cid:105) − , (C3)where (cid:104) χ ( (cid:126)k ) (cid:105) = (cid:80) (cid:15) ⊂ ∆ ( L ) pr( (cid:15) ) χ ( (cid:126)k ), (cid:126)k is the wavevector and (cid:126)k = (2 π/L, , χ ( (cid:126)k ) = (cid:88) { s v } N (cid:32)(cid:88) u ∈ U s u e i(cid:126)k · (cid:126)r u (cid:33) e − βH X(cid:15) ( { s v } ) Z (cid:15) ( β ) . (C4)where (cid:126)r u denotes the position of the vertex spin s u in a sublattice U ⊂ ∆ ( L ) of single-color vertices. Near a phasetransition at temperature T c , the normalized correlation length is expected to scale as ξ L ( T ) L ∼ f ( L /ν ( T − T c )) , (C5)where f is a dimensionless scaling function and ν is the correlation length critical exponent. We can estimate T c byplotting ξ L ( T ) /L as a function of temperature T for different system sizes L and finding their crossing point. If thereis no crossing, then we conclude that there is no phase transition.
3. Wilson loop operator
When the system under consideration has a local (gauge) symmetry, one cannot use a local order parameter todetect a phase transition. Rather, one needs to consider gauge-invariant quantities, such as the Wilson loop operator W ( γ ) in Eq. (16). Suppose γ is a square loop. We denote by P ( γ ) and A ( γ ) the perimeter of γ and the minimalarea enclosed by γ , respectively. The scaling of the averaged Wilson loop operator (cid:104) W ( γ ) (cid:105) in the limit of large loopschanges between the ordered (Higgs) and disordered (confinement) phases. Namely, • in the disordered phase: (cid:104) W ( γ ) (cid:105) ∼ exp( − const · A ( γ )), • in the ordered phase: (cid:104) W ( γ ) (cid:105) ∼ exp( − const · P ( γ )).We consider a system of finite size L and denote by γ l a square loop of linear size l ≤ L/
2. Since A ( γ l ) ∝ l and P ( γ l ) ∝ l , then log (cid:104) W ( γ l ) (cid:105) should scale either quadratically or linearly in l , depending on the phase of the system.Due to finite-size effects, there are some corrections to the area and perimeter scaling. In particular, we numericallyfind that − log (cid:104) W ( γ l ) (cid:105) l ∼ al + b + c log l, (C6)where a, b, c are some constants. We identify the disordered phase as the region where the fitting parameter a ispositive, a > Appendix D: Classical Ising gauge theory
As an example of using specific heat and the scaling of the Wilson loop operator to identify a phase transitionwe study a known model, the three-dimensional random plaquette Ising model (RPIM); see Fig. 5. The RPIM is ageneralization of the Z Ising gauge theory, which is relevant for studying the optimal error correction threshold for1D string-like operators in the 3D toric code [3]. The RPIM is a statistical-mechanical model with classical spins s e = ± e ∈ ∆ ( C ) of the cubic lattice C and disorder (cid:15) ⊂ ∆ ( C ). The Hamiltonian describing theRPIM H RPIM (cid:15) ( { s e } ) = − (cid:88) f ∈ ∆ ( C ) ( − [ (cid:15) ] f (D1)contains 4-body terms, which are products of four edge spins around every square face f ∈ ∆ ( C ) of the lattice C . Weset [ (cid:15) ] f = 1 if f ∈ (cid:15) , otherwise [ (cid:15) ] f = 0. We observe that H RPIM (cid:15) ( { s e } ) has a local Z symmetry, generated by flips ofspins on all edges incident on any vertex v ∈ ∆ ( C ). The Wilson loop operator W ( γ l ) is a gauge-invariant quantity,where γ l is a square loop of linear size l . The disorder-temperature phase diagram of the 3D RPIM is shown in Fig. 6. -3 (d) (e) (f)(a) (b) (c) FIG. 5. Results for the 3D RPIM. (a) For p = 0 we can estimate the critical temperature T c = 1 . c L for different system sizes L and exploiting finite-size scaling. In (b)-(d) wecheck for p = 0, p = 0 .
031 and p = 0 .
035 whether the Wilson loop operator W ( γ ) satisfies the perimeter law by plotting − log (cid:104) W ( γ ) (cid:105) /P ( γ ) as a function of perimeter P ( γ ) of the square loop γ for different temperatures T . (e) For fixed temperature T = 0 .
45 we analyze scaling of − log (cid:104) W ( γ ) (cid:105) /P ( γ ) for different disorder values p . (f) We find a fit − log (cid:104) W ( γ l ) (cid:105) /l ∼ al + b + c log l to the data in (c) and plot the fit coefficient a as a function of temperature T . We identify the critical temperature T c = 0 . a = 0. In (b),(c) and (e) we see a change of scaling as the system undergoes aphase transition at T c = 1 . T c = 0 . p c = 0 . FIG. 6. The disorder-temperature ( p, T )-phase diagram of the 3D random plaquette Ising model on the cubic lattice. Theintersection of the Nishimori line (blue) with the anticipated phase boundary (red dotted line) gives the 3D toric code threshold p (1)3DTC (cid:39) .
3% for optimal error correction associated with 1D string-like logical operators (and point-like excitations). Notethat the location of a phase transition for T = 0 was found in [21]. Appendix E: Numerical simulation details
The numerical complexity of simulating the statistical-mechanical models, such as the 4- and 6-body RCIM and theRPIM, increases with the disorder strength p , which is reminiscent of a spin glass behavior. To speed up simulations weuse the parallel tempering technique. The parallel tempering technique requires simultaneous simulation of multiplecopies k = 1 , . . . , n of the system with the same disorder (cid:15) but different spin configurations { s i } k and temperatures T < . . . < T n . After performing single-spin Metropolis updates for all spins in every copy of the system, swaps ofspin configurations { s i } k ↔ { s i } k +1 of copies at neighboring temperatures T k and T k +1 are allowed with probabilitypr( k ↔ k + 1) = exp (cid:18) ( E k − E k +1 ) (cid:18) T k − T k +1 (cid:19)(cid:19) , (E1)where E k and E k +1 denote energies of spin configurations { s i } k and { s i } k +1 . We choose temperatures T < . . . < T n is such a way that the exchange rate { s i } k ↔ { s i } k +1 is approximately flat; for more in-depth discussions see e.g. [47].Equilibration of the system is tested by a logarithmic binning of data. Numerical simulation details for the 4-bodyRCIM, the 6-body RCIM and the RPIM are presented in Table I.To estimate statistical error bars of quantities analyzed in the simulation we use the bootstrap technique. The mainidea behind the bootstrap technique is to repeat sampling from the existing data set D = { d , . . . , d N } and evaluatinga quantity of interest q = q ( D ). In particular, for i = 1 , . . . , n we perform the following steps1. from the data set D randomly choose N data points d i ( j ) , where i ( j ) ∈ { , . . . , N } ,2. evaluate the quantity q i = q ( D i ) from the data set D i = { d i (1) , . . . , d i ( N ) } .Note that in step 1 we allow to choose the same data point multiple times. The relevant quantity q is estimated to be q = ¯ q ± (cid:118)(cid:117)(cid:117)(cid:116) n (cid:88) i =1 (¯ q − q i ) n − , (E2)where ¯ q = n (cid:80) ni =1 q i .1 p L max N (cid:15) τ N T T min T max L max and N (cid:15) denote the linear size of the biggest simulated system and the number of randomly chosen disorder samples. N T denotes the number of temperatures in the range [ T min , T max ] chosen in a way that the exchange rate of spin configurations isapproximately constant. 2 τ is the number of equilibration steps, where one equilibration step consists of an update of everyspin in all N T copies of the system followed by swaps { s i } k ↔ { s i } k +1 of spin configurations. [1] A. Kitaev, Annals of Physics , 1 (2003), arXiv:9707021v1 [arXiv:quant-ph].[2] S. B. Bravyi and A. Y. Kitaev, , 6 (1998), arXiv:9811052 [quant-ph].[3] E. Dennis, A. Kitaev, A. Landahl, and J. Preskill, Journal of Mathematical Physics , 4452 (2002), arXiv:0110143[quant-ph].[4] S. Bravyi and A. Kitaev, Physical Review A , 022316 (2005).[5] A. G. Fowler, M. Mariantoni, J. M. Martinis, and A. N. Cleland, Physical Review A , 032324 (2012).[6] M. Levin and X.-G. Wen, Physical Review B , 045110 (2005).[7] H. Bombin and M. A. Martin-Delgado, Physical Review Letters , 180501 (2006), arXiv:0605138 [quant-ph].[8] H. Bombin, in Topological Codes , edited by D. A. Lidar and T. A. Brun (Cambridge University Press, 2013) arXiv:1311.0277.[9] B. Zeng, A. Cross, and I. L. Chuang, IEEE Transactions on Information Theory , 6272 (2011), arXiv:0706.1382.[10] B. Eastin and E. Knill, Physical review letters (2009).[11] H. Bombin and M. Martin-Delgado, Physical Review B , 075103 (2007), arXiv:0607736 [cond-mat].[12] H. Bombin, New Journal of Physics (2015), 10.1088/1367-2630/17/8/083002, arXiv:1311.0879.[13] A. Kubica, B. Yoshida, and F. Pastawski, New Journal of Physics , 083026 (2015).[14] H. Bomb´ın, Physical Review X , 031043 (2015).[15] B. J. Brown, N. H. Nickerson, and D. E. Browne, Nature Communications , 4 (2015), arXiv:1503.08217.[16] A. Kubica and M. E. Beverland, Physical Review A , 032330 (2015).[17] Y. Ozeki and N. Ito, Journal of Physics A: Mathematical and General , 5451 (1998).[18] M. Hasenbusch, F. P. Toldin, A. Pelissetto, and E. Vicari, Physical Review B , 184202 (2007).[19] T. Ohno, G. Arakawa, I. Ichinose, and T. Matsui, Nuclear Physics B , 462 (2004), arXiv:0401101 [quant-ph].[20] A. Honecker, M. Picco, and P. Pujol, Phys. Rev. Lett. , 047201 (2001), arXiv:0010143 [cond-mat].[21] C. Wang, J. Harrington, and J. Preskill, “Confinement-Higgs transition in a disordered gauge theory and the accuracythreshold for quantum memory,” (2003), arXiv:0207088 [quant-ph].[22] A. A. Kovalev and L. P. Pryadko, , 15 (2013), arXiv:1311.7688.[23] H. G. Katzgraber, H. Bombin, and M. A. Martin-Delgado, Physical Review Letters , 090501 (2009).[24] H. Bombin, R. S. Andrist, M. Ohzeki, H. G. Katzgraber, and M. A. Martin-Delgado, Physical Review X , 021004 (2012).[25] R. S. Andrist, Understanding topological quantum error-correction codes using classical spin models , Ph.D. thesis, ETHZurich (2012).[26] K. Hukushima and K. Nemoto, Journal of the Physical Society of Japan , 1604 (1996), arXiv:9512035 [cond-mat].[27] The 3D subsystem color code has the gauge group G = (cid:104) (cid:81) ∆ ( L ) (cid:51) t ⊃ e X ( t ) , (cid:81) ∆ ( L ) (cid:51) t ⊃ e Z ( t ) | e ∈ ∆ ( L ) (cid:105) and the stabilizergroup S = (cid:104) (cid:81) ∆ ( L ) (cid:51) t ⊃ v X ( t ) , (cid:81) ∆ ( L ) (cid:51) t ⊃ v Z ( t ) | v ∈ ∆ ( L ) (cid:105) .[28] D. Gottesman, Physical Review A , 1862 (1996), arXiv:9604038 [quant-ph].[29] A. Calderbank, E. Rains, P. Shor, and N. Sloane, Physical Review Letters , 405 (1997), arXiv:9605005 [quant-ph].[30] For correction of Z -errors, exchange X ↔ Z and redefine B = ∆ ( L ), B = ∆ ( L ), B = ∆ ( L ), ∂ e = (cid:80) ∆ ( L ) (cid:51) t ⊃ e t and ∂ t = (cid:80) ∆ ( L ) (cid:51) v ⊂ t v for any e ∈ ∆ ( L ) and t ∈ ∆ ( L ).[31] M. H. Freedman and D. A. Meyer, Found. Comp. Math. , 325 (2001), arXiv:9810055 [quant-ph].[32] N. Delfosse, Physical Review A , 012317 (2014), arXiv:1308.6207.[33] F. J. Wegner, Journal of Mathematical Physics , 2259 (1971).[34] H. Nishimori, Progress Of Theoretical Physics , 1169 (1981).[35] M. Palassini and S. Caracciolo, Physical Review Letters , 5128 (1999).[36] N. Goldenfeld, Lectures on phase transitions and critical phenomena (Addison-Wesley Publishing Company, 1992).[37] S. Elitzur, Physical Review D , 3978 (1975), arXiv:0305013 [hep-lat].[38] K. G. Wilson, Physical Review D , 2445 (1974).[39] J. Kogut, Reviews of Modern Physics (1979).[40] M. Creutz, L. Jacobs, and C. Rebbi, Physical Review D , 1915 (1979).[41] Note that this threshold was estimated for the code on a different 3D lattice and thus the comparison should only beconsidered as a guide to intuition.[42] K. Binder and A. P. Young, Reviews of Modern Physics , 801 (1986).[43] M. S. S. Challa, D. P. Landau, and K. Binder, Physical Review B , 1841 (1986).[44] K. Binder, Rep. Prog. Phys. , 783 (1987).[45] J. Lee and J. M. Kosterlitz, Physical Review B , 3265 (1991).[46] C. Borgs and R. Kotecky, Physical Review Letters , 1734 (1992).[47] H. G. Katzgraber, S. Trebst, D. A. Huse, and M. Troyer, Journal of Statistical Mechanics: Theory and Experiment2006