Surface codes, quantum circuits, and entanglement phases

Surface codes—leading candidates for quantum error correction (QEC)—and entanglement phases—a key notion for many-body quantum dynamics—have heretofore been unrelated. Here, we establish a link between the two. We map two-dimensional (2D) surface codes under a class of incoherent or coherent errors (bit flips or uniaxial rotations) to (1 + 1)D free-fermion quantum circuits via Ising models. We show that the error-correcting phase implies a topologically nontrivial area law for the circuit’s 1D long-time state | Ψ ∞ ⟩ . Above the error threshold, we find a topologically trivial area law for incoherent errors and logarithmic entanglement in the coherent case. In establishing our results, we formulate 1D parent Hamiltonians for | Ψ ∞ ⟩ via linking Ising models and 2D scattering networks, the latter displaying respective insulating and metallic phases and setting the 1D fermion gap and topology via their localization length and topological invariant. We expect our results to generalize to a duality between the error-correcting phase of ( d + 1)D topological codes and d -dimensional area laws; this can facilitate assessing code performance under various errors. The approach of combining Ising models, scattering networks, and parent Hamiltonians can be generalized to other fermionic circuits and may be of independent interest.


I. INTRODUCTION
The entanglement of quantum states characterizes many-body phases.For example, zero modes appear in the entanglement spectrum of topologically ordered [1,2] and symmetry-protected topological phases [3][4][5], including topological insulators and superconductors [6].The entanglement entropy is a dynamical probe: for example, starting from local product states, it grows linearly in generic many-body systems, but only logarithmically in many-body localized phases [7,8].Similarly, ground states of gapped local Hamiltonians [9,10], short-range correlated states [11][12][13], and almost all many-body localized eigenstates [14] exhibit an area law, i.e., the entanglement entropy grows with a subsystem's area, while generic (random) states follow a volume law [15].
The entanglement entropy and entanglement spectrum can also characterize purely dynamical phases without an underlying Hamiltonian.While the long-time evolution of a density matrix with a unitary random circuit will generally yield a volume-law [16,17], non-unitary elements change this picture: When following the quantum trajectory behind a density matrix, i.e., post-selecting measurement outcomes, hybrid circuits that consist of unitary gates and measurements exhibit a transition between area-law and volume-law phases as a function of measurement rate [18][19][20][21].This transition can also occur in measurement-only dynamics [22][23][24][25], and similar area-law to logarithmic-law transitions occur for weak measurements [26][27][28].Due to the post-selection, directly studying these phases experimentally requires a number of runs that is exponential in the system size and circuit depth [29,30].This difficulty may however be overcome via local probes of the entanglement transition [29], rotating space and time directions in the circuits [30,31], or correlating with classical simulations [32][33][34][35].
In this work, we show that entanglement features also usefully characterize quantum error correction (QEC) [36][37][38].Specifically, we establish a link between entanglement phases in hybrid quantum circuits that we derive from the surface code [39][40][41][42][43], and the phases of QEC in the surface code, the latter being a leading candidate for QEC with recent proof-of-principle experiments [44,45].Our results are summarized in Fig. 1.
The link we describe is distinct from recent entanglement-QEC relations via the scrambling of quantum information in hybrid circuits [32,[46][47][48][49][50][51][52][53].There, the counter-intuitive robustness of the volume-law phase against a small but non-zero rate of local measurements is explained by this phase supporting emergent QEC code spaces generated by the scrambling dynamics [46][47][48].While this may have potential quantum applications, its uses for fault-tolerant quantum computing have both practical and fundamental limitations [54,55].By focusing on the surface code, the links we describe pertain to codes and error models that are explicitly defined (instead of being emergent), and relate entanglement phases and practically motivated schemes for QEC.Our work is also conceptually distinct from topological order emerging in 2+1-dimensional quantum circuits with surfacecode ingredients [22,56]: Both our concept of interest (i.e., the phases of surface-code QEC) and the quantum circuits this relates to (1+1 instead of 2+1-dimensional circuits as we next note) are different.
Our main result is to map the two-dimensional (2D) surface code under incoherent X errors (i.e., bit-flips) or coherent exp[iϕX] errors (with angle ϕ) to (1 + 1)D free-fermionic hybrid quantum circuits and to embed the phases of QEC in the entanglement phases of these circuits' long-time 1D states.This links QEC phases to entanglement phases.Interestingly, unlike the volumelaw-QEC relation in scrambling, we find that the errorcorrecting phase (QEC phase for short) maps to a 1D area law, namely a topologically nontrivial 1D phase.
The overture to establishing these links is a mapping from the 2D surface code to 2D random-bond Ising mod-els (RBIMs) [42,57].This opens a direct route to (1+1)D dynamics upon viewing the RBIM transfer matrix [58][59][60] as a quantum circuit.While for bit flips, yielding real Ising couplings [42], this is the familiar (d + 1)D classical to d-dimensional quantum duality, additional considerations are needed for the coherent case where Ising couplings are complex [42,57].This is provided by a further mapping between 2D Ising models and 2D scattering networks [57,59,60].
The entanglement phases of QEC, and the broader entanglement (and Ising) phases they are embedded in, are sketched in Fig. 1: For bit-flips [real Ising couplings, panel (a)], we find area-law phases both below and above the QEC threshold.These phases correspond to an insulating network (see also Refs. 42,57,60,and 61).The nontrivial long-time-state topology below threshold is signaled by a zero mode in the entanglement spectrum and by the (interrelated) topological invariants for the 1D state and the 2D network [57,60].For coherent errors [complex Ising couplings, panel (b)], the QEC phase corresponds to the same entanglement phase (and networkmodel phase, cf.Ref. 57) as the incoherent QEC phase.Above threshold, however, we find a phase with entanglement entropy increasing logarithmically with system size.Here the network is metallic (see also Ref. 57).
While these results build on two existing links, from surface code QEC to Ising and network models on the one hand [42,57,60] and between network models and entanglement phases on the other [61], they establish a conceptually novel link between surface code QEC and entanglement phases in transfer matrix space, a connection we expect to exemplify a broader correspondence with intriguing implications [62].Firstly, our construction naturally generalizes to dualities between (d + 1)D codes with a local structure (such as topological codes [38,63]) and d-dimensional entanglement phases in (d+1)D quantum circuits.Secondly, by having found it to emerge for qualitatively different error types (bit flips and coherent errors), we expect a general correspondence between the QEC phase and transfer matrix area laws.This opens the door to using the area laws' classical simulability to chart the QEC phase for various codes and errors, including settings with coherent errors where a free-fermion description is unavailable, thus tackling a key challenge in QEC [64][65][66][67][68][69].Furthermore, by mapping error-corrupted codes to the entanglement structure arising from the long-time (i.e., infrared) dynamics of a system one dimension lower, our results anticipate a deep connection to the characterization of error-corrupted topologically ordered states via boundary phases and their transitions [70][71][72][73].
On a technical level, our use of Ising and network models together allows us to analytically establish the correspondence between Anderson insulating (i.e., disordered) networks and free-fermion area-laws, which we further characterize via quasi-local 1D parent Hamiltonians.These advances, extending clean-system analytic results and disorder numerics [61], may be of independent interest for free-fermion hybrid quantum circuits.

II. ISING MODELS FOR RANDOM BIT FLIPS AND FOR COHERENT ERRORS
The link between surface-code QEC and hybrid quantum circuits will be through a generalized 2D RBIM [57][58][59][60] on the square lattice, with Hamiltonian and partition function Z = {σv} exp(−H) (the inverse temperature is absorbed into the couplings).The Ising spins σ v = ±1 have nearest-neighbor coupling constant J with random signs η vv ′ ; the latter are drawn from an uncorrelated random distribution where η vv ′ = 1 with probability 1 − p and η vv ′ = −1 with probability p.We consider two choices of J: either purely real or complex

A. Surface code basics
As we explain below, both choices for the couplings in Eq. ( 1) originate in surface-code QEC [42,57], cf.Fig. 2. We consider the 2D toric code [40,74] on the square lattice.This is a topological stabilizer code [38] with qubits on the lattice links and with two types of stabilizers that (in the bulk) each act on four neighboring qubits: X-stabilizers S X v = i∈v X i are assigned to vertices v and Z-stabilizers S Z w = i∈w Z i to plaquettes w, where X i and Z i are Pauli operators [40].The states |ψ⟩ that for all v and w satisfy S X v |ψ⟩ = |ψ⟩ and S Z w |ψ⟩ = |ψ⟩ constitute the logical subspace.The number of states that satisfy these conditions depends on the boundary conditions [40].For concreteness, we focus on a cylinder geometry with boundary conditions yielding a two-dimensional computational space, i.e., one logical qubit.(Our considerations, however, are more general, cf.Ref. 57 for details on a planar geometry.)In particular we use "smooth boundaries" [42] so that one of the logical operators is X = i∈γ X i with the product being over qubits i on the shortest sequence γ of vertical links (in terms of Fig. 2) along the length of the cylinder, while the conjugate logical operator Z = i∈γ ′ Z i is the product of Z i along its circumference (with the shortest path γ ′ of vertical links).Equivalent logical operators arise from these upon stabilizer multiplication.
The toric code can correct X-errors and Z-errors independently; the considerations for the two are analogous.In what follows, we focus on X-errors, considering incoherent bit flips, i.e., X i being applied with probability p, or coherent errors from the application of exp[iϕX i ] on each qubit.(The latter arise from unwanted gate rotations-ubiquitous in quantum devices.)A string of X i being applied, whether from bit-flips or as a contribution from i exp[iϕX i ], can be detected by syndrome measurements: S Z w = −1 mark the end points of applied X i -strings.The set of S Z w eigenvalues is called the syndrome s.Given syndrome s, applying an X i string C s with the same end points returns the state to the computational space [42].While the end points are fixed, the strings C s themselves can vary: Applying S X v to C s = i∈s X i adds or removes loops of X i operators, and thus changes the strings contained in C s but not their end points.Furthermore, by S X v C s = C s S X v and S X v |ψ⟩ = |ψ⟩, this leaves C s |ψ⟩ invariant.Applying X, however, also leaves the end points invariant, but XC s |ψ⟩ ̸ = C s |ψ⟩.Hence there are two inequivalent classes (homology classes [42]) of error: those equivalent to C s and those to XC s .In QEC, given syndrome s, a decoder must decide which homology class the error is in and hence whether to apply C s or C s X to return the state back to the logical subspace.We denote both cases by C s Xq with q = 0, 1.

B. Ising mappings
We now relate surface code QEC to Eq. (1), starting with random bit flips (i.e., incoherent errors [75,76]).For this case, we follow Ref. 42.On each qubit j, a bit flip X j occurs with probability p; the qubit stays intact with probability 1 − p.Thus the probability of an X-string C s Xq occurring is the product over all qubits where η (Cs,q) j = −1 if X j occurs (i.e., contained in C s Xq ) and η (Cs,q) j = 1 otherwise.Henceforth we suppress the superscripts in η (Cs,q) j .As we noted above, a syndrome s does not determine a unique string C s , but only its end points.To obtain the probability P s,q that syndrome s occurs, and does so via an error in the homology class q of C s Xq , we thus need to sum over P C ′ s ,q for all other strings C ′ s Xq with the same s and q.Fixing a reference string C s Xq and multiplying it by X-stabilizers generates another such string, Xq , where n v ∈ {0, 1}.The set of all {n v } configurations generates all such homologically equivalent strings, i.e., all strings given s and q.
On each vertex, we now introduce Ising spins σ v = (−1) nv = ±1 (valued −1 when S X v is contained in the stabilizer product and 1 otherwise), cf.Fig. 2. Each qubit has two neighboring X-stabilizers.(This is evident in the bulk from Fig. 2; we also use cylinder termination with this property [57].)Each qubit thus corresponds to a bond between nearest neighbor σ v .When exactly one X-stabilizer neighboring qubit j is contained in the stabilizer product, we must swap p ↔ 1 − p in the corresponding factor in P C ′ s ,q ; in order words, the exponent e ηj J → e −ηj J .We can express this swap via the σ v : For each qubit, i.e., Ising bond, we include the product of the two neighboring σ v and thus write the product over all qubits as one over all Ising bonds Here, we also relabeled η j → η vv ′ , using again that each qubit is located at bonds between nearest-neighbor σ v .Summing over all possible strings with a given s and q, or equivalently over all Ising spin configurations {σ v }, we obtain P s,q = ( j p(1 − p))Z s,q , where Z s,q = {σv} exp(−H s,q ) is the RBIM partition function.Here, H s,q defined in Eq. ( 1) with J = (1/2) log((1 − p)/p); the subscripts s, q denote the reference string C s Xq that sets the configuration {η vv ′ } for H s,q [implicit in Eq. ( 1)].Since random bit flips occur with probability p, the signs η vv ′ = −1 with probability p and η vv ′ = 1 with probability 1 − p.This choice of p and J defines the Nishimori line [42,77].Our exploration of entanglement phases includes both this line, but we are also interested in the phase diagram in the broader p − J space [cf.Fig. 1(a)].
We next review, following Ref.57, the Ising mapping for coherent errors of the form U = i U i with U i = exp(iϕX i ) [78].The probability P s,q now arises from an overlap, P s,q = |⟨ψ|C s Xq U |ψ⟩| 2 .Here, we take |ψ⟩ to be the +1 eigenstate of Z so that |ψ⟩ and X|ψ⟩ are orthogonal and hence P s,q are probabilities.(This P s,q is also related a QEC fidelity under a suitable Bloch-sphere average [57,79].)The amplitude ⟨ψ|C s Xq U |ψ⟩ can be evaluated similarly to how P s,q was in the incoherent case, but now instead of a sum over the probabilities P C ′ s ,q of various X-strings, the expansion of U involves their coherent

FIG. 2. (a)
A bulk patch of the toric code, with physical qubits marked by black dots, while X-and Z-stabilizers by white and gray disks, respectively.(b) The code is mapped to a RBIM with real couplings for incoherent and complex couplings for coherent errors; X-stabilizers map to Ising spins σv.The nearest-neighbor couplings have sign η vv ′ .sum.To get the amplitude, we must thus replace p → i sin(ϕ) and 1 − p → cos(ϕ) in our previous derivation.As a result, P s,q = ( j | sin ϕ cos ϕ|)|Z s,q | 2 with J = −(1/2) log(i tan ϕ).As before, Z s,q = {σv} exp(−H s,q ) with the H of Eq. ( 1).
For the coherent-error QEC problem, now ϕ sets the syndrome distribution and hence η vv ′ , in a coherent generalization of the Nishimori line [57].Sampling η vv ′ according to this is more difficult than for bit flips: instead of sampling independently for each qubit (i.e., Ising bond), one must now sample bonds in certain sequence [67,79] to sample from P s = P s,0 + P s,1 [67,79].While this is needed for quantitative accuracy (e.g., for the error threshold or for critical properties), here we use a simplified model where we draw η vv ′ from an uncorrelated distribution with η vv ′ = −1 occurring with probability p and η vv ′ = 1 with probability 1 − p.This model thus has a ϕ − p phase diagram [Fig.1(b)].
Taking p = sin 2 ϕ [shown dashed in Fig. 1(b)] in this ϕ − p space mimics the QEC problem in a manner reminiscent of the Pauli twirl approximation [80,81] which replaces each U i by a bit flip occurring with probability p = sin 2 ϕ.Pauli twirling would however make this replacement from the outset, yielding the incoherent RBIM at this p, in contrast to using p = sin 2 ϕ with the complex RBIM.The latter goes qualitatively beyond Pauli twirling: the partition function (and hence the quantum circuit below) accounts for the coherent sum over X-string amplitudes-the key feature distinguishing coherent from incoherent errors.For this reason, we call p = sin 2 ϕ in the complex RBIM a "partial Pauli twirl".
Along the p = sin 2 ϕ partial Pauli twirl line we expect the qualitative structure of the QEC phase diagram to be the same for the full coherent error model and our simplified model; we shall further substantiate this expectation in Sec.IV using the scattering network description.

III. QUANTUM CIRCUIT
To relate the Ising models to quantum circuits, we express the partition function using the transfer matrix [58,60].Following standard steps [58][59][60][82][83][84] agnostic to whether couplings are real or complex, the partition function for a system on a cylinder of length L and circumference M is Z = ⟨α L |M|α 0 ⟩ where |α r ⟩ encodes boundary conditions at the x = 0, L ends of the cylinder and M is the transfer matrix where the hat distinguishes this many-body operator from its single-particle counterpart M in Sec.IV.The two kinds of transfer matrix layers are with where the labels h and v distinguish horizontal (η n,i ).The Pauli X i and Z i in Eq. ( 5) act in an M -site 1D transfer matrix space.The V i and H i involve complementary terms from the transverse field Ising model.The layers thus commute with the Z 2 symmetry P = j X j shared with this model.
Since the individual terms in the exponentials in Eqs.(5) mutually commute, we write ).The transfer matrix thus consists of a successive application of layers of onebody and two-body gates: it is a quantum circuit.The gates are not unitary, but depending on whether J is real or complex, they can yield both real and imaginary time evolution (see below and Fig. 3).We define the entanglement phases of M as of the long-time state |Ψ ∞ ⟩ obtained by time-evolution with M, starting a generic definite-parity initial state |Ψ 0 ⟩ (not necessarily |α 0 ⟩).
It will be beneficial to explore this in a fermionic setting.This will allow us to show that M is (essentially) a free-fermion circuit and that |Ψ 0 ⟩ can be taken as a fermionic Gaussian state, without loss of generality.(Gaussian states are ground states or thermal states of free-fermion Hamiltonians [6,85].)To construct a fermionic quantum circuit from M, we switch to a Majorana basis [60] via a Jordan-Wigner transformation [86].The Majorana fermions γ † i = γ i (i = 1 . . .2M ) satisfy the canonical anticommutation relations {γ i , γ j } = 2δ ij [87], and allow one to express the parity as n,i < 0, the Hn,i gates involve unitary doublebraids (blue) and imaginary time evolution (dark red).(b) For complex Ising couplings J = −(1/2) log(i tan ϕ), from coherent errors, the gates Hn,i are always unitary (blue), while the Vn,i consist of unitary braids (blue) and imaginary time evolution (dark red).
The appearance of P in V n,M is due to the nonlocality of the Jordan-Wigner transformation; it arises from describing a bosonic, i.e., qubit-based, quantum circuit with fermions.As swapping P swaps the sign of a fermionic hopping around the cylinder, changing the fermion parity changes between periodic and antiperiodic boundary conditions (pbc and apbc, respectively) on fermions.[Note that an η n,i string corresponding to X achieves the same, so changing q also changes between pbc and apbc.]While retaining P in V n,M , and hence the intertwined parity and fermion boundary conditions, is important for establishing surface code QEC features from the circuit [57], it is less crucial for establishing the circuit's entanglement phases, provided we consider both pbc and apbc and both parities for fermions (see Secs.III A and III B).In this way, we can also use V n,M = exp (−iκ n,M γ 2M γ 1 ).Henceforth we call circuits with this V n,M "purely fermionic", to distinguish from the fermionized transfer matrix M (henceforth called "bosonic").
These quantum circuits are different from previous fermionic mappings of the toric code [67,88]: these involve Abrikosov pseudo-fermions [89] that have a parity constraint.While the pseudo-fermion mapping can be used to sample from the syndrome probabilities in the coherent case [67], the circuit that arises from a statisticalmechanics mapping computes the probabilities for the different homology classes, and incorporates both incoherent and coherent errors on a unified footing.
We now discuss how the circuit combines real and imaginary time evolution, cf.Fig. 3. Real couplings J (from incoherent bit flips) correspond to purely imaginary time evolution up to double braiding of Majoranas: The V n,i are matrix exponentials of Hermitian operators, but for the H n,i this is true only when η where γ 2i−1 γ 2i is the unitary double braiding of Majoranas, cf.Fig. 3(a) [90][91][92].Complex couplings J = −(1/2) log(i tan ϕ) (from coherent errors) correspond to a mixed real-and imaginary time evolution: Here, the H n,i are unitary operators since κn,i = i(ϕ 2 and imaginary time evolution exp(−i Re κ n,i γ 2i γ 2i+1 ), as illustrated in Fig. 3(b).The surface code, in particular with coherent errors, thus provides a concrete physical motivation for fermionic quantum circuits alternating real and imaginary time evolution, studied in relation to emergent conformal symmetries [27] and classifications of fermionic quantum circuits and tensor networks [61].

A. Final states as 1D ground states
We now further specify the settings for defining the entanglement phases of M. We consider the properties of a long-time state |Ψ ∞ ⟩.That is, we consider the large L limit of the evolution where normalizing |Ψ L ⟩ (as required by the evolution not being unitary) is left implicit.As M is nonunitary, a useful view on its features can be obtained from its singular value decomposition.We write [57,60] where the left singular vectors |φ n ⟩are the eigenvectors of M M † and the right singular vectors | φ n ⟩ are the eigenvectors of M † M. The energies E n can be interpreted as those of a 1D Hamiltonian H, defined by M M † = e −LH , which has |φ n ⟩ as its eigenvectors.
We next define the large L limit more carefully: Considering that M is parity conserving, and denoting by δε the gap between the lowest and second-to-lowest energies of H eigenstates with the same parity as that of |Ψ 0 ⟩, we define the large L limit by Lδε ≫ 1. (The energy levels of H become increasingly non-random upon increasing L, cf.Sec.V A.) For |Ψ ∞ ⟩, this implies hence |Ψ ∞ ⟩ ∝ |φ min ⟩, the lowest-energy state (with energy E min ) of H with the same parity as that of |Ψ 0 ⟩.This is the ground state of H (or a ground state if there is a degenerate ground space) only if a ground state exists with this parity.This distinction is important when H is gapped.In this case, we consider states |Ψ 0 ⟩ with each parity and, depending on whether we deal with the fermionized bosonic M or its purely fermionic version, we also consider both pbc and apbc, i.e., q = 0, For the characterization of |Ψ ∞ ⟩, a further key feature is that the gates V n,i and H n,i are quadratic, and hence H is a 1D free fermion Hamiltonian.(This holds as is for the purely fermionic version of M; for the bosonic M it holds for each parity.)This implies that |Ψ ∞ ⟩ is a free-fermion state for any definite-parity initial state |Ψ 0 ⟩.Viewing |Ψ ∞ ⟩ as a free-fermion ground state is particularly useful in establishing its topological and entanglement features.
To establish a topological characterization, we will use that gapped free-fermion Hamiltonians in 1D are distinguished by the response of their ground-state fermion parity to a change between pbc and apbc.Specifically, in a topologically nontrivial system, the respective groundstate parities satisfy P pbc GS = −P apbc GS [87].In a topologically trivial system we have P pbc GS = P apbc GS .This allows one to define a topological invariant I = P pbc GS P apbc GS with I = −1 in a topological phase [87].
The topological aspects and boundary conditions are thus intertwined.In particular, while the ground state is always unique for the purely fermionic gapped H, subtleties arise in the bosonic problem when I = −1.This is because in this problem, changing P changes the fermion boundary conditions and for I = −1, changing these boundary conditions changes P GS , where P GS is the ground-state parity of the purely fermionic H. Hence, either P = P GS for both P = ±1 or for neither.That is, when the purely fermionic H is gapped and has I = −1, the state |φ min ⟩ in Eq. ( 8) for the bosonic problem is either a ground state of this purely fermionic H for both P = ±1, or it is its lowest excited state for both P = ±1.The two-fold ground space degeneracy in the former case, of course, just corresponds to spontaneous symmetry breaking in the spin-chain, generalizing that in the transverse-field Ising chain.
In this I = −1 case, one can switch between |φ min ⟩ being a ground or excited fermionic state by switching P GS (without changing P ).This is achieved by changing boundary conditions via changing q, i.e., changing η along the cylinder, corresponding to the application of X. (For I = 1, both q-values work because P pbc GS = P apbc GS .)The above considerations highlight that, depending on whether we use a purely fermionic or the bosonic form of the circuit M, exploring the entanglement phases requires considering both parities and boundary conditions.(In practice, choices exist that work for most disorder realizations.For example, q = 0 for the bosonic M works because I = −1, as we will show, arises for small p where long η (v) n,i = −1 chains, effecting a spurious boundary con-dition change to be undone by q = 1, have probability exponentially suppressed in L.)

C. Characterizing free-fermion entanglement
The quantum circuit having quadratic gates also enables both the single-particle characterization and the efficient numerical evaluation of entanglement properties.Using that |Ψ ∞ ⟩ is the same Gaussian state regardless of the details of the definite-parity initial state |Ψ 0 ⟩, we can choose |Ψ 0 ⟩ to be Gaussian as well.We can then use that any Gaussian state evolved by H n and V n remains Gaussian [85], with the same parity as that of |Ψ 0 ⟩.This implies that all many-body quantities can be computed using fermionic linear optics [85].The central object of this approach is the correlation matrix from which all higher correlators follow [85].Following Bravyi, we evolve the matrix directly [85] (with n denoting the time step) instead of the (exponentially large) density matrix ρ To calculate the entanglement entropy and entanglement spectrum for a subsystem A, the correlation matrix C(n) A of the corresponding reduced density matrix ρ(n) A (also a Gaussian state [6,85]) can be obtained from C (n) by keeping only those indices contained in A.

Since C(n)
A is real and antisymmetric, it can be block-diagonalized via a Youla decomposition [93] C(n) A = QΣQ T where Σ = diag({iλ r Y }) and Y is the 2 nd Pauli matrix.The set of λ r is the single-particle entanglement spectrum [6]; the matrix i

C(n)
A is the single-particle "entanglement Hamiltonian".It determines the entanglement entropy as (10) We will from now on consider the entanglement spectrum and entropy only for bipartitions of the system into two halves of size M/2, and denote the entanglement entropy by S M/2 .As we shall see in the Sec.VI, the entanglement spectrum and the entanglement entropy are key characteristics of the circuit M and hence also characterize surface-code QEC.

IV. NETWORK MODEL
We now turn to the network model (cf.Fig. 4).For a many-body operator T = exp(i ij γ i q ij γ j ), single Majorana operators transform as [87] T γ i T −1 = j t ji γ j , t = exp(4iq).At each Ising bond in Fig. 2(b) we now have junction matrices Hn,i and Vn,i.They scatter directed Majorana modes residing on the networks' links.Translating Hn,i and Vn,i into junction scattering matrices requires different link direction layouts for real and complex Ising couplings: counter-propagating pairs of links for real couplings [60] and co-propagating pairs for complex Ising couplings [57].
We can thus switch to single-particle matrices h n,i = exp(2κ ′ n,i Y ) and v n,i = exp(2κ n,i Y ) instead of the respective many-body operators H n,i and V n,i [61], where the 2 × 2 Pauli matrix Y acts on the (2i − 1, 2i) th (for h n,i ) and (2i, 2i + 1) th degrees of freedom (for v n,i ).We denote the resulting 2M × 2M transfer matrix by M.
For real J, the single-particle operators are pseudounitary, Zt −1 Z = t † with t = v n,i or t = h n,i , and can thus readily be interpreted as single-particle transfer matrices [60] that, when acting on a pair of counterpropagating modes c = (c n , c n+1 ), conserve their current c † Zc [94], cf.Fig. 4(a).In this way, each matrix t = v n,i or t = h n,i describes the scattering at a "junction", and the junctions form a scattering network.
For complex J, neither the v n,i nor the h n,i are pseudounitary [57].However, since κn,i = i(ϕ − (1 − η (h) n,i )π/4) is purely imaginary, h n,i is always unitary and we can interpret it as a scattering matrix connecting co-propagating modes.While v n,i is not unitary, the product Xv n,i is pseudounitary: That is, when acting on two counterpropagating modes, v n,i swaps them and conserves current [57], cf.Fig. 4(b).
The single-particle transfer matrices for both real and complex J imply real scattering matrices for each junction [57,60,82].This places the networks into Altland and Zirnbauer's symmetry class D [95].Hence, the links of the networks can be interpreted as describing directed 1D Majorana modes.(Upon taking the networks together with their time-reversed partners, the real and complex J cases realize two limits of the time-reversal symmetric class DIII Majorana network in Ref. 96.) Networks in symmetry class D often include disorder in the form of randomly placed "vortices" [59,60,82,97,98].(A vortex is a point-defect such that a mode encircling it picks up an extra π phase.)In our case the disorder is via the η n,i and this indeed introduces vortices.In the incoherent case, as is well known from the RBIM [60,82,97], η n,i = −1 imprints a pair of vortices adjacent to junction n, i in one of the sublattices.In terms of the surface code, a vortex appears at the adjacent S Z w , i.e., where S Z w = −1 due to the bit flip represented by η n,i = −1.In the coherent case, η n,i = −1 has the same effect, however, the manner in which a vortex can be encircled is different than in the incoherent case due to the propagation directions being laid out differently in the coherent errors' scattering network.
The networks, together with the vortex distribution, determine the phase of QEC [57].Using this fact, we can further substantiate why the p = sin 2 (ϕ) partial Pauli twirl line of our simplified ϕ − p model is expected to capture the qualitative phase structure of the full coherent error model: The network itself is the same for the two models since, apart from the bond signs, they originate from the same complex-J Ising model.(As we noted in Sec.II B, this captures the coherent summation over X-strings, the key feature of the coherent error problem.)By sampling η n,i differently, the two models differ in their vortex distribution.However, for both models, the rarity of η n,i = −1 for small ϕ implies tightly bound vortex pairs, while for sufficiently large ϕ vortices proliferate.These basic features dictate [97] that the qualitative phase structure along p = sin 2 (ϕ) is the same for the two models, albeit the quantitative details such as the phase boundary ϕ c or the critical properties may differ.
Our characterization of network models will include transport properties, specifically the dimensionless conductivity g = (L/M )⟨tr[T † T ]⟩ dis .Here T denotes the transmission matrix from the transmission-reflection grading of the total scattering matrix S = R T ′ T R ′ [94] and ⟨. . .⟩ dis denotes the disorder average.In an insulator, i.e., a localized network, the conductivity satisfies g ∝ e −2L/ξ where ξ is the localization length [98].A metallic network, in contrast, displays g ∝ ln(L).Both expressions hold in the large L limit, understood to be taken with fixed aspect ratio L/M .

V. ENTANGLEMENT PHASES VIA 2D ISING MODELS, NETWORKS, AND 1D FERMIONS
We next discuss how 2D Ising considerations combined with links between 2D scattering networks and 1D free-fermion parent Hamiltonians illuminate the entanglement phases of |Ψ ∞ ⟩.Our approach in this Section can be generalized to other fermionic quantum circuits, beyond our motivating surface-code problems, and hence may be of independent interest.The Ising model and parent Hamiltonian perspectives complement recent tensor-network-and scattering-network-based approaches [61] to entanglement phases in free-fermion circuits.In Sec.VI, we shall numerically confirm the insights we obtain here, returning our focus to the entanglement phases in the quantum circuits dual to the surface code with bit flips and coherent errors.

A. 2D networks and 1D Hamiltonians
We now link some features of 2D networks and of H from M M † = e −LH .We follow Ref. 57, where we noted that the links we describe bridge between the approach of Ref. 99 relating 1D and 2D topological phases via scattering matrices (the 1D Hamiltonians there, however, arise differently than here) and Ref. 60's pioneering insights linking topology in 2D networks and 1D systems.We focus on purely fermionic M.
The first key observation is that an insulating (i.e., localized) network implies that H = i 2 ij a ij γ i γ j is gapped.(Here we introduced the single-particle Hamiltonian ia with a real antisymmetric matrix a.) To see this, we note that Eq. ( 11) implies for the matrix M for M.This links the single-particle energies ε j ≥ 0 of ia to transport properties [60].In particular, one can show that the conductivity satisfies In an insulator, the g ∝ e −2L/ξ large-L asymptotics (with fixed L/M ) implies lim M →∞ ε 1 > 0 for the smallest energy ε 1 .(The energies ε j , and as such ε 1 , become increasingly non-random upon increasing the system size [60,94].)Hence, H is gapped, with gap lim M →∞ ε 1 = αξ −1 (with α > 0 order of unity accounting for the difference between average and typical ξ [98]).
In what follows, we refer to gapped H and insulating networks interchangeably.
The second key link, also implied by Eq. ( 12), is between the 1D topological invariant I and the reflection matrix R ′ of the 2D scattering network.Specifically, one can show that for a gapped H, i.e., an insulating network, we have . (Replacing R ′ by R gives the same result [99].)

B. Area law phases
We now show that when the purely fermionic H is gapped, i.e., the corresponding network is insulating, then |Ψ ∞ ⟩ satisfies the entanglement area law.If we knew that H is a local Hamiltonian, this would be an immediate consequence of its gap [9,10].However, from M M † = e −LH the locality is not obvious, even if L ∝ the circuit depth of M makes it plausible.To establish the area law, we will show that the correlations ⟨Ψ ∞ |γ a γ b |Ψ ∞ ⟩ decay exponentially with |a − b| (for ξ ≪ |a − b| ≪ M ); this provides a sufficient condition for |Ψ ∞ ⟩ to display an area law [11][12][13].Our approach does not assume the absence of disorder from η n,i ; in this way it complements the analytical arguments in Ref. [61] based on disorder-free networks.
We start by noting that for large L −iC where we take the trace in terms of the bosonic problem (i.e., use P -dependent boundary conditions) [100].This allows us to view ⟨Ψ ∞ |γ a γ b |Ψ ∞ ⟩ as an Ising correlator on the torus.This enables the use of space-time duality [30,31] to evaluate the correlation function.
The corresponding Ising model is defined by M † M; it consists of two coupled Ising patches, one for M † and one for M.These two Ising patches are in the same phase: the Hamiltonian for M and M † have identical spectra, so both of them are gapped; the corresponding phases are labeled by I for |Ψ ∞ ⟩.By I being defined by the fermion parity, and the parity being the same for |φ min ⟩ and | φ min ⟩, the value of I for M is the same as for M † .
The correlation function is thus that of γ a γ b embedded in the bulk in the transfer matrix of this 2D Ising model.To interpret this in the Ising language, we take a < b without loss of generality, and implement γ a γ b = γ a γ a+1 γ a+1 γ a+2 . . .γ b−1 γ b , up to an overall phase, by κ → κ + iπ/2, κ → κ + iπ/2 in the last layers of M, while leaving M † unchanged.This introduces Jη n,i along the line from a to b.In the 2D Ising language, the former yields σ a σ b , while the latter yields a seam of flipped horizontal bonds from a to b.The corresponding correlator is that of products of Ising spins and disorder operators: an Ising fermion correlator [101].This decays exponentially for both I = ±1 due to either the disorder or the Ising correlators decaying exponentially while the other being constant [57,60,101].Using space-time duality to orient the fermion string for γ a γ b along the temporal direction, one can show that C (∞) ab ∝ e −α|a−b|/2ξ , with ξ the localization length in the scattering network.This holds both typically and on average because the Ising model, or network, for M † M has η = −1 strings appear in pairs thus the rare long η = −1 strings (cf.Sec.III B) in the dual-temporal direction are inoperative.
This establishes the I = ±1 gapped phases of M, and the respective insulating phases of the scattering networks, as yielding an area-law |Ψ ∞ ⟩.The exponentially decaying correlations also imply that iC (∞) is a quasilocal (i.e., with couplings exponentially decaying with distance) single-particle Hamiltonian; it has eigenvalues ±1 and hence defines a gapped quasilocal parent Hamiltonian for |Ψ ∞ ⟩ [6,102,103].This results in the following signatures for the single-particle entanglement spectrum [5,6]: For both I = ±1, the entanglement Hamiltonian i C(∞) M/2 has a bulk "entanglement gap".When I = 1, the entire single-particle entanglement spectrum is gapped.When I = −1, however, the nontrivial topology implies entanglement zero modes (analogous to Majorana end states at physical boundaries).For finite M , the zero modes are split, yielding an entanglement energy level λ 0 satisfying λ 0 ∝ e −M/c with c > 0 increasing with ξ.
A gapless H can also arise; this happens if the network is metallic.While this is ruled out for an Ising model with real couplings [82], a metallic phase is generically part of the phase diagram when the couplings are complex [97].In this case, from |Ψ ∞ ⟩ being the ground state of a gapless 1D H, by analogy to the logarithmic S M/2 at criticality [18][19][20][104][105][106][107][108] we expect S M/2 ∝ ln M , i.e., a logarithmic entanglement phase.(See also Ref. 61 for linking metallic networks to logarithmic entanglement phases.) A prediction on scaling beyond these asymptotics can also be made if we note that the physics of metallic 2D networks is described by a nonlinear σ model [82,98].(We numerically verify this in Sec.VI for the coherenterror network.)In this model, the conductivity g is the only coupling; as a consequence, it follows singleparameter scaling g(L; p, J) = g[L/ℓ p,J ] tracing out the renormalization group flow of g. (Here ℓ p,J is an effective length scale.)This is a characteristic feature of the metallic phase that holds beyond the asymptotic g ∝ ln L regime.(The key requirement is diffusive transport, setting in for L much larger than the mean free path, i.e., the short-distance cutoff for the nonlinear σ model.)Based on this, we similarly expect single-parameter scaling for the entanglement, M/2 (p, J) = S(M/m p,J ), providing an fingerprint of the nonlinear σ model.

VI. NUMERICAL RESULTS
We now return to the link between the phases of QEC in the surface code and the entanglement phases of their dual quantum circuits.

A. Real Ising couplings
For the real-coupling RBIM, and thus the surface code with random bit flips, the transport properties of the network model have been extensively discussed in the literature [60,97,98,109,110].Here, we highlight one key observation: The phases on both sides of the transition are insulating, but characterized by different topological invariants [57,60]: we have I = −1 in the ordered Ising phase (including the error-correcting part of the Nishimori line) and I = 1 otherwise, as shown Fig. 1(a).
Turning to entanglement, in Fig. 5 we show the entanglement spectrum and entropy S M/2 along the Nishimori line (i.e., for surface-code QEC).Our initial state |Ψ 0 ⟩ is a random half-filled state [defined in terms of fermions c j = (γ 2j−1 +iγ 2j )/2], which we evolve for long cylinders, L = 5M .We find that the entanglement spectrum and entropy converge, indicating that |Ψ ∞ ⟩ has been reached.
In the entanglement spectrum, we observe the following features: Below the error threshold, p c ≈ 0.1093 [42,

60]
, where I = −1, the single-particle entanglement spectrum is gapped and has a zero mode whose energy decays exponentially with system size; the many-body entanglement spectrum is thus degenerate in the infinite-system limit.These features confirm expectations from Sec. V B for a topologically nontrivial phase.
The smallest entanglement eigenvalue λ 1 of the bulk is minimal close to the transition.With increasing M , the p where λ 1 is minimal shifts towards p c and the minimum itself λ 1 | min decreases as a power law with M , consistent with a critical phase at the transition [110].
On both sides of the transition, the entanglement entropy scales as an area law, i.e., it does not increase with the system width M .This again confirms expectations from Sec. V B for |Ψ ∞ ⟩ associated to insulating networks.
For I = −1, the entanglement entropy is bound from below by log 2, which reflects the presence of a zero mode.For I = 1, the entropy S M/2 goes to zero for large M and sufficiently large p. Near the transition, the entanglement entropy grows with M .Consistently with the area law away from p c , this is expected saturate unless p = p c .This is consistent with the p for which S M/2 is maximal shifting towards p c with increasing M .

B. Complex Ising couplings
In Fig. 6, we show the conductivity g for the complex RBIM motivated by coherent errors, focusing on the p = sin 2 ϕ partial Pauli twirl line [shown dashed in Fig. 1(b)].To probe the bulk value of g we work with a wide cylinder, M = 5L [96,111].The results for different ϕ are shown with different colors.When rescaling the length to a dimensionless L/ℓ(ϕ) with an appropriately chosen function ℓ(ϕ), the conductivity data collapses onto one of two scaling curves, depending on ϕ. (For completeness, we show the unscaled data in Appendix B.) FIG. 6. Dimensionless conductivity g for complex couplings at p = sin 2 ϕ as a function of the rescaled system length L/ℓ(ϕ) for wide systems with M = 5L, averaged over 100-10 4 configurations of η; error bars (2×standard error) are imperceptible.For angles above a critical ϕc, the conductivity increases with L/ℓ(ϕ) (metallic phase; dashed black line shows g ∝ (1/π) log[L/ℓ(ϕ)]).Below the transition, it decreases exponentially to zero (localized phase, dashed gray line shows the exponential tail).
For angles ϕ > ϕ c , the system is metallic: g increases with L and for sufficiently large systems it approaches the universal class-D result [98] g ∝ (1/π) ln L (dashed black line).For ϕ < ϕ c the system is in an insulating phase: for large systems, g decreases exponentially with L (dashed gray line).In this phase, ℓ(ϕ) is the localization length; it diverges close to the transition.The metal-insulator transition occurs at ϕ c = (0.095 ± 0.005)π-note that this value is significantly smaller than the coherent error threshold ϕ th = (0.14 ± 0.005)π we found in Ref. 57 by sampling the syndromes according to P s instead of sampling each η independently as we do here.
In the insulating phase, we find I = −1, as in the ordered Ising insulator for real J. (Our results are also consistent with the I = −1 insulator for vortices sampled according to P s [57].)On leaving the asymptotics, the scaling curves we find are qualitatively similar to previous results for class-D metal-insulator transitions [57,111,112].Furthermore, the scaling in the metallic regime follows closely the nonlinear σ model renormalization group flow for g [98].This excellent agreement with nonlinear σ model predictions is in contrast to the results for η being sampled according to P s for coherent errors [57].
We now discuss the signatures of these phases in the entanglement spectrum and entropy.In Fig. 7, we show these quantities, continuing to focus on the partial Pauli twirl line p = sin 2 ϕ.We again start the evolution from a random half filled state and converge to |Ψ ∞ ⟩ using long cylinders with L = 5M .
In the insulating phase, the entanglement spectrum displays a zero mode and has a bulk gap.The entanglement entropy displays an area law and it slowly decreases with ϕ to the ϕ = 0 value S M/2 = log 2. The entanglement zero mode λ 0 decays exponentially with M , shown in the inset of Fig. 8 for various angles ϕ < ϕ c .These FIG. 8. Entanglement entropy S M/2 in the metallic phase (ϕ > ϕc) for complex couplings at p = sin 2 ϕ as a function of the rescaled circumference M/m(ϕ).We averaged over 2 5 -2 8 configurations of η; error bars are 2×standard error and the gray dashed line serves as a guide for the eye showing a logarithmic increase.The inset shows λ0, the exponentially decaying zero mode; the decay length increases with ϕ.
features agree with the behavior expected for an I = −1 insulator, i.e., a topological area-law |Ψ ∞ ⟩, cf.Sec.V B.
In the metallic phase, the entanglement spectrum gap decreases as a power law in M (with a ϕ-dependent power) and the entanglement entropy increases with M .The large-M asymptotic is S M/2 ∝ ln M (shown dashed in Fig. 8), indicating a logarithmic entanglement phase.(The data fit S M/2 ∝ ln 2 M , derived in a related context [113], similarly well.)Similarly to g, rescaling M → M/m(ϕ) by a ϕ-dependent length m(ϕ) collapses data points a smooth curve, shown in Fig. 8.This confirms the expectations from Sec. V C: S M/2 in the logarithmic entanglement phase inherits single-parameter scaling from g. [The function m(ϕ), however, does not equal ℓ(ϕ) used for g in Fig. 6.]

VII. CONCLUSION
In this work, we related the phases of surface-code QEC for coherent and incoherent errors to entanglement phases.In particular, using a mapping to a RBIM with real couplings for incoherent [42] and complex couplings for coherent errors [57], we could interpret the RBIM transfer matrix as a quantum circuit for mixed realimaginary time Gaussian evolution that converges to a long-time Gaussian state |Ψ ∞ ⟩ from generic (e.g., random) definite-parity initial states.
For both error types, the QEC phase is dual to phase where |Ψ ∞ ⟩ satisfies the entanglement area law.This phase is topologically nontrivial (I = −1), which implies that its gapped single-particle entanglement spectrum supports a zero mode.Consequently, the entanglement entropy is bounded from below by log 2. Above threshold, and for incoherent errors, we find an I = 1 area law.The state |Ψ ∞ ⟩ again has gapped entanglement spectrum but without a zero mode, and the entanglement entropy approaches zero away from the transition between the two area-law phases.For coherent errors, we find a logarithmic entanglement phase above the threshold.
The duality between QEC codes and entanglement phases provides a new perspective from which to study the dynamics of hybrid quantum circuits that is entirely distinct from previously considered emergent QEC in hybrid circuits [46][47][48].In particular, the surface code with coherent errors provides a natural physical system in terms of which to interpret hybrid dynamics alternating real-time and imaginary-time evolution and the associated transitions between area-law and logarithmic entanglement phases.In this sense, it is tempting to think of the logical error rate-a direct indicator of which phase of QEC the system is in-as an indirect fingerprint of entanglement phases and transitions, albeit via the dual system: the QEC code.
Our results not only show that such hybrid circuits can be motivated by QEC, but the entanglement phases also offer a novel characterization for the phases of QEC.The area law for the QEC phase is especially important in this regard [114].While we demonstrated this area law only for the specific error models we studied, previous results on more general incoherent errors suggest [62] that the entanglement entropy continues to exhibit an area law in the QEC phase for a broader class of errors.This suggests that the quantum circuits dual to the QEC problemwhich has more complicated statistical mechanics models for more general errors [115]-can be efficiently simulated in the QEC phase using matrix product states [116,117].Using this, and generalizing our approach to deriving statistical mechanics models for coherent errors, one may chart out the QEC phase for a broad class of errors, including the important open problem of coherent errors with generic SU(2) rotations [64][65][66][67][68][69].
Our analysis using scattering networks also offers new perspectives on the relations between such network models and entanglement [61].Our Ising considerations link insulating networks and area-law phases explicitly via correlations; this complements existing arguments [61] based on disorder-free networks.The link to quasilocal parent Hamiltonians and their topological invariants are, to our knowledge, also new aspects connecting scattering networks and entanglement phases [61].The entanglement gap and the presence or absence of entanglement zero modes emerge directly and naturally in this approach.The link between metals and logarithmic entanglement phases we find agrees with Ref. 61.To elucidate this link further, we showed that the entanglement entropy follows single-parameter scaling, similarly to the conductivity.The very good agreement we found for the latter with nonlinear σ model predictions suggests that a σ-model theory may be developed also for the entanglement entropy.(See Refs.113 and 118, that appeared independently of this work, for such σ-model theories.) Our results can also be viewed as pertaining to error-corrupted topological quantum memories.Unlike Refs.70-73, that appeared independently of this work, we focus on the state post stabilizer measurement, and encode all stabilizer information in a (1 + 1)D circuit.The circuit can be interpreted as the boundary theory of the error-corrupted state and the phase transition of this boundary theory expresses the loss of topologically encoded information in the bulk [70,71].
The quantum circuit duals for the surface code problems we study display some analogies to a family of Gaussian fermionic circuits studied recently [22,24,25,[119][120][121][122].It would be interesting to generalize our approach to construct network and Ising models for these circuits and thereby to characterize the "gapped" (area law) and "Goldstone" (logarithmic entanglement) phases found in their hybrid [119][120][121][122] and measurement-only [22,24,25] variants.Network models may shine light on the classification of these area-law phases, including explicitly establishing the topological origin of the log 2 entanglement entropy arising in one of these phases [24,25,122].

FIG. 1 .
FIG. 1. Phase diagrams of surface-code QEC, 1D entanglement, and 2D scattering networks.Bit-flip errors give real Ising couplings J [panel (a)]; coherent exp[iϕX] errors yield complex couplings J = −(1/2) log(i tan ϕ) [panel (b)].The RBIM bond signs, related to X error configurations, are swapped with probability p.The solid lines sketch the phase boundaries, the black dots mark computed phase transition points.The dashed line shows (a) the Nishimori line representing incoherent-error QEC and (b) the p = sin 2 ϕ "partial Pauli twirl" line for coherent errors .The error-correcting phase (QEC✓) is dual to a circuit yielding a topologically nontrivial 1D area law with entanglement-spectrum zero modes (ALtop), and to an insulating network with topological invariant I = −1 (I−1).Above threshold (QEC✗) we find (a) a topologically trivial 1D area law (no zero modes, ALtriv) and a I = 1 insulator (I1), or (b) a logarithmic entanglement phase (LL) and a metallic network.
FIG.3.Fermionic quantum circuits acting on Majorana fermion lines.The "time" direction, along the cylinder, is upwards.(a) Circuit for real Ising couplings J (from incoherent bit flips).If η (h)

FIG. 4 .
FIG.4.Network model for the (a) real and (b) complex Ising couplings from incoherent and coherent errors, respectively.At each Ising bond in Fig.2(b) we now have junction matrices Hn,i and Vn,i.They scatter directed Majorana modes residing on the networks' links.Translating Hn,i and Vn,i into junction scattering matrices requires different link direction layouts for real and complex Ising couplings: counter-propagating pairs of links for real couplings[60] and co-propagating pairs for complex Ising couplings[57].
data into the insulating regime [panel (b)] and metallic regime [panel (c)].Note panel (b) uses a log-scale (for better visibility of the exponential decay) and (c) a loglog scale (for better visibility of power laws at large ϕ).
exponential convergence due to the gap.)That is, when H is gapped, by the entanglement phases of M we mean those of this ground-state-converged |Ψ ∞ ⟩. (When H is gapless we do not need such qualification because |Ψ ∞ ⟩ is similar for either parity.) 1 (see Sec. III B).In this way, when H is gapped, we can take |Ψ ∞ ⟩ to be a ground state.(This is easily identifiable by the fast B. Gaps, ground states, boundary conditions