Symmetry-enforced many-body separability transitions

We study quantum many-body mixed states with a symmetry from the perspective of separability, i.e., whether a mixed state can be expressed as an ensemble of short-range entangled (SRE) symmetric pure states. We provide evidence for 'symmetry-enforced separability transitions' in a variety of states, where in one regime the mixed state is expressible as a convex sum of symmetric SRE pure states, while in the other regime, such a representation is not feasible. We first discuss Gibbs state of Hamiltonians that exhibit spontaneous breaking of a discrete symmetry, and argue that the associated thermal phase transition can be thought of as a symmetry-enforced separability transition. Next, we study cluster states in various dimensions subjected to local decoherence, and identify several distinct mixed-state phases and associated separability phase transitions, which also provides an alternate perspective on recently discussed 'average SPT order'. We also study decohered p+ip superconductors, and find that if the decoherence breaks the fermion parity explicitly, then the resulting mixed state can be expressed as a convex sum of non-chiral states, while a fermion-parity preserving decoherence results in a phase transition at a non-zero threshold that corresponds to spontaneous breaking of fermion parity. Finally, we briefly discuss systems that satisfy NLTS (no low-energy trivial state) property, such as the recently discovered good LDPC codes, and argue that the Gibbs state of such systems exhibits a temperature-tuned separability transition.

Suppose one has the ability to apply unitary gates that act in a geometrically local fashion on a many-body system.Starting from a product state, a specific circuit composed of such gates results in a specific pure state, and an ensemble of such circuits can therefore be associated with the mixed state ρ = i p i |ψ i ⟩⟨ψ i | where the pure state |ψ i ⟩ is prepared with probability p i .If one is limited to only constant depth unitary circuits, then the corresponding mixed state can be regarded as 'shortrange entangled' or 'trivial' [1,2], which generalizes the notion of short-range entangled pure state [3][4][5][6][7][8].In parallel with the notion of symmetry protected topological phases for pure states [9][10][11][12], it is then natural to define a trivial/short-ranged entangled symmetric mixed state (a 'sym-SRE' state) as one that can be obtained from an ensemble of pure states, where each element of the ensemble is prepared with only a constant depth circuit consisting of local, symmetric gates under some given symmetry.Motivated from experimental progress in controllable quantum devices where both unitary quantum dynamics and decoherence play an important role [13][14][15][16], in this paper we will explore mixed state phase diagrams where in one regime a mixed state is sym-SRE, and in the other regime, it is not.We will call such phase transitions 'symmetry-enforced separability transitions', since a sym-SRE state is essentially separable [1] (i.e. a convex sum of unentangled states) upto short-distance correlations generated by constant-depth unitaries.In the absence of any symmetry constraint, analogs of such transitions were recently studied in Ref. [17] in the context of decohered topologically ordered mixed states [18][19][20][21][22][23][24].
To make progress, we will try to leverage our understanding of the complexity of preparing pure many-body states using unitaries.Some of the questions that will motivate our discussion are: Do there exist separability phase transitions when pure-state symmetry protected topological (SPT) phases are subjected to decoherence, and if yes, what is the universality class of such transition?When a 2d chiral pure state (e.g. the ground state of an integer quantum Hall phase) is subjected to local decoherence, can the resulting density matrix be expressed as a convex sum of non-chiral states?Can the conventional, finite temperature phase transitions corresponding to the spontaneous breaking of a global symmetry be also thought of as separability transitions?
As an example, consider the transverse field Ising model on a square lattice.We provide an argument (Sec.III) that the Gibbs state for this model can be prepared using an ensemble of finite-depth local unitary circuit at all temperatures, including at T ≤ T c , where T c is the critical temperature for spontaneous symmetry breaking.It is crucial here that one is not imposing any symmetry constraint on the unitaries.This is consistent with previous works [25][26][27] where evidence was provided that the mixed-state entanglement corresponding to a Gibbs state that exhibits spontaneous symmetry breaking remains short-ranged at all non-zero temperatures, including at the finite temperature critical point (assuming absence of any co-existing finite-temperature topological order).However, if one only allows access to an ensemble of short-depth unitary circuits composed of Ising symmetric local gates, then using results of Ref. [21], we provide a rigorous argument that the Gibbs state can not be prepared for any T ≤ T c .We expect similar results to hold for other symmetry broken Gibbs states as well.Therefore, the conventional, finite-temperature symmetry breaking phase transition in a transverse-field Ising model can be thought of as a symmetry-enforced separability transition.This statement is true even when the transverse field is zero (i.e. for a classical Ising model) -the quantum mechanics still plays a role since the imposition of symmetry implies that one is forced to work with 'cat' (GHZ) states, which are long-range entangled.
In the context of pure states, a well-known example of symmetry enforced complexity is an SPT phase whose ground state can not be prepared using a finite depth circuit composed of symmetric local gates [9][10][11][12].Recent works have provided a detailed classification of SPT phases protected by zero-form symmetries that are being subjected to decoherence using spectral sequences and obstruction to a short-ranged entangled (SRE) purification [28,29].Progress has also been made in understanding non-trivial decohered SPTs using string operators [30] and 'strange correlators' [31,32], concepts that were originally introduced to characterize pure SPT states [10,33,34].Here we will be interested in understanding decohered SPT states from the viewpoint of separability, which, as we discuss in Sec.II, is a different notion of entanglement of mixed states than that based on SRE purification considered in Ref. [28,29].As hinted above, we define a sym-LRE (symmetric, longrange entangled) state as one which does not admit a decomposition as a convex sum of pure states which are all preparable via a finite-depth circuit made of symmetric local gates.If so, it is interesting to ask if there exist separability transitions between sym-LRE and sym-SRE states as a function of the decoherence rate, analogous to the phase transitions in mixed states with instrinsic topological order [17].We will not consider a general SPT state, and focus primarily on cluster states in various dimensions to illustrate the broad idea.A key step in our analysis will be the following result that was also briefly mentioned in Ref. [17] and which we discuss in detail in Sec.IV: for a large class of SPTs, including the cluster states in various dimensions, Kitaev chain in 1d, and several 2d topological phases protected by zero-form Z 2 symmetry, one can find local, finite-depth channels that map the pure state to a Gibbs state.We will discuss decoherence induced separability transitions due to such channels in Sec.IV.
When trying to understand complexity of mixed SPT states, we will often find the following line of inquiry helpful.One first asks: Does assuming that a mixed state is trivial (i.e.decomposable as a convex sum of SRE pure states) lead to an obvious contradiction?If the answer to this question is 'yes', then we already know that the mixed state is necessarily non-trivial.In this case, there may still exist interesting transitions between two different kinds of non-trivial mixed states, and we will consider a couple of such examples as well.On the other hand, if the answer to this question is 'no', we will attempt to find an explicit decomposition of the mixed state as a convex sum of SRE states.The aforementioned relation between local and thermal decoherence will again be instrumental in making analytical progress.
As an example, consider the ground state of the 2d cluster state Hamiltonian H subjected to a local channel that locally anticommutes with the terms in the Hamiltonian.One can show that resulting decohered state ρ d takes the Gibbs form: ρ ∝ e −βH where tanh(β) = 1 − 2p and p is the decoherence rate.In this example, H has both a zero-form and a one-form Ising symmetry.We will provide arguments that this system undergoes a separa-bility transition as a function of p: for 0 < p < p c , ρ cannot be decomposed as a sum of pure states that respect the aforementioned two symmetries, while for p > p c , such a decomposition is feasible.Moreover, for p > p c we will express ρ d explicitly as m p m |ψ m ⟩⟨ψ m |, where |ψ m ⟩ are pure, symmetric states that are statistically SRE.More precisely, one can define an ensemble averaged string-correlation, [⟨S C ⟩ 2 ] ≡ m p m ⟨S C ⟩ 2 m , where ⟨S C ⟩ m = ⟨ψ m |S C |ψ m ⟩/⟨ψ m |ψ m ⟩, and S C is a stringoperator whose non-zero expectation value implies longrange entanglement.We will show that [⟨S C ⟩ 2 ] precisely corresponds to a disorder-averaged correlation function in the 2d random-bond Ising model along the Nishimori line [35].Therefore, in this example, the separability transition maps to the ferromagnetic transition in the random bond Ising model.For the 3d cluster state, we will find an analogous relation between separability and 3d random-plaquette Ising gauge theory.We note that similar order parameters and connections to statistical mechanics models also appear in the setting of measurement protocols to prepare long-range entangled SPT states [36,37].We briefly discuss connection to these works.
As another byproduct of the relation between local decoherence and Gibbs states, we also study recently introduced non-trivial class of mixed states which are protected by a tensor product of 'exact' and 'average' symmetries [28][29][30][31].One says that a density matrix ρ has an 'exact symmetry' if U E ρ ∝ ρ for some unitary U E , while it has an 'average symmetry' if U † A ρU A = ρ for some unitary U A .Refs. [28][29][30][31] have provided several non-trivial examples of such mixed-state SPTs by showing that they possess non-trivial correlation functions, and/or cannot be purified to a short-ranged entangled (SRE) pure state.Here we will focus on examples of such states that are based on cluster states in various dimensions, and using locality/Lieb-Robinson bound [38][39][40], show that the corresponding mixed states cannot be written as a convex sum of symmetric, pure states.For 1d cluster state, we also provide an alternative proof of non-separability by using the result from Ref. [41] that in one-dimension if a state has an average Z 2 symmetry, and its connected correlation functions are short-ranged, then the corresponding 'order parameter' and the 'disorder parameter' can't be both zero or non-zero at the same time.
Next, in Sec.V we consider fermionic chiral states subjected to local decoherence.We primarily focus on the ground state of a 2d p x + ip y superconductor (p + ip SC) as our initial state (we expect integer quantum Hall states to have qualitative similar behavior).We first consider subjecting this pure state to a finitedepth channel with Kraus operators that are linear in fermion creation/annihilation operators, so that the decoherence breaks the fermion parity symmetry.In the pure state classification of topological superconductors, fermion parity is precisely the symmetry responsible for the non-trivial topological character of the p + ip SC [42,43].Therefore, it is natural to wonder about the fate of the mixed state obtained by breaking this symmetry from exact down to average.One potential path to make progress on this problem is to map the mixed state to a pure state in the doubled Hilbert-space using the Choi-Jamiolkowski (C-J) map [44,45] (we will call such a state the 'double state', similar to the nomenclature in Ref. [20]).There are interesting subtleties in applying the C-J map to fermionic Kraus operators that we clarify.Following the ideas in Refs.[18,20,22,31], one may then map the double state to a 1+1-D theory of counter-propagating free CFTs coupled via a fermion bilinear term, which is clearly relevant and gaps out the edge states in these doubled picture.However, a shortdepth channel cannot qualitatively change the expectation value of state-independent operators (i.e.tr(ρO) where O is independent of ρ) [18,19], and it is not obvious what does the gapping of edge modes imply for the actual mixed state.We conjecture that the physical implication of the gapping of the edge states in the doubled formulation is that the actual mixed state can now be expressed as a convex sum of SRE states with zero Chern number, which is equivalent to the statement that they can be obtained as a Slater determinant of Wannier states, unlike the pure p + ip state where such a representation is not possible [46][47][48].Therefore, the transition from the pure state to the mixed state can be thought of as a 'Wannierizability transition'.We consider an explicit ansatz of such a decomposition, and provide numerical support of our conjecture by calculating the entanglement spectrum and modular commutator of the pure states whose convex sum corresponds to the decohered density matrix.
A more interesting channel that acts on the 2d p + ip SC corresponds to Kraus operators that are bilinear in fermion creation/annihilation operators.To make progress on this problem, we use the C-J map to obtain a field theoretic description for this problem in terms of two counter-propagating chiral Majorana CFTs interacting via a four-fermion interaction, where the strength of the interaction is related to the strength of the interacting decohering channel.This theory admits a phase transition at a critical interaction strength in the supersymmeteric tricritical Ising universality class, which can be thought of as corresponding to spontaneous breaking of the fermion parity.Although we don't have an understanding of this transition directly in terms of the mixed state in the non-doubled (i.e.original) Hilbert space, it seems reasonable to conjecture that at weak decoherence, the density matrix can not be expressed as a convex sum of area-law entangled non-chiral states, while at strong decoherence, it is most naturally expressible as a convex sum of states with GHZ like character that originates from the aforementioned spontaneous breaking of the fermion parity.
Incidentally, the kind of arguments we consider to rule out sym-SRE mixed states in the context of symmetry broken phases or SPT phases also find an application in an exotic separability transition where symmetry plays no role.In particular, we consider separability as-pects of Gibbs state of Hamiltonians that satisfy the 'no low-energy trivial state' (NLTS) condition introduced by Freedman and Hastings in Ref. [49].Colloquially, if a Hamiltonian satisfies the NLTS condition, then any pure state with energy density less than a critical non-zero threshold cannot be prepared by a constant depth circuit.Recently, Ref. [50] showed that the 'good LDPC code' constructed in Ref. [51] satisfies the NLTS condition (we note that 'Good LDPC codes' [51][52][53] have the remarkable property that both the code distance and the number of logical qubits scale linearly with the number of physical qubits).Ref. [50] already showed that the NLTS condition holds also for mixed states, if one defines the circuit depth of a mixed state as the minimum depth of unitary needed to prepare it by acting on system⊗ancillae, both initially in a product state, where the ancillae are traced out afterwards [54].Under such a definition of a non-trivial mixed state (namely, a mixed state that can not be prepared by a constant depth circuit under the aforementioned protocol), even mixed states with long-range classical correlations (e.g. the Gibbs state of 3d classical Ising model) would be considered non-trivial.In contrast, under our definition of a non-trivial mixed state, such classical states will be trivial since they can be written as a convex sum of SRE states.Therefore we ask: assuming that one defines a trivial (non-trivial) mixed state as one which can (can't) be expressed as a convex sum of SRE states, is the Gibbs state of a Hamiltonian that satisfies the NLTS property non-trivial at a low but non-zero temperature?Under reasonable assumptions, in Sec.VI we provide a short argument that this is indeed the case.This implies that one should expect a non-zero temperature separability transition in such Gibbs states.
In Sec.VII we briefly discuss connections between separability criteria and other measures of the complexity of a mixed state such as the ability to purify a mixed state to an SRE pure state, entanglement of the doubled state using C-J map, and strange correlators.
Finally, in Sec.VIII we summarize our results and discuss a few open questions.
For convenience, below we list some of the acronyms that appear frequently in this work: -SRE: Short-range entangled.

II. SEPARABILITY CRITERIA WITH AND WITHOUT SYMMETRY
Motivated from Werner and Hastings [1,2], we call a mixed state ρ short-range entangled (SRE) if and only if it can be decomposed as a convex sum of pure states where each |ψ m ⟩ is short-range entangled (SRE), i.e., it can be prepared by applying a constant-depth local unitary circuit to some product state.The physical motivation for this definition is rather transparent: if a mixed state can be expressed as Eq.( 1), only then it can be prepared using an ensemble of unitary-circuits (acting on the Hilbert space of ρ) whose depth does not scale with the system size.We note that this definition of an SRE mixed state has been employed to understand phase transitions in systems with intrinsic topological order subjected to thermal or local decoherence [17,55].
One can generalize the notion of an SRE mixed state in the presence of a symmetry.Specifically, we say that a mixed state ρ satisfying U (g)ρU † (g) = ρ, ∀g ∈ G is a 'symmetric SRE' (sym-SRE in short) if and only if one can decompose it as a convex sum of pure states, where each of these pure states can be prepared by applying a finite-depth quantum circuit made of local gates that all commute with U , to a symmetric product state.
Several comments follow.
1.The 'only if' clause in our definition for a (sym-) SRE state is a bit subtle.For example, consider a density matrix where there exists no decomposition that satisfies Eq.( 1), but there exists q m |ϕ m ⟩⟨ϕ m | such that the relative weight of the non-SRE states is zero in thermodynamic limit (i.e.m q m /( m (p m + q m ) → 0 in the thermodynamic limit).In this case, it might seem reasonable to regard ρ as SRE.One may also define an average circuit complexity of a density matrix as ⟨C⟩ = inf{ m p m C(ψ m )}, where C(|ψ m ⟩) is the minimum depth of a circuit composed of local gates to prepare the state |ψ m ⟩ and the infimum is taken over all possible decompositions of the mixed state ρ.One may then consider calling a mixed state ρ as SRE if and only if ⟨C⟩ does not scale with the system size.But even then, there may be special cases where the average behavior is not representative of a typical behavior.We will not dwell on this subtlety further at this point, and use physical intuition to quantify the separability of a density matrix, were we to encounter such a situation.
where U is a finite-depth circuit and Z = tr e −H cl .We are unable to show that the definition in Eq.( 1) is equivalent to Eq.( 2).Although we will primarily use the former definition (Eq.( 1)), in Sec.VII we will briefly discuss potential connections between the two definitions, and also relation with other diagnostics of mixed-state entanglement.
3. The symmetry (U (g)ρU † (g) = ρ, ∀g ∈ G) we consider is called weak symmetry (average symmetry) in Ref. [30] (Ref.[28]), which highlights its difference with the stronger symmetry U (g)ρ = ρU (g) = e iθ(g) ρ, ∀g ∈ G, termed strong symmetry (exact symmetry) in Ref. [30] (Ref.[28]).Physically, exact symmetry enforces the constraint that the density matrix must be written as an incoherent sum of pure states, where each of them is an eigenstate of U (g) with the same eigenvalue e iθ(g) .On the other hand, while the mixed state ρ with only average symmetry can be written as a convex sum of symmetric pure states having different charge under G, one may as well express ρ as a convex sum of nonsymmetric pure states.Therefore, our requirement that each of the pure states respects the symmetry puts a further constraint on a mixed state with only average symmetry.
On that note, Ref. [56] defined a symmetric-SRE state for a symmetry U as one which satisfies Eq.( 2) where e −H cl is replaced by P θ(g) e −H cl where P θ(g) is a projector onto a given symmetry charge θ(g).Therefore, in this definition one is always working with a density matrix that has an exact symmetry.As already mentioned, we will instead only impose the average symmetry in our definition of a sym-SRE state (of course, there may be special quantum channels that happen to preserve an exact symmetry).
4. An alternative definition of an SRE mixed state was considered in Refs.[28,29,54] whereby a mixed density matrix is considered SRE if it can be obtained from a pure product state in a system⊗ancillae Hilbert space via a finite-depth unitary followed by tracing out ancillae.In contrast, as already mentioned above in comment #2, Ref. [2] defines a mixed density matrix as SRE if it can be obtained from the 'classical mixed state' ρ cl ∝ e −H cl of system⊗ancillae via a finite-depth local quantum channel.Therefore, a mixed state can be trivial/SRE using the definition of Ref. [2] while remaining non-trivial/LRE using the definition of Refs.[28,29,54].The physical distinction between these two definitions is most apparent when one considers a mixed state for qubits of the form This state is clearly separable (unentangled).However, any short-depth purification of this state must be long-ranged entangled.This is because tr(ρZ i Z j ) − tr(ρZ i ) tr(ρZ j ) is non-zero and the purified state can't change this correlation function due to the Lieb-Robinson bound [38,39] (this is also related to the fact the entanglement of purification [57] is sensitive to both quantum and classical correlations, and therefore is not a good mixedstate entanglement measure).Thus the aforementioned ρ will be SRE using definition of Ref. [2], and LRE using the definition of Ref. [28,29,54].Of course, it will also be SRE via Eq.( 1), which is the definition we will use throughout this paper.

III. AN ILLUSTRATIVE EXAMPLE: SEPARABILITY TRANSITION IN THE GIBBS STATE OF THE 2D QUANTUM ISING MODEL
Let us consider an example to illustrate the difference between an SRE mixed state and a sym-SRE mixed state, that will also provide one of the simplest examples of a separability transition.Consider the density matrix ρ for qubits (i.e.objects transforming in the spin-1/2 representation of SU (2)) given by ρ(β) = e −βH /Z where H is a local Hamiltonian that satisfies U † HU = H with U = i X i being the generator of the Ising symmetry, and Z = tr e −βH is the partition function.Let us further assume that ρ(β) exhibits spontaneous symmetry breaking (SSB) for β > β c where 0 < β c < ∞ (for a range of other parameters that specify the Hamiltonian).For concreteness, one may choose H as the nearest neighbor transverse-field Ising model on the square lattice, i.e., H = − ⟨i,j⟩ Z i Z j − h i X i although the only aspect that will matter in the following discussion is that H is local with a zero-form Ising symmetry, and the order parameter in the symmetry breaking phase is a real scalar (e.g. one may as well consider a transverse-field Ising model on a cubic lattice).Therefore, for a range of the transverse-field h and β > β c (where β c depends on h), the two-point correlation function tr (ρZ i Z j ) is non-zero for |i − j| → ∞.We will argue that ρ is SRE for all nonzero temperatures, while it is sym-SRE only for β < β c .Partial support for ρ being an SRE at all non-zero temperatures was provided in Refs.[25][26][27] and we will argue below for an explicit decomposition of ρ in terms of SRE states.
The statement that ρ is not sym-SRE for β ≥ β c was also hinted in [56], and intuitively follows from the fact that for β > β c , SSB implies that if one decomposes ρ as a convex sum of symmetric, pure states, those pure states must have GHZ-like entanglement.Let us first consider a rigorous argument for this statement which, upto small modifications, essentially follow the one in Ref. [21] for a closely related problem of non-triviality of a density matrix with an exact symmetry and long-range order.
To show that for β > β c , ρ can't be a sym-SRE state, let us first decompose ρ as ρ = ρ + + ρ − where ρ ± = ( 1±U 2 )ρ are the projections of ρ onto even and odd charge of the Ising symmetry.ρ + and ρ − are valid density matrices with an exact Ising symmetry, that is, they satisfy, U ρ ± = ±ρ ± .Now let us make the assumption that for β > β c , ρ is a sym-SRE state.We will show that this assumption leads to a contradiction.Therefore, we write ρ ± = α p α,± |ψ α,± ⟩⟨ψ α,± | where p α,± are positive numbers, and |ψ α,± ⟩ are SRE states ∀α that satisfy U |ψ α,± ⟩ = ±|ψ α,± ⟩.Since U anti-commutes with Z i , ⟨ψ α,± |Z i |ψ α,± ⟩ = 0. Further, since |ψ α,± ⟩ are all SRE states, correlation functions of all local operators decay exponentially (notably, we assume that the associated correlation length is bounded by a systemsize independent constant for all |ψ α,± ⟩), and therefore, However, this leads to a contradiction, because this implies that tr (ρZ j Z k ) = ± α p α,± ⟨ψ α,± |Z j Z k |ψ α,± ⟩ itself vanishes, which we know can't be true since as mentioned above, for β > β c , the system is in an SSB phase with long-range order.Therefore, our assumption that ρ is a sym-SRE state for β > β c must be incorrect.The same conclusion also holds for β = β c since the correlations at the critical point decay as a power-law.
As mentioned in the introduction, our general approach would be to first look for general constraints that lead to a mixed state being necessarily non-trivial.If we are unable to find such a constraint, we will attempt to find an explicit decomposition of the density matrix as a convex sum of SRE states.For example, above, we noted that ρ cannot be a sym-SRE state for β ≥ β c , and we also claimed that ρ is an SRE state for all nonzero temperatures.Let us therefore try to find an explicit decomposition of ρ as a convex sum of SRE pure states for any non-zero temperature, and as a convex sum of symmetric, pure SRE states for β < β c .The key player in our argument will be a particular convex decomposition ansatz (CDA in short) that is motivated from "minimally entangled typical thermal states" (METTS) construction introduced in Ref. [58], and which was employed in Ref. [55] to show that the Gibbs state of 2d and 3d toric code is SRE for all non-zero temperatures.Note that despite the nomenclature, METTS construction as introduced in [58] does not involve minimization of entanglement over all possible decompositions, and is simply an ansatz that is physically motivated (which is why we prefer the nomenclature CDA over METTS for our discussion).
First, let us specialize to zero transverse field.In this case, ρ is clearly an SRE state at any temperature since ρ ∝ m e −βEm |z m ⟩⟨z m | where |z m ⟩ denotes a product state in the Z-basis and E m = ⟨z m |H|z m ⟩.To obtain a symmetric convex decomposition, we write: where the set {|x m ⟩} corresponds to the complete set of states in the X basis, |ψ m ⟩ = e −βH/2 |x m ⟩/ √ Z is the unnormalized wave function.The states |ψ m ⟩ are clearly symmetric under the Ising symmetry, and their symmetry charge (= ±1) is determined by the parity of the number of sites in the product state |x m ⟩ where spins point along the negative-x direction.We will now argue that the states |ψ m ⟩ are SRE for β < β c and LRE for β ≥ β c .To see this, we first consider the "partition function with respect to |ψ m ⟩" defined as Z m = ⟨ψ m |ψ m ⟩ and study its analyticity as a function of β.In this specific example, since transverse field is set to zero, one finds that for all m, Z m is simply proportional to the partition function of the 2d classical Ising model at inverse temperature β, and therefore is non-analytic across the phase transition.Similarly, the two-point correlation function ⟨ψ m |Z i Z j |ψ m ⟩/⟨ψ m |ψ m ⟩ is just the two-point spin-spin correlation function in the 2d classical Ising model, which is long-ranged for β ≥ β c and exponentially decaying for β < β c .These observations strongly indicate that |ψ m ⟩ is SRE (and correspondingly, ρ sym-SRE) if and only if β < β c .Note that the states |ψ m ⟩ are expected to be area-law entangled for all β.This is because one may represent the imaginary time evolution e −βH |m⟩ as a tensor network of depth β acting on |m⟩ (which is a product state), which can only generate an area-law worth of entanglement.Further, even the state at β = ∞ is area-law entangled (= the ground state of H).Therefore, shortrange correlations are strongly suggestive of short-range entanglement.
Now, let's consider non-zero transverse field.To argue that ρ is SRE for any non-zero temperature, we again decompose it as ρ = m |ψ m ⟩⟨ψ m | where |ψ m ⟩ = e −βH/2 |z m ⟩/ √ Z.The corresponding Z m = ⟨ψ m |ψ m ⟩ can now be expressed in the continuum limit as an imaginarytime path integral x,y r|ϕ| 2 + u|ϕ| 4 , ω n = 2πn/β are the Matsubara frequencies, and, crucially, the Dirichlet boundary conditions ϕ(x, y, τ = 0) = ϕ(x, y, τ = β) = ϕ 0 (x, y) are imposed by the 'initial' state z m ∼ ϕ 0 (x, y).Since β ̸ = ∞, the discrete sum over the Matsubara frequencies will be dominated by ω n = 0, which corresponds to space-time configurations that are translationally invariant along the imaginary-time-direction.Furthermore, the Dirichlet boundary conditions imply that there is just one such configuration, namely, ϕ(x, y, τ ) = ϕ 0 (x, y) such that Z m ∼ e S[ϕ0(x,y)] , and thus the fluctuations of ϕ will be completely suppressed at all non-zero temperatures (including at the finite temperature critical point which corresponds to renormalized r = 0) Therefore, we ex-pect that Z m will not exhibit singularity across the finite temperature critical point, which indicates that the states |ϕ m ⟩ are SRE. To Crucially, since the initial state is now a product state in the X basis, the fields at the two boundaries τ = 0, β are being integrated over all possible configurations Again, the path integral will be dominated by ω n = 0 which only implies that the dominant contribution comes from configurations ϕ(τ, x, y) = ϕ(x, y).Therefore, unlike the aforementioned case when the CDA states corresponded to e −βH/2 |z m ⟩/ √ Z, here dominant contribution to Z m precisely corresponds to the partition function of the 2d classical Ising model, which is in the paramagnetic phase for β < β c .The correspondence with 2d classical Ising model makes physical sense since the universality class of the phase transition at any non-zero temperature is indeed that of the 2d classical Ising model.Therefore, we expect that the states |ψ m ⟩ = e −βH/2 |x m ⟩/ √ Z are SRE for β < β c and LRE for β ≥ β c .Correspondingly, we expect that the Gibbs state is sym-SRE (sym-LRE) for β < β c (β > β c ).
To summarize, we provided arguments that the Gibbs state of a transverse field Ising model is an SRE state at any non-zero temperature, and a sym-SRE state only for β < β c .Therefore, we expect that it undergoes a separability transition as a function of temperature if one is only allowed to expand the density matrix as a convex sum of symmetric states.We expect similar statements for other models that exhibit a finite temperature zeroform symmetry breaking phase transition.In the following sections, we will employ broadly similar logic as in this example, with primary focus on topological phases of matter subjected to local decoherence.Specifically, we write ρ = ΓΓ † and employ the following CDA: where |ψ m ⟩ = Γ|m⟩, p m = ⟨m|Γ † Γ|m⟩ = ⟨ψ m |ψ m ⟩, and | ψm ⟩ = |ψ m ⟩/ ⟨ψ m |ψ m ⟩ are normalized versions of |ψ m ⟩.We note that here Γ is not unique (note that Γ is not restricted to be a square matrix, see e.g.Ref. [17]), and CDA in Eq.(3) corresponds to choosing Γ = ρ 1/2 for the Gibbs state ρ.We will sometimes call states {|ψ m ⟩} that enter a particular CDA as 'CDA states'.We further note that in general we do not know how to find the matrix Γ that is 'optimal' i.e. a matrix Γ that guarantees that the states Γ|m⟩ are SRE whenever ρ is SRE.However, as we will see in the rest of the paper, for a large class of problems, after a judicious choice of the basis {|m⟩}, the decomposition in the form of Eq.4 turns out to be optimal.

IV. SEPARABILITY TRANSITIONS IN SPT STATES
The fundamental property of a non-trivial SPT phase is that it cannot be prepared using a short-depth circuit consisting of local, symmetric, unitary gates [9][10][11][12].Therefore, it is natural to ask: if an SPT phase is subjected to local decoherence, is the resulting mixed state sym-SRE, i.e., can it be expressed as a convex sum of symmetric, SRE pure states?This is clearly a very challenging question for many-body mixed states since to our knowledge, there does not exist an easily calculable measure of mixed-state entanglement that is non-zero if and only if the mixed state is unentangled [59] (if such a measure did exist, then it would be useful to study its universal, long-distance component, similar to topological part of negativity [19,55,60]).As already hinted in the introduction, our general scheme will be to first seek sufficient conditions that make a given mixed state sym-LRE (i.e.not sym-SRE).We will do this by decomposing the decohered state into its distinct symmetry sectors as ρ = Q ρ Q , with ρ Q the projection of the density matrix onto symmetry charge Q, and then examining whether the assumption of each ρ Q being an SRE leads to a contradiction.If we are unable to find an obvious contradiction, we will then attempt to use the decomposition outlined in Eq.( 4) to express ρ as a convex sum of sym-SRE states.In either of these steps, we will exploit the connection between local and thermal decoherence for cluster states that was briefly mentioned in Ref. [17], and which is described in the next subsection in detail.

A relation between local and thermal decoherence
Systems with intrinsic topological order typically behave rather differently when they are coupled to a thermal bath, compared to when they are subjected to decoherence induced by a short-depth quantum channel.For example, when 2d and 3d toric codes are embedded in a thermal bath, so that the mixed state is described by a Gibbs state, the topological order is lost at any non-zero temperature [2,55,61,62].In contrast, when 2d or 3d toric codes are subjected to local decoherence, then the error-threshold theorems [63][64][65][66][67][68] imply that the mixedstate topological order is stable upto a non-zero decoherence rate [17-20, 61, 69].Given this, it is interesting to ask if there exist situations where a local short-depth channel maps a ground state to a Gibbs state.Here we show that this is indeed the case if the corresponding Hamiltonian satisfies the following properties: (1) It can be written as a sum of local commuting terms where each of them squares to identity: (2) There exists a local unitary O j which anticommutes (commutes) with Specifically, denoting the total system size as N , the channel maps the ground state density matrix ρ 0 to a Gibbs state for H.
To verify the claim, we first note that Eq.( 5) implies that ρ 0 can be written as the product of the projectors on all sites ρ 0 = 1 2 N j (I − h j ).Now, using Eq.( 6), it is straightforward to show that Since h 2 j = I, which implies e −βhj = cosh(β)I − sinh(β)h j , one may now exponentiate Eq.( 8) to obtain E[ρ 0 ] = 1 Z e −βH where tanh β = (1 − 2p), Z = tr e −βH .In Sec.VIII, we also discuss a Z N generalization of this construction.For the rest of the paper, the aforementioned Z 2 version will suffice.In the following we will exploit the connection between local and thermal decoherence to study decoherence-induced separability transitions for the cluster states in various dimensions (Secs.IV A,IV B,IV C).We will also briefly discuss a couple examples where the pure state is protected by a single zero-form symmetry (Sec.IV D).

A. 1d cluster state
The Hamiltonian for the 1d cluster state is where a and b denote the two sublattices of the 1d chain, see Fig. 1(a).H has a global Z 2 ×Z 2 symmetry generated by We assume periodic boundary conditions, so that their is a unique, symmetric, ground state of H which is separated from the rest of the spectrum with a finite gap.
It is obvious that H satisfies Eq.( 5).To satisfy Eq. ( 6), we choose Kraus operators O a/b,j = Z a/b,j .Therefore, under the composition of the channel E a/b,j [ρ] = (1 − p a/b )ρ + p a/b Z a/b,j ρZ a/b,j on all sites, the pure state density matrix becomes with tanh ) and Z a/b = tr e −β a/b j h a/b,j .In the following, we will suppress the arguments p a , p b in ρ a (p a ), ρ b (p b ) if there is no ambiguity.Note that ρ a and ρ b commute with each other.To decompose ρ as a convex sum of symmetric states, we write ρ = Qa,Q b ρ Qa,Q b , where each ρ Qa,Q b is an unnormalized density matrix that carries exact symmetry: Note that the probability for a given sector (Q a , Q b ) is given by tr (ρ Qa,Q b ) which can be used to obtain the normalized density matrix ρQa, To discuss whether the decohered mixed state ρ is trivial based on our definition of sym-SRE mixed state, we start from considering the special case p a > 0, p b = 0, i.e., the mixed state obtained by applying the aforementioned quantum channel only on sublattice a.This case was studied in detail in Ref. [29] from a different perspective and is an example of an 'average-SPT' phase [28,29,31,32].In particular, it was shown in Ref. [29] that this mixed state can not be purified to an SRE pure state using a finite-depth local quantum channel.As discussed in Sec.II, our definition of SRE mixed state is a bit different (namely, whether a mixed state can be written as a convex sum of SRE pure states), and therefore, it is worth examining whether this state continues to remain an LRE mixed state using our definition.
When p a > 0, p b = 0, only the sector corresponding to Q b = 0 survives, and in this sector, ρ Qa,Q b ∝ j (I − h b,j )e −βa j ha,j P Qa .We will now provide two separate arguments that show that ρ Qa,Q b is a sym-LRE (i.e.not sym-SRE) mixed state when p a > 0, p b = 0. m p m |ψ m ⟩⟨ψ m | where |ψ m ⟩ are SRE states that can be prepared via a short-depth circuit consisting of symmetric, local gates.We utilize the result in Ref. [41], which shows that for an area-law entangled state in 1D (which we will take to be |ψ m ⟩) which is symmetric under an Ising symmetry (which we will take here to be U a = j X a,j ), both order and disorder parameters cannot vanish simultaneously.Note that we are assuming that |ψ m ⟩ has an area-law entanglement, as otherwise, it is certainly not SRE and there is nothing more to prove.
Therefore, following the results in Ref. [41], |ψ m ⟩ must either (a) have a non-zero order parameter corresponding to the symmetry U a i.e. ⟨ψ m | Zj Zk |ψ m ⟩ ̸ = 0 where |j − k| ≫ 1 and Z is an operator that is odd under U a e.g.Zi = Z a,i , or, (b) it must have a non-zero 'disor-der parameter' corresponding to the symmetry = 0 where |j − k| ≫ 1, and O L , O R are operators localized close to site j and k respectively that are either both even or both odd under U a .In case (a), the system has a long-range GHZ type order since the state |ψ m ⟩ is symmetric under U a .In case (b), we now argue that the system has an SPT order.
As a consequence, for ⟨ψ m |O L k l=j X a,l O R |ψ m ⟩ to be non-zero, O L and O R must both be odd under U b .If so, then the disorder order parameter precisely corresponds to one of the two SPT string order parameters, namely, S a (j, k) = Z b,j−1 k l=j X a,l Z b,k+1 up to finite depth symmetric unitary transformation.At the same time, the other SPT string order parameter )), and therefore, we arrive at the conclusion that in case (b), |ψ m ⟩ must possess non-trivial SPT order since string order parameters on both sublattices are non-zero.Therefore, in either case (a) or (b), |ψ m ⟩ cannot be prepared by a short-depth circuit composed of local gates that respect both U a and U b , starting with a symmetric product state.
Second argument: This argument is essentially the same as the one introduced in Ref. [40] to show that the circuit depth of various states with a non-trivial string order parameter cannot be a system-size independent constant due to locality/Lieb-Robinson bound [38,39].Again, recall that we want to show that ρ Qa,Q b ∝ j (I − h b,j )e −βa j ha,j P Qa can not be written as m p m |ψ m ⟩⟨ψ m | where |ψ m ⟩ are SRE.Since ρ Qa,Q b carries an exact symmetry charge of U a , U b , so do each of the pure states |ψ m ⟩.As already discussed above, the expectation value of the string order parameter S b (j, k) = k l=j (−h b,l ) = Z a,j k l=j X b,l Z a,k+1 is unity with respect to ρ Qa,Q b , which implies that its expectation value is also unity with respect to each of the states |ψ m ⟩.Let us assume that |ψ m ⟩ can be obtained from a symmetric product state (i.e. an eigenstate of Pauli X on all sites) which we denote as |x a,b ⟩ = ⊗ j |x a,j , x b,j ⟩, i.e., |ψ m ⟩ = V |x a,b ⟩ (here x a/b,j = ±1 are chosen so as to satisfy the symmetry Moreover, since V is a finite-depth unitary, the operator V † S b (j, k)V is still a string operator with each 'end point operator' V † Z a,j/k+1 V a sum of local operators (due to the locality of V ) that are charged under U a (due to V being a symmetric unitary, i.e., [V, which is unity, as discussed above.Therefore, we arrive at a contradiction.This implies that our assumption that |ψ m ⟩ is a symmetric SRE state must be incorrect. We now discuss the general case of both p a and p b being non-zero.Based on our discussion above, it is instructive to evaluate the string order parameter with respect to each ρ Qa,Q b , i.e., tr ρ Qa,Q b S a/b / tr(ρ Qa,Q b ).One finds (see Appendix A) that both string order parameters can be mapped to two-point correlation functions of spins in the 1d classical Ising model at nonzero temperature and hence decay exponentially with the length of the strings.This result merely implies that the corresponding mixed state ρ = Qa,Q b ρ Qa,Q b doesn't satisfy the aforementioned sufficient condition for non-trivial sym-SRE, and does not guarantee that ρ must be trivial.We now use the CDA in Eq.( 4) to argue that ρ is indeed sym-SRE.In particular, we choose One can also compute the string order parameters S a (S b ) for |ψ m ⟩ and show its equivalence to ⟨z j z k ⟩ 1D Ising at inverse temperature β a (β b ).Therefore, |ψ m ⟩ does not develop string order as long as β a/b < ∞.The triviality of |ψ m ⟩ is also manifested by the non-zero expectation value of disorder operator U a/b (k, j) = k l=j X a/b,l .For example, consider the expectation value of the disorder operator on a sublattice: Using the fact that the only terms in e −(βa j ha,j +β b j h b,j )/2 that anticommutes with U a (k, j) are h b,j−1 and h b,k , we find that ⟨U a (k, j)⟩ m = ( k l=j x l ) sech 2 (β a ), which is nonvanishing except for β a = ∞.This is of course expected based on the result of Ref. [41], since |ψ m ⟩ does not have any GHZ type order.The result for U b (k, j) is similar.
It is also instructive to apply the aforementioned convex decomposition to the case β b = ∞, β a ̸ = ∞, i.e., the above discussed case of 'average SPT order'.In this case we find that the corresponding state |ψ m ⟩ develops GHZ type long-range entanglement.To see this, one can rewrite |ψ m ⟩ as |ψ m ⟩ ∼ e −βa k ha,j /2 |χ m ⟩, where |χ m ⟩ ∼ Using the fact that the only terms in e −βa j ha,j /2 that anticommute with Z a,j Z a,k are h a,j and h a,k , one finds that To summarize the results in this subsection, the decohered state ρ as a function of p a and p b can be divided into four regimes (see Fig. 1): . This is just the pure state SPT.
(ii) p a > 0 and p b = 0: tr(ρ Qa,Q b S a (j, k)) decays exponentially with |j − k| and tr(ρ This regime is sym-LRE i.e. a non-trivial mixed state, in agreement with the non-trivial 'average SPT' discussed in Ref. [29]. (iii) p a = 0 and p b > 0: this is similar to the case (ii) with a ↔ b and is again a sym-LRE state.
(iv) p a , p b > 0: both tr(ρ Qa,Q b S a (j, k)) and tr(ρ Qa,Q b S b (j, k)) decay exponentially with |j − k|.This is a sym-SRE state.
Based on our discussion above, we also provide one possible 'phase diagram' to express ρ as a convex sum of symmetric states using CDA states |ψ m ⟩ = ρ 1/2 |x a , x b ⟩, as summarized in the third column of Fig. 1(a).Note that the boundary of the phase diagram using the employed CDA matches the boundary of regimes (i)-(iv), and therefore, the CDA is optimal in this sense.However, it's worth noting that the decomposition we chose is just one possible choice, and the label 'GHZ' on the x and y axis in the third column of Fig. 1(a) is tied to this choice.One may also chose to expand ρ as a convex sum of SPT states.Therefore, the result that is independent of any specific choice of CDA is that the regime (iv) is sym-SRE, while the regimes (i), (ii) and (iii) are sym-LRE.

B. 2d cluster state
The 2d cluster state Hamiltonian H 2d Cluster is: Here the Hilbert space consists of qubits residing on both the vertices v and the edges e of a 2D square lattice, see with the corresponding generators where p labels the plaquette on the lattice and ∂p is the boundary of p.We assume periodic boundary conditions, so that H has a unique, symmetric, gapped ground state.Using Eqs.( 5),( 6), if one subjects the ground state of H 2d Cluster to Kraus operators O v/e = Z v/e with respective probabilities p v/e , the resulting decohered density v +βe e h (1)   e ) with tanh Let us decompose ρ as a convex sum of symmetric states by writing ρ = , where each ρ Q (0) ,Q (1) carries the exact symmetry: . Here, the one-form symmetry charge is labeled by the set p = 0, 1 defined on each plaquette p. Crucially, the number of one-form symmetry sectors grows exponentially as a function of the system size, and this implies that the probability for a given sector (Q (0) , Q (1) ), i.e., tr ρ Q (0) ,Q (1) , is exponentially small in general.It follows that even if there exists some ρ that is not sym-SRE, the decohered state ρ may still be well approximated by a sym-SRE mixed state as long as the total probability corresponding to the non-trivial sectors is exponentially small.Therefore, the notion of ρ being sym-SRE must take into account the probability for each symmetry sector, and can only be made precise in a statistical sense (a similar situation arises for a certain non-optimal decomposition for decohered toric code [17].We will return to this point in detail below.For now, let's focus on the physical observables in each symmetry sector.
The observables that characterize the 2d cluster ground state are the expectation value of the membrane operator M S = v∈S (−h v ) with S a surface (for simplicity, we will assume that the boundary ∂S of this surface is contractible), and the string operator e ) with C a curve (the expectation value of either of these operators equals unity in the 2d cluster ground state).To detect whether ρ Q (0) ,Q (1) is sym-SRE, i.e., it can be expanded as a convex sum of pure SRE states that each carries a definite symmetry charge (Q (0) , Q (1) ), it is instructive to calculate the expectation value of these operators with respect to ρ To proceed, we first compute the denominator in these expressions, i.e., tr ρ Q (0) ,Q (1) .
Similar to the 1d cluster state, this can be easily done by inserting the complete basis {|x e,v ⟩} and {|z e,v ⟩}, where |x e,v ⟩ = ⊗ e,v |x e , x v ⟩ and |z e,v ⟩ = ⊗ e,v |z e , z v ⟩ denote the product state in Pauli-X and Z basis, respectively.Following a calculation quite similar to that in the 1D cluster state, one finds that tr ρ Here, Z 2D gauge,xv = ze e βv v xv( e∋v ze) is the partition function of the 2D Ising gauge theory with the sign of interaction on each vertex given by x v while Z 2D Ising,xe = zv e βe e xe( v∈e zv) is the partition function of the 2D Ising model with the sign of Ising interaction given by x e .
In the summation, the notation x v ∈ Q (0) denotes all possible x v which satisfy v x v = (−1) Q (0) while x e ∈ Q (1) denotes all possible x e which satisfy e∈∂p x e = (−1) Q (1) p , ∀p.For a system with periodic boundary conditions, all possible x v ∈ Q (0) (x e ∈ Q (1) ) can be reached by the transformation One may verify that Z 2D gauge,xv (Z 2D Ising,xe ) is invariant under the aforementioned transformation by changing the dummy variables z e → σ e z e (z v → s v z v ).It follows that Z 2D gauge,xv (Z 2D Ising,xe ) is only a function of the charge Q (0) (Q (1) ), and therefore we will label it as One may similarly compute tr ρ Q (0) ,Q (1) M S and tr ρ Q (0) ,Q (1) S C , the numerators in the expectation value for the membrane and the string operators.Let us first consider the membrane order parameter in the sector (Q 0 , Q 1 ) which we denote as ⟨M S ⟩ Q0,Q1 .One finds = ze ( v∈S x v e∈∂S z e )e βv v xv( e∋v ze) Z 2D gauge,Q (0) where ⟨W ∂S ⟩ 2D gauge,xv is the expectation value of the Wilson loop operator along the curve ∂S for the 2D Ising gauge theory with interaction x v while Area(S) is the area enclosed by the surface S. The area law follows because the 2d Ising gauge theory is confining at any non-zero temperature.We conclude that ρ Q (0) ,Q (1) has no membrane order as long as p v > 0.
On the other hand, the string order parameter = zv ( e∈C x e )z v1 z v2 e βe e xe( v∈e zv) Z 2D Ising,Q xe∈Q (1)   = ( where v 1 and v 2 label the end points of the curve C and ⟨z v1 z v2 ⟩ 2D Ising,xe is the spin-spin correlation function of the 2D Ising model with the sign of the Ising interaction determined by x e .Clearly, ⟨S C ⟩ Q0,Q1 can show longrange order at low-temperature, and following the same argument as that for the 1d cluster state, long-range order for a given sector implies that the (unnormalized) density matrix ρ Q (0) ,Q (1) is sym-LRE.For example, in the sector corresponding to all x e = 1, the long range order sets in below 2d Ising critical temperature.However, since the ordering temperature clearly depends on the sector Q (1) , to understand whether the full density ma- is sym-LRE, one needs to statistically quantify the string order as a function of the error rate.To do so, we introduce the following 'average string order parameter': Eq.( 17) is equivalent to the disorder averaged spin-spin correlation function of RBIM along the Nishimori line [35].It follows that [⟨S C ⟩ 2 ] decays exponentially as a function of |C| when p e > p c ≈ 0.109 [71].
Based on above analysis, the decohered state ρ as a function of p e and p v can be divided into four regimes using the qualitative behavior of membrane and average string orders [see Fig. 1(b)]: (i) p v = 0 and p c > p e ≥ 0: ⟨M S ⟩ Q0,Q1 = 1 (in the sector In this regime ρ must be sym-LRE.
(iv) p v > 0 and p e > p c : ⟨M S ⟩ Q0,Q1 ∼ e −Area(S) and [⟨S C ⟩ 2 ] ∼ e −|C| .This is suggestive that in this regime ρ is (statistically) sym-SRE and we provide an argument in favor of this conclusion below using an explicit convex decomposition.
We now use the CDA in Eq.( 4) with Γ = √ ρ to argue that the regime (iv) above, namely p v > 0 and p e > p c , is indeed sym-SRE.To ensure that each CDA state |ψ m ⟩ satisfies the symmetry, we choose {|m⟩ = |x v , x e ⟩}.Similar to the 1D case, we consider the singularity of 'partition function' Z m = ⟨ψ m |ψ m ⟩ as a diagnostic for transition from SRE to LRE as β is increased from zero.Since Z m = ⟨x e , x v |ρ|x e , x v ⟩, a calculation similar to that for tr ρ Q (0) ,Q (1) shows that Z m is proportional to Z 2D Ising gauge,xv Z 2D Ising,xe .One can also compute the expectation values of membrane and average string order operators with respect to |ψ m ⟩ and obtain ⟨M S ⟩ m is proportional to the expression in Eq.( 15) while [⟨S C ⟩ 2 m ] is proportional to the expression in Eq.( 17), and therefore both vanish when p v > 0 and p e > p c .
Alternatively, one may define an 'average free energy' ) with respect to |ψ m ⟩ to detect whether the ensemble {ψ m ⟩} encounters a phase transition as a function of the error rate.When β = 0, |ψ m ⟩ = |x a , x b ⟩ m is the trivial product state.On the other hand, |ψ m ⟩ becomes the 2D cluster state when β → ∞.One expects that the phase transition point can be located by the singular behavior of [log Z].Since [log Z] is proportional to the disorderaveraged free energy of the 2d RBIM along the Nishimori line, it is singular at p e ≈ 0.109.This leads to the same conclusion that {|ψ m ⟩} remains SRE in the regime (iv) above.
Interestingly, if one adopts the aforementioned CDA in regimes (ii) and (iii), then |ψ m ⟩ hosts intrinsic topological order and GHZ order, respectively.This can be argued by first considering the extreme case (p v , p e ) = (0, 0.5) in regime (ii) and (p v , p e ) = (0.5, 0) in regime (iii).When (p v , p e ) = (0, 0.5), which is an eigenstate of toric code.On the other hand, when (p v , p e ) = (0.5, 0), |ψ m ⟩ ∝ e (I + h e )|m⟩ ∝ (|x e ⟩ ⊗ e (I + x e v∈e Z v )|x v ⟩) is the 2D GHZ state.The argument based on the analyticity of the average free energy [log Z] then indicates that regimes (ii) and (iii) continue to host topological order and GHZ order, respectively.The phase diagram using the current decomposition is summarized in Fig. 1(b).
Finally, we note that order parameters similar to [⟨S C ⟩ 2 ] (Eq.( 17)) and the connections between the decohered cluster states and RBIM have also appeared in the context of preparing long-range entangled states using measurement protocols in Refs.[36,37].In paricular, our phase diagram (Fig. 1(b)) along the line p v = 0.5 is similar to the finite-time measurement induced phase transtions in Ref. [36,37].However, one crucial difference is that the mixed states in Refs.[36,37] do not respect the Z (1) 2 symmetry and therefore the corresponding transitions can not be interpreted as separability transitions protected by symmetry between a sym-LRE phase and a sym-SRE phase.Instead, the role of differ-ent sectors corresponding to the 2 symmetry is played by the flux f p = e∈p s e through a plaquette p, where s e is the measurement outcome.One may then regard the transition in Ref. [36,37] as a separability transition where in the non-trivial phase it is impossible to decompose the density matrix as a convex sum of SRE states which carry both definite Z (0) 2 charge and flux f p .Similar statements hold true for the case of 3D cluster state, which we discuss next.

C. 3d cluster state
The 3d cluster state Hamiltonian H 3d Cluster is: The Hilbert space consists of qubits residing at both the faces f and the edges e of a cubic lattice, see Fig. 1(c), or equivalently, at the edges of a cubic lattice, and the edges of its dual lattice (recall that each edge (plaquette) of the original lattice is in one-to-one correspondence with a plaquatte (edge) of the dual lattice).We assume periodic boundary conditions.This model has a symmetry whose generators are given by where c(c) specifies the cube in the lattice (dual lattice) and ∂c(∂c) denotes the faces on the boundary of c(c).Choosing Kraus operators O e/f = Z e/f with respective probabilities p e/f , using Eqs.( 5), (6), one obtains the decohered state ρ = 1 Z e −βe e h (1)   e −β f f h We now decompose ρ as a convex sum of symmetric states by writing ρ . Here, two one-form symmetry charges are labeled by = 0, 1 defined on each cube c and Q (1) = {Q e ) with S a noncontractible surface on the dual lattice.Thus, we want to compute tr ρ Similar to the cases in previous sections, we first compute the denominator tr ρ Q (1 ′ ) ,Q (1) in these expressions by inserting the complete complete basis {|x f ,e ⟩} and {|z f ,e ⟩}, and obtain tr ρ Here, Z 3D gauge,x f = ze e β f f x f ( e∈f ze) is the partition function of the 3D Ising guage theory with the sign of the interaction on each face labeled by x f , and For a system with periodic boundary condition, all possible x f ∈ Q (1 ′ ) can be reached by the transformation x f → x f e∋f σ e , σ e = ±1.Further, one may verify that Z 3D gauge,x f is invariant under the aforementioned transformation by changing the dummy variables z e → σ e z e .It follows that . Analogous statements hold true for Z 3D gauge,xe .Therefore, we write One may similarly compute tr ρ M S , and obtain the following expressions: where ⟨W ∂S ⟩ 3D gauge,x f is the expectation value of the Wilson loop operator (= z e along a closed curve) along the boundary of S for the 3d classical Ising gauge theory whose Hamiltonian is defined by the term that multiplies β f in the exponential in the second line of Eq.( 21).Since the plaquette interaction term in this Ising gauge theory depends on x f ∈ Q (1 ′ ) , similar to the discussion for 2d cluster state, we introduce an average membrane order parameter Eq.( 22 Using an argument similar to Ref. [21], and also similar to those already used in previous subsections for 1d and 2d cluster states, one can show that in regimes (i)-(iii), ρ cannot be a convex sum of symmetric pure states where membrane operators only exhibit an area-law.This suggests that these three regimes are sym-LRE.In regime (iv), ρ does not develop any average membrane orders, which strongly suggests that it is a sym-SRE state.We now use a CDA to support this expectation.
We again choose a CDA (Eq.( 4)) with Γ = √ ρ.To ensure that each |ψ m ⟩ that enters the CDA satisfies  22), one concludes that both of these order parameters vanish in regime (iv).
One may also conclude that the aforementioned decomposition in regimes (ii) and (iii) correspond to topologically ordered phases.This can be argued by first considering the extreme case (p f , p e ) = (0, 0.5) in (ii) and (p f , p e ) = (0.5, 0) in (iii).When (p f , p e ) = (0, 0.5), which is an eigenstate of the 3D toric code.The argument based on the singularity of average free energy [log Z] then indicates that in regime (ii) CDA states are topologically ordered.Similar arguments hold for regime (iii).The phase diagram using such a convex decomposition is summarized in the third column of Fig. 1(c).
It is interesting to compare our results with Ref. [72] where the Gibbs state of 3d cluster Hamiltonian was studied.The main difference between the decohered state we study, which is also takes the Gibbs form, with the state studied in Ref. [62] is that in Ref. [62], the Gibbs state is projected to a single charge sector of both 1form symmetries (and therefore possesses an exact symmetry, see comment #3 in Sec.II), which results in a phase transition as a function of temperature that is in the 3d Ising universality.In contrast, the decoherence we are considering leads only to an average (instead of exact) symmetry, and therefore, we obtain an ensemble of density matrices ρ Q (1 ′ ) ,Q (1) labeled by the symmetry charges Q (1 ′ ) , Q (1) .As discussed above, this implies that the universality class of the transition is related to the 3d random plaquette gauge model (and not 3d Ising transition).

D. 1d and 2d topological phases protected by a Z
(0) 2 symmetry Aside from the cluster states in several dimensions, Eq.( 5) and Eq.( 6) also holds for various stabilizer models realizing 1d and 2d SPT phases protected by a Z (0) 2 symmetry, which we now discuss briefly.An example in 1d is the non-trivial phase of the Kitaev chain [73]: where γ j denotes the majorana operator satisfying {γ j , γ k } = 2δ ij .It is straightforward to see that the Hamiltonian satisfies Eq.( 5) and one can choose O j as γ 2j−1 or γ 2j such that Eq.( 6) is satisfied.Therefore, under the composition of the channel E j [ρ] = (1 − p)ρ + pγ 2j−1 ργ 2j−1 , the pure state density matrix becomes the finite temperature Gibbs state with tanh β = 1 − 2p.A 2d example is the Levin-Gu state [74], where the Hamiltonian is defined on the triangular lattice and can be written as where the product runs over the six triangles ⟨pqq ′ ⟩ containing the site p.The ground state has non-trivial SPT order for the Z (0) 2 symmetry generated by U = p X p .One can verify [B p , B p ′ ] = 0 and B 2 p = 1 by straightforward algebra, and thus Eq.( 5) is satisfied.Besides, one can choose O j = Z j such that Eq.( 6) is satisfied.Therefore, under the composition of the channel E j [ρ] = (1 − p)ρ + pZ j ρZ j , the pure state density matrix becomes the finite temperature Gibbs states with tanh β = 1 − 2p.Using the CDA in Eq.( 4), one may then argue that both the decohered Kitaev chain and Levin-Gu state are sym-SRE for any non-zero p (we assume periodic boundary conditions so that there are no boundary modes).

V. SEPARABILITY TRANSITIONS FOR 2D CHIRAL TOPOLOGICAL STATES A. Setup and motivation
In this subsection, we consider subjecting chiral fermions in 2d to local decoherence.The starting pure state we consider is the ground state of a p x + ip y superconductor (p + ip SC in short), although we expect that the results will qualitatively carry over to other noninteracting chiral states.
Our motivation is as follows: it is generally believed that the 2d p + ip SC cannot be prepared from a product state using a constant-depth unitary circuit (as suggested by the fact the thermal Hall conductance of a p + ip SC is non-zero while that for a trivial, gapped paramagnet is zero).Indeed, one may think of a p + ip SC as an SPT phase protected by the conservation of fermion parity [42].Therefore, it is natural to ask what happens if one applies a quantum channel to this system where Kraus operators anticommute with the fermion parity.This is conceptually similar to our discussion in Sec.IV where we subjected a non-trivial SPT ground state to Kraus operators odd under the symmetry responsible for the existence of (pure) SPT ground state.An example of such a Kraus operator is the fermion creation/annihilation operator, and we will study this case in detail.Alternatively, one may consider subjecting p + ip ground state to decoherence with Kraus operators bilinear in fermion creation/annihilation operators.In this latter case, the fermion parity remains an exact symmetry.Based on our discussion in Sec.IV, one may expect a qualitative difference in these two cases, namely, Kraus operators linear Vs bilinear in fermion creation/annihilation operators.Let us briefly outline such a qualitative difference as suggested by field-theoretic considerations whose details are presented in Sec.V D.
Let us first consider Kraus operators linear in fermion operators.This is equivalent to bringing in ancillae fermions and entangling them with the fermions of the p + ip SC by a finite-depth unitary.Since this is a finite depth unitary operation on the enlarged Hilbert space (= ancillae + original p + ip SC), the expectation value of any observable, including non-local ones that detect chiral topological order [75,76], cannot become zero.At the same time, intuitively, the resulting mixed state for the electrons belonging to the original p + ip SC must somehow "lose its chirality" at infinitesimal coupling to the ancillae.This is indicated by treating the density matrix as a pure state in the doubled Hilbert space using C-J isomorphism, which we discuss below in detail, where we also clarify subtleties pertinent to the mapping of Kraus operators linear in fermion operators.Under the C-J map, the effect of the channel becomes a coupling bilinear in fermion operators between two chiral Ising CFTs with opposite chirality, and which, therefore, gaps out the counter-propagating chiral CFTs.The gapping out of the edge states in the double state is also manifested in the entanglement spectrum of the double state, which we also study.In particular, we show that infinitesimal decoherence leads to a gap in the entanglement spectrum.
Although working with the double state using C-J map is insightful, it does not directly tell us the nature of the decohered mixed state.One of our central aims is to understand the difference between the original pure (nondecohered) state and the decohered state not in terms of the double state obtained via the C-J map, or in terms of non-linear functions of density matrix, but directly in terms of the separability properties of the mixed state.Our main result is that the resulting mixed state can be expressed as a convex sum of non-chiral states, and in this sense, is non-chiral (i.e. it can be prepared using an ensemble of finite-depth unitaries that commute with fermion parity).
Let us next consider Kraus operators bilinear in the fermion operators.We study this problem only using the double-state formalism (i.e. the aforementioned C-J map), and obtain an effective action consisting of two counter-propagating free, chiral Majorana CFTs coupled via a four-fermion interaction.Such a Hamiltonian has already been studied in the past (see e.g.Refs.[77,78]), and we simply borrow the previous results to conclude that unlike the case for Kraus operators linear in Majo-rana operators, this system is stable against infinitesimal decoherence.Furthermore, the field-theory corresponding to the double state indicates that this system undergoes a spontaneous symmetry breaking where the gapless modes corresponding to the CFT are gapped out.The university class for this transition lies in the (supersymmetric) c = 7/10 tricritical Ising model.We discuss this below in detail in Sec.V D. We note that recently, Ref. [22] studied chiral topological phases subjected to decoherence using a generalization of strange correlator [34] to mixed states [31,32].Although Ref. [22] did not study the problem of our interest (namely, p + ip SC subjected to Kraus operators bilinear in Majorana fermions), the overall structure of the field theories obtained in Ref. [22] using strange correlator bears resemblance to the one we motivate using entanglement spectrum in Sec.V D.

B. Separability of p + ip SC subjected to fermionic Kraus operators
Our starting point is the ground state of the p + ip superconductor [47] described by the following Hamiltonian on a square lattice When t = ∆ = 1/2 and the chemical potential µ = 1, the system is in the topologically non-trivial phase.This can be diagnosed, for example, by studying the entanglement spectrum which will exhibit chiral propagating modes [79,80], or, by studying the modular commutator [81][82][83][84] which is proportional to the chiral central charge of the edge modes that appear if the system had boundaries.Relatedly, in the topological phase, the ground state cannot be written as a Slater determinant of exponentially localized Wannier single-particle states [46][47][48].
In our discussion, we assume periodic boundary conditions, so that there are no physical edge modes.We are interested in subjecting the ground state of Eq.( 25) to the composition of the following singlemajorana channel on all sites: Is the chiral nature of the ground state ρ 0 stable under the channel?More precisely, can we express the decohered density matrix as a convex sum of pure states, where each of these pure states now does not exhibit chiral states in its entanglement spectrum, and relatedly, has a vanishing modular commutator in the thermodynamic limit?Under the aforementioned channel (Eq.( 26)), the den-sity matrix will continue to remain Gaussian, and is fully determined by the covariance matrix M defined as M jk = −i tr(ρ(γ j γ k − δ jk )).As shown in Appendix B 1, under the channel in Eq.26, M evolves as E(M ) = (1 − 2p) 2 M .We write the decohered density matrix ρ as ρ(p) = e −Hρ (p) , where H ρ (p) can be determined explicitly in terms of E(M ) = (1 − 2p) 2 M as detailed in the Appendix B 1.
To write the decohered mixed state ρ as a convex sum of pure states, we consider the decomposition in Eq.( 4), and write where |0⟩ is a state with no fermions.One can analytically show at any non-zero decoherence, the realspace wavefunction for this state is a Slater determinant of localized Wannier orbitals, unlike the (undecohered) ground state of p + ip SC [46][47][48].The argument is as This expression may then be exponentiated to obtain the standard BCS-like form for As p → 0, β → ∞ (recall tanh(β) = (1 − 2p) 2 ), and one recovers the p + ip ground state where h(k) ∼ v k /u k diverges as 1/(k x + ik y ) and results in a power-law decay of Wannier orbitals [47].In contrast, at any non-infinite β (i.e.non-zero decoherence rate p), h(k) is non-infinite for any k (since |u k | 2 +|v k | 2 = 1), and therefore, the Wannier orbitals corresponding the state |ψ 0 ⟩ are exponentially localized.As an aside, this same argument also applies to the decohered 1d Kitaev chain (Sec.IV D), and more generally, to other decohered non-interacting fermionic topological superconductors.
The above argument only applies to the translationally invariant state |ψ 0 ⟩ that enters the convex decomposition in Eq.( 27).To make progress for general |ψ m ⟩, we found it more helpful to consider diagnostics which directly access the topological character (or lack thereof) of a wavefunction, and which are also more amenable to finite-size scaling.In particular, we employ the 'modular commutator' introduced in Refs.[81][82][83][84].Modular commutator is a multipartite entanglement measure that quantifies the chiral central charge for a pure state, and can be completely determined by the manybody wavefuntion [81][82][83][84].Specifically, it is defined as J ABC := i tr(ρ ABC [ln ρ AC , ln ρ BC ]) with ρ X the reduced density matrix in region X obtained from a pure state |ψ⟩ (i.e.ρ X = tr X |ψ⟩⟨ψ|).
In the absence of decoherence, the modular commutator of |ψ m ⟩ for this setup is J 0,ABC = πc/3 = π/6, as the chiral central charge c = 1/2 for the p + ip superconductor.Fig. 2 shows the modular commutator J ABC /J 0,ABC on a L × L torus as a function of L. We choose the error rate p = 0.04 and several different initial states, including |m⟩ = |0, ..., 0⟩ (uniform), |0, 1, 0, 1, .., 0, 1⟩ (staggered), and also a random bit string in the occupational number basis.We find that in all cases, J ABC vanishes in the thermodynamic limit.We also studied other values of p, and our results are again consistent with the claim that at any non-zero p, the modular commutator for the states |ψ m ⟩ vanishes in the thermodynamic limit.This provides numerical evidence that at any non-zero error rate, the  27)), as well as the modular commutator of the decohered mixed state itself.We choose error rate p = 0.04, and the following initial states |m⟩ in Eq.( 27): |m⟩ = |0, ..., 0⟩ (uniform), |0, 1, 0, 1, .., 0, 1⟩ (staggered), and |m⟩ = a random bit string in the occupational number basis.The inset shows the geometry of regions A, B, C used to define the modular commutator.We use anti-periodic boundary conditions along both directions so that the ground state is unique.
decohered mixed state can be expressed as a convex sum of states that do not have any chiral topological order, and hence must be representable as Slater determinants of single-particle localized Wannier states [46] (note that all states |ψ m ⟩ are area-law entangled).
It is important to note that in contrast to the pure states |ψ m ⟩, the modular commutator for the decohered mixed state ρ does not show any abrupt behavior change at p = 0 (dashed plot in Fig. 2).This is consistent with the fact that the arguments relating modular commutator to the chiral central charge rely on the overall state being pure [81][82][83][84], and therefore, we don't expect that modular commutator for the mixed-state ρ captures the separability transition at p = 0.This again highlights the utility of the convex decomposition of ρ into pure states.
In addition, we also numerically compute the entanglement spectrum of |ψ m ⟩ with |m⟩ the uniform product state (so that momentum along the entanglement bipartition is a good quantum number).For a chiral topological state, one expects that the edge spectrum of a physical edge will be imprinted on the entanglement spectrum of a subregion [79].Since |ψ m ⟩ is Gaussian, the entanglement spectrum is encoded in the spectrum of the matrix iM A , where M A is the restriction of the covariance matrix M to the region A in the inset of Fig. 3. Fig. 3 shows the spectrum of iM ABC (denoted as ν) as a function of the momentum k y with error rate p = 0 and p = 0.04.The geometry is again chosen as a torus, with length L x = 60, and height L y = 30.In the absence of error (p = 0), all  27)) as a function of the momentum ky for error rates p = 0 (i.e.non-decohered) and p = 0.04 (i.e.decohered).Here, we put the system on a Lx × Ly torus with Lx = 60 and Ly = 30.
states |ψ m ⟩ are projected to the p + ip ground state, and thus the spectrum shows chirality, resembling the edge states of the p + ip SC (note that we have two entanglement boundaries resulting in counter-propagating chiral states in the entanglement spectrum).After the decoherence is introduced, one finds that the chiral mode in the entanglement spectrum is gapped out, see Fig. 3.We also confirmed that the gap between the two 'bands' of the entanglement spectrum increases with the system size (not shown).Overall, both the modular commutator and the entanglement spectrum provide numerical evidence that the decohered density matrix can be written as a convex sum of free-fermion, pure states that have no chiral topological order.

C. Double-state formalism for fermions
Previous subsection focused on the single-Majorana channel that breaks the fermion parity symmetry of the initial density matrix from exact (U ρ = ρU = ρ) down to average (U † ρU = ρ).As briefly mentioned above, if one instead uses a channel where Kraus operators are bilinear in Majorana operators (so that the fermion parity remains an exact symmetry), one might expect a more interesting behavior, in particular the possibility of a phase transition between different non-trivial mixed states.One way to make progress on this case is to study appropriate non-linear functions of the density matrix [18-20, 22, 31, 85].Relatedly, one may employ the double state obtained using C-J map, which has been used in [18,20,[86][87][88][89] to study decoherence in bosonic problems.Specifically, given a density matrix ρ H acting on the Hilbert space H, one can define a state vector |ρ⟩ H⊗ H in the doubled Hilbert space H ⊗ H (with H having the same dimension as H) using the C-J map [44,45] (see footnote [90]): Here I H denotes the identity in H and |Φ⟩ H⊗ H is the product of (unnormalized) maximally entangled pairs connecting H and H, i.e., |Φ⟩ H⊗ H = ⊗ j |ϕ⟩ j,H⊗ H with |ϕ⟩ j,H⊗ H = ⊗ j ( d p=1 |p H , p H⟩ j ) and d the Hilbert space dimension on a single site.Henceforth, for notational simplicity we omit the subscript labeling the Hilbert space if there is no confusion.For bosons, it is straightforward to see that under Eq.( 30), the density matrix ρ = p,q ρ p q |p⟩⟨q| is mapped to |ρ⟩ = p,q ρ p q |p, q⟩.On the other hand, the channel α is mapped to the operator This can be derived by expressing |E[ρ]⟩ as an operator acting on |ρ⟩, i.e., |E[ρ]⟩ = N E |ρ⟩.See App.B 2 for details.However, a similar correspondence for fermions is a bit subtle.For example, naively applying Eq.( 31) to the single-majorana channel in Eq.( 26) gives where we denote η = γ as the Majorana operators in the Hilbert space H. Eq.( 32) suggests that the channel generates a real time evolution for the double state, which contradicts our intuition that the channel instead gives rise to an imaginary time evolution.Another hint that Eq.( 32) is incorrect comes from setting p = 1/2, where the relation holds.However, Eq.( 32) gives E j E j |ρ⟩ = γ j η j |ρ⟩/2, which is not equal to E j |ρ⟩.Therefore, to find the correct correspondence between E[•] and N E for fermions, one should begin with the more fundamental definition of the double state, i.e, |ρ⟩ = ρ ⊗ I|Φ⟩.
Due to the linearlity of Eq.( 30), one can consider each α ⟩, one finds K α is unchanged under Eq.( 30).On the other hand, since one can always write K † α as a function of c and c † , it suffices to consider how to express |ρc⟩ and |ρc † ⟩ as an operator applying to |ρ⟩.In Appendix B 2, we find that One can then use Eq.( 33) to derive N E given E[•].For example, for the Kraus operator given by K = γ 1 ≡ (c + c † ), one finds This implies the C-J transformed operator D. Phase transition induced by an interacting channel in a p + ip SC Being equipped with the correspondence between E[•] and N E , we now return to our discussion of decoherence induced transitions in chiral topological states of fermions.We first revisit the problem discussed in V B, and then consider a more interesting problem where the Kraus operators are bilinear in fermions so that the decohered density matrix is not Gaussian.
There are different ways to employ the double state to probe the effect of decoherence.For example, one may consider non-linear functions such as the normalization of the double state [18,20,31,85].Here we will motivate the entanglement spectrum of a state obtained from the double state |ρ⟩ (after space-time rotation) as a probe of the decoherence-induced phase transitions.
To begin with, consider the normalization of the double state If the bulk action describing |ρ 0 ⟩ = |Ψ 0 , Ψ * 0 ⟩ is rotationally invariant, one can map ⟨ρ|ρ⟩ to the path integral of the (1 + 1)-D boundary fields following Ref.[20]: Here, ψ L and ψ * L denote the low-energy field variables in H and H, respectively.S 0,L is the partition function on the left side of the spatial interface x = 0 − (the meaning of (ψ R , ψ * R ) and S 0,R are similar with left ↔ right).S int describes the effect of the channel E † E and has two contributions: Here, S 1 denotes the action that exists even in the absences of decoherence.In particular, S 1 strongly couples the fields in the absence of decoherence.On the other hand, S E describes the action that merely comes from the decoherence and vanishes when the error rate p = 0. We note that a similar field theory has also been discussed in Ref. [91] in a different context when evaluating the system-environment entanglement in the (1+1)D sys-FIG.4. The spectrum of iML for the double state |ρ⟩, where ML is the restriction of the covariance matrix M to the region L, as a function of the momentum ky for different error rates p.Here, we put the system on a cylinder with circumference Lx = 60, and height Ly = 16.
tem.In general, the exact form of S E involves four fields (ψ L , ψ * L , ψ R , ψ * R ) and may be schematically captured by the following Hamiltonian: where H 1 strongly couples the L and R fields.One may then consider the reduced density matrix for L fields that is obtained after tracing out the R fields.One expects [92,93] that the corresponding entanglement Hamiltonian (= logarithm of the reduced density matrix) will essentially correspond to H 0,L + H E,L .Working with entanglement Hamiltonian has the advantage that the number of fields one needs to keep track of are now halved.Similar simplification occurs if one considers the fidelity tr(ρ d ρ 0 ) between the decohered density matrix ρ d and the non-decohered density matrix ρ 0 , see Ref. [22].Since we are now working only with the L fields, in the following we will omit the subscript L for notational simplicity.
As an example, let us first revisit the case of p + ip superconductor perturbed by a channel that is linear in Majorana fermions (Sec.V B).Recall that here the Kraus map corresponds to the composition of the following map on all sites: E x [ρ] = (1 − p)ρ + pγ x ργ x .Based on our discussion above on the C-J map for fermions, this translates to a term of the form H E = ig dy γ η where p ∼ g, and γ, η respectively denote the fields corresponding to H and H of the L fields.In the absence of any decoherence, the spatial boundary of the p + ip superconductor has a simple description in terms of a chiral Majorana fermion.The entanglement Hamiltonian in the doubled Hilbert space then corresponds to stacking the boundary of p + ip and p − ip superconductors, and is given by H 0 = i dy(γ∂ y γ − η∂ y η).Therefore, one expects that the entanglement Hamiltonian for the L fields in the presence of decoherence takes the form: The counter-propagating edge modes are gapped out for any non-zero g (∝ p), in line with our earlier discussion where we provided evidence that at any non-zero p, the density matrix can be written as a convex sum of pure states that are SRE.The gapping out of the edge modes can also be seen by numerically evaluating the entanglement spectrum of the double state obtained via C-J map.Fig. 4 shows the spectrum of iM L (denoted as ν) as a function of the momentum k y with different error rate p.
Here, we put the system on the cylinder with circumference L x = 60, and height L y = 16.In the absence of error (p = 0), there are two counter-propagating modes, resembling the edge states of the initial double state |ρ 0 ⟩.After the decoherence is introduced, one can clearly see from Fig. 4 that these counter-propagating modes are gapped out for arbitrary small error rate.Note that we did not perform any space-time rotation to obtain Fig. 4.This is suggestive that the entanglement Hamiltonian of the double state |ρ⟩ may already have the same qualitative behavior as the one obtained after space-time rotation.We leave further investigation of this point to the future.Let us return to the problem of our main interest in this subsection, namely that of Kraus operators that commute with the fermion parity operator.The simplest possibility is the composition of the following Kraus map on all nearest-neighbor bonds ⟨x, y⟩ of the square lattice: The interaction term H E induced by such a Kraus map in the double state should respect the following Z 2 × Z 2 symmetries: γ → −γ and/or η → −η.Since Majorana fermions square to identity, the simplest term that is bilinear in both γ and η and respects all the symmetries involves derivatives: where g ∝ p.Therefore, the full entanglement Hamiltonian for the L fields in the presence of decoherence is given by: This field theory has been studied earlier in Refs.[77,78].At a particular g = g c , the system undergoes a phase transition in the tricitical Ising university class with central charge c = 7/10.For g < g c , the interaction term is irrelevant, while above g c , the system spontaneously breaks the Z 2 × Z 2 symmetry down to the diagonal Z 2 symmetry.Physically, this means that the exact fermionparity symmetry (i.e.U ρ = ρ where U is the generator of the fermion parity), has been spontaneously broken down to an average symmetry (i.e.U ρU † = ρ).We note that a class of 2d chiral topological phases subjected to decoherence with fermion-bilinear Kraus operators was also studied in Ref. [22].One notable difference between Ref. [22] and our problem is that in the examples considered in Ref. [22], the decoherence always reduces the effective central charge of the action corresponding to the double state.In contrast, in our problem, the effective central charge increases from c = 1/2 to c = 7/10.It is interesting to contemplate the implications of the above phase transition in terms of the separability properties of the original mixed state ρ (instead of the double state |ρ⟩).We conjecture that for p ≪ p c , there exists no decomposition of the density matrix as a convex sum of area-law entangled pure states without any chirality, while p > p c the density matrix is expressible as a convex sum of area-law entangled pure states with GHZ-like entanglement (due to spontaneous breaking of fermion parity).Similar to the case of intrinsic topological orders subjected to local decoherence [17,18,20,84], we anticipate that the universality class as well as the location of the critical point obtained from the double-state formalism will differ from that of the 'intrinsic' mixed-state transition for the density matrix, e.g., when viewed from the perspective of separability.We do not know the universality for the latter transition and we leave it as an open question.

VI. SEPARABILITY TRANSITION IN GIBBS STATES OF NLTS HAMILTONIAN
In this section we will consider an exotic separability transition in a Gibbs state relevant to certain quantum codes.Although this transition does not require any symmetry, which has been a main ingredient in the rest of this work, the argument below to deduce the existence of a separability transition is broadly similar in spirit to that in Secs.III and IV.
Recently, there has been discovery of 'good LDPC codes' where the code distance as well as the number of logical qubits scale with the total number of qubits [51][52][53].Moreover, Ref. [50] showed that the construction of a good LDPC code in Ref. [51] satisfies Freedman-Hastings' No Low-Energy Trivial States (NLTS) conjecture [49] which, when satisfied by a Hamiltonian, means that any state |ψ⟩ with energy density less than a nonzero value e c can not be prepared by a constant depth unitary circuit (the energy density e of a state |ψ⟩ is defined as e = lim N →∞ (⟨ψ|H|ψ⟩ − E 0 ) /N where E 0 is the groundspace energy of H).Here we ask: does the Gibbs state of a NLTS Hamiltonian show a separability transition at a non-zero temperature?That is, does there exist a T c > 0 so that for T < T c , the Gibbs state can not be written as a convex sum of SRE pure states?Firstly, we note Ref. [50] already proved that any mixed state whose energy density is less than an e c > 0 cannot be purified to a pure SRE state by a short-depth channel, i.e., it cannot be prepared by first enlarging the Hilbert space to include ancillae, which are initially all in a product state, followed by a finite-depth unitary that entangles the 'system' qubits (which are also initially in a product state) with the ancillae qubits, and eventually integrating out the ancillae.However, as already discussed in Sec.II, inability to purify to an SRE state via a short-depth channel does not imply that a mixed state is SRE using our definition (= expressibility of a mixed state as a convex sum of SRE pure states).Just to briefly re-iterate the example discussed in Sec.II that illustrates these two different notions of mixed-state entanglement (see comment #4 in Sec.II), any Gibbs state that exhibits spontaneous symmetry breaking (which therefore has long-range correlations for the operator corresponding to the order parameter), can not be purified to an SRE pure state via a short-depth channel.Therefore such a mixed state will be SRE using our definition and LRE using the definition in Refs.[50,54].Here we provide a simple argument that the Gibbs state of an NLTS staisfying Hamiltonian shows a separability transition at a non-zero temperature.
Let us assume that the Gibbs state of an NLTS satisfying Hamiltonian H can be expressed as a convex sum of SRE pure states for any temperature T > 0, i.e., ρ(T ) = e −H/T /Z = i p i |ψ i ⟩⟨ψ i | where each |ψ i ⟩ can be prepared via a unitary whose depth is independent of the number of qubits N .For simplicity of notation, we set the groundspace energy E 0 to zero (this can always be achieved by adding a constant N c to the Hamiltonian, where c is a constant).We will show that this assumption leads to a contradiction.Since all pure states |ψ i ⟩ are SRE, by NLTS condition, they must all satisfy ⟨ψ i |H|ψ i ⟩/N > e c as N → ∞.Therefore tr(ρ(T )H)/N = i p i ⟨ψ i |H|ψ i ⟩/N > i p i e c = e c .This implies that if the Gibbs state can be expressed as a convex sum of SRE pure states, then its energy density is non-zero.However, non-zero energy density necessarily implies non-zero temperature.This is equivalent to showing that as T → 0, tr(ρ(T )H)/N → 0. This is indeed the case because as T → 0, tr(ρ(T )H)/N ≈ E 1 e −E1/T /N which indeed vanishes as T → 0 (E 1 denotes the energy of the first-excited state, which is a constant independent of N since the LDPC code Hamiltonian under discussion is a sum of commuting projectors).Therefore, if we assume that the Gibbs state is separable for all non-zero temperatures, we arrive at a contradiction.Hence the Gibbs state must be long-range entangled upto a nonzero temperature T .It seems reasonable to assume that at sufficiently high temperature, the Gibbs state is SRE.Therefore, one expects a separability transition at some temperature T c that satisfies 0 < T c < ∞.
It is important to note that the above-argued separability transition does not necessarily imply that the Gibbs state has a thermodynamic phase transition, i.e., it need not be accompanied by a singularity of the partition function.

VII. SOME CONNECTIONS BETWEEN SEPARABILITY AND OTHER MEASURES OF MIXED-STATE COMPLEXITY
In this section, we comment on some connections among the separability criteria, purifications, double states, and strange correlators.

Connections among separability, purification, and double state
In Sec.V C, we used the double-state formalism to probe decoherence-induced transitions.However, the connection between ρ being sym-SRE and the double state |ρ⟩ being trivial remains unclear.In this subsection, we will attempt to bridge the gap between them using purification of the mixed state.
We first recall the idea of purification: given a mixed state ρ H in the Hilbert space H, there exists a purification in the double Hilbert space H⊗ H with H having the same dimension as H: where |Φ⟩ H⊗ H is a maximally entangled state between H and H.It is straightforward to see that tr H(|ρ 1/2 ⟩⟨ρ 1/2 |) = ρ H . Besides, we note that Eq.( 43) is somtimes called 'standard purification' [57], and all possible purifications are equivalent up to an isometry applied merely in H.If one uses Eq.( 2) as a definition of an SRE mixed state ρ, then |ρ 1/2 ⟩ being an SRE implies that ρ is SRE.However, it is not obvious to us how to show that this implies that ρ can be written a convex sum of SRE states (Eq.( 1)).Instead, we are only able to show that if |ρ with m = (m 1 , m 2 , ...) and n = (n 1 , n 2 , ...).Since |Φ m,n ⟩ are obtained by applying local unitary in H to a maximally entangled state, they are all also maximally entangled.Now, we use the same idea as we used to de-fine CDA states (Eq.( 4)) by writing ρ ⊗ I/dim(H) as In the second line, we use the property that j (Z j, H) mj (X j, H) nj and ρ 1/2 commute, as they act on different Hilbert space.Since |ρ 1/2 m,n ⟩ is related to |ρ 1/2 ⟩ by a unitary acting solely in H , if |ρ 1/2 ⟩ is SRE, then so is |ρ 1/2 m,n ⟩.Therefore, if there exists an SRE purification for ρ (Eq.( 43)), then ρ ⊗ I/dim( H) can be written as a convex sum of SRE pure states (Eq.( 45)).
However, we emphasize that the converse is not true: if |ρ 1/2 ⟩ is not an SRE, it does not rule out the possibility that the mixed state ρ is still an SRE.This can be most easily seen by considering the following counterexample that also appeared in Sec.II.Let ρ be the convex sum of two product state |0⟩ N and |1⟩ N , i.e., It follows that the purified state is the GHZ state: which is clearly long-range entangled.This implies that |ρ 1/2 ⟩ being trivial is a sufficient but not necessary condition for ρ ⊗ I/dim( H) being trivial.
The advantage of studying ρ using its purification is obvious: instead of finding the decomposition in Eq.( 4), one only needs to deal with a single pure state |ρ 1/2 ⟩.However, it is in general difficult to compute |ρ 1/2 ⟩, as taking a square root of ρ is non-trivial if H ρ = − log(ρ) doesn't admit a simple compact form.An alternative is to consider the double state |ρ⟩ = ρ ⊗ I|Φ⟩ in Eq. (30) (note that if the original density is pure, i.e., ρ 2 = ρ, then the double state |ρ⟩ is equivalent to the purified state |ρ 1/2 ⟩).Heuristically, since the coefficient in front of H ρ for |ρ⟩ is higher than the coefficient for |ρ 1/2 ⟩, we expect that if |ρ⟩ is SRE, then |ρ 1/2 ⟩ is SRE as well, but we do not know how to prove this.This is consistent with the result in Ref. [20], where the critical error rate for |ρ⟩ being trivial is higher than the error rate that the topological entanglement negativity drops to zero, and also consistent with the results in Ref. [17] on the error threshold for separability for topologically ordered mixed states.
Finally, as a small application of the state |ρ 1/2 ⟩, let us briefly consider the 2d toric code under the phase-flip channel.We will now argue that under such a decoherence process, |ρ 1/2 ⟩ undergoes a transition from being double topologically ordered to single topologically ordered precisely when the decoherence rate corresponds to the multicritical Nishimori point in RBIM, which further coincides with the separability transition for ρ [17], as well as the decoding transition and other metrics of mixed state entanglement [18][19][20]61].The calculation is almost identical to that for a closely related pure state studied in Ref. [17], and therefore we will only sketch the calculation, and refer the interested readers to Ref. [17] for details.The density matrix of the 2d toric code (on a torus of size L×L) under the phase-flip channel can be written as [17,18]: p } defines the flux configuration through all plaquettes p of the square lattice via the equation Q (1) p = e∈∂p x e ( e∈∂p means that the product is taken around the edges belonging to the plaquette p).Further, we choose L x = e∈ℓ,e∥x x e = L y = e∈ℓ,e∥ŷ x e = 1 with ℓ a non-contractible loop along x/ŷ direction.Therefore, one finds |ρ 1/2 ⟩ ∝ Q (1)   Z 2d Ising,Q (1) |Ω Q (1) ΩQ (1) ⟩.The Higgs transition for this pure state can be detected by the t'Hooft loop T l = e∈ l Z e Ze where l denotes a homologically non-contractible loop on the dual lattice.
Specifically, a straightforward calculation shows that ⟨ρ 1/2 |T l|ρ 1/2 ⟩/⟨ρ where ∆F (Q l ) is the free energy cost of inserting a domain wall of size l ∼ L. It follows that ⟨T l⟩ = 0 (a non-zero constant) in the thermodynamic limit when p < p c = p 2d RBIM (p > p c ), indicating the transition of |ρ 1/2 ⟩ from being double topologically ordered to single topologically ordered.

Connections between convex decomposition and strange correlator
In Sec.IV, we studied separability transitions for cluster state SPTs in various dimensions using CDA (Eq.( 4)) with the initial basis {|m⟩} as product states satisfying the corresponding symmetry of the cluster state SPT (which was Pauli-X basis in all the cases we considered).Fortuitously, as we discussed, the threshold for the CDA states being sym-SRE exactly corresponded to the error rate beyond which ρ must be sym-LRE using general arguments, indicating that our choice of CDA is optimal.
Intriguingly, the symmetric product state basis to generate CDA has an apparently close connection with the strange correlator [34], which was originally devised as a diagnosis of the SPT pure states and has recently been used to probe the non-trivial SPT mixed states [31,32].To see the connection between them, we briefly review the original strange correlator for SPT pure states and two types of strange correlator introduced in Ref. [31].Choosing |m⟩ as the disordered product state respecting the symmetry group G, the strange correlator for a pure state |ψ⟩ is defined as [34] where O is some operator that transforms non-trivially under G.The basic idea of strange correlator is that the temporal edge of SPT pure state (when the many-body wavefunction is expressed as an imaginary time path integral) mimics its spatial edge.Since at least 2d SPT possess nontrivial spatial edge states (in 3d, there also exists a possibility of boundary topological order), one may also use the temporal correlation defined in Eq.( 48) to probe non-trivial SPT.To generalize the strange correlator from pure states to mixed states, two types of strange correlator were introduced in Ref. [31].The type-I strange correlator is defined as In the pure state limit ρ = |ψ⟩⟨ψ|, the type-I strange correlator reduces to Eq. (48).Therefore, in the case of subjecting local decoherence to an SPT pure state, C I m can be intuitively regarded as asking whether the local decoherence destroys the temporal edge states.However, it has been shown in Ref. [31] that the type-I strange correlator is unable to detect the average SPT order mentioned in Ref. [28].Instead, it was argued that the non-triviality of such an SPT order should be detected by the type-II strange correlator, defined as On the other side, the necessary condition for the mixed state to be non-trivial using separability criteria is the non-triviality of CDA states |ψ m ⟩, which may be probed by several physical observables S as discussed in Sec.IV: Comparing Eq.( 49), Eq.( 50) and Eq.( 51), one finds that the denominator is always the fidelity between a symmetric product state and the mixed state of interest: For the numerator, Eq.( 51) involves inserting an operator between ⟨m|ρ 1/2 and ρ 1/2 |m⟩, while the strange correlator involves inserting an operator between ⟨m|ρ and |m⟩.

VIII. SUMMARY AND DISCUSSION
In this work we explored the interplay of complexity and symmetry for many-body mixed states.Specifically, we asked whether a given mixed state can be expressed as a convex sum of symmetric short-ranged entangled pure states, which we took as a definition of an SRE mixed state subject to a given symmetry (a 'sym-SRE' mixed state, Sec.II).Our primary aim was to identify 'many-body separability transitions' as a function of an appropriate tuning parameter (e.g decoherence rate or temperature) across which the nature of the mixed state changes qualitatively -on one side of transition the mixed state is sym-SRE, and on the other side it is sym-LRE (= not sym-SRE).Analogous phase diagrams for intrinsic topological orders subject to local decoherence [18][19][20][21][22] were recently studied in Ref. [17].Our general approach was to first seek constraints that imply that a mixed state is necessarily long-range entangled, and absent such constraints, we developed tools to find the regime where a mixed state can be shown to be sym-SRE.One of the tools that allowed us to make progress was that local decoherence converts ground states of several SPTs, e.g.cluster states in various dimensions, to a Gibbs state.
In the context of SPTs subjected to local decoherence, we focussed on cluster states in various dimensions and obtained their 'separability phase diagram' as shown in Fig. 1.As evident from the figure, the phase diagram gets progressively richer as one moves up in spatial dimensionality.The paths solely along the x and y axes in these phase diagrams correspond to the special case of 'average SPT' mixed states where one of the symme-tries is exact while the other is average [28][29][30][31][32].It is crucial to note that although the decohered mixed state takes a Gibbs form, the corresponding partition function is not singular at any non-zero temperature for any of these cluster states.This is because any local channel can be Stinespring dilated as a local unitary circuit in the enlarged system, and thus any physical observables tr(ρO) with O acting on a large but finite region must be a smooth function of the error rate [18,19,91].Therefore, the different phases in Fig. 1 arise only because we are requiring that the density matrix be expressible as a convex sum of pure, symmetric states.As a consequence, these transitions are conceptually distinct from thermal phase transitions, and are more akin to 'complexity phase transitions' for the mixed state, when a symmetry is enforced.We briefly discussed relation with other approaches to classifying mixed-state SPTs [28][29][30][31][32] in Sec.VII.
It is also interesting to contrast the symmetry-enforced separability transitions in decohered 2d and 3d cluster states with decoherence induced separability transitions in 2d and 3d toric codes, studied in Ref. [17].In both cases, one finds the appearance of same statistical mechanics models (e.g.RBIM in 2d).This similarity can be traced to the fact that the ground state of toric codes can be obtained from the ground state of the cluster states by performing appropriate projective measurements [94][95][96][97], along with the equivalence between local and thermal decoherence for cluster states (this statement holds true also for the fractonic X-cube model [98] and its parent cluster state [96]).
We also studied non-stabilizer topological states subjected to local decoherence.In particular, for free fermion chiral states corresponding to a p+ip superconductor, we argued that if the quantum channel responsible for decoherence breaks the fermion parity, the resulting Gibbs state can be expressed as a convex sum of non-chiral states, and is therefore SRE at any non-zero decoherence rate (Sec.V).We also studied a case where the channel respects the fermion parity and identified a mixed state phase transition as a function of the decoherence rate using the double-state formalism.This transition can be thought of as corresponding to spontaneous breaking of the fermion parity, and as far as we know, does not have a pure-state counterpart.Intuitively, in a purestate context, breaking fermion parity spontaneously essentially requires assigning a non-zero expectation value to fermionic operators, which is unphysical.In contrast, in the context of a mixed state, breaking fermion parity spontaneously means that the environment can exchange fermions with the system 'spontaneously', which is not unphysical (in the double-state formalism, this corresponds to assigning non-zero expectation value to the bosonic order-parameter η γ where η and γ respectively denote the fields corresponding to bra and ket Hilbert spaces).
We also analyzed separability transitions in the Gibbs state of the quantum Ising model and argued that the Gibbs state is SRE at any non-zero temperature, and sym-SRE only for T > T c , where T c is the critical temperature corresponding to the spontaneous symmetry breaking (Sec.III).We expect similar results to hold for other models whose Gibbs state shows a spontaneous breaking of zero-form discrete symmetry.
Finally, in Sec.VI, we provided a short argument that the Gibbs states of Hamiltonians that satisfy NLTS (no low-energy state) condition [49] must exhibit a separability transition at a non-zero temperature.
In the rest of this section, we discuss various aspects of our results and motivate questions for further exploration.
1. SPT and chiral states: The technique we used to study phase diagrams of various cluster states relied on the fact the quantum channel resulted in Gibbs states (Eqs.(5),( 6)).It is not obvious how to generalize it to other SPT states.On that note, the following Z N generalization may be helpful to study Z N cluster states and topological orders produced by partial measurement of such states.Let us consider a commuting projector Hamiltonian of the form H = P i , where P i are projector (P 2 i = 1) written as Let us now introduce the following set of Kraus operators on each site i, K one may verify that the application of this channel on all sites again results in a Gibbs state for H.
It might also be interesting to study 'intrinsically mixed' SPT states introduced in Refs.[28,29] from the point of view of separability.These are SPT states that can exist only in the presence of decoherence.Conversely, it will be interesting to understand our results on nontrivial mixed SPTs protected by higher-form symmetries, such as 2d and 3d cluster state, using the techniques in Refs.[28,29] which primarily focussed on zero-form symmetry SPTs.
In the context of chiral states, we studied a phase transition driven by a channel where the Kraus operators were Majorana bilinears (Sec.V C).We analyzed this problem only using the double-state formalism.As suggested by the problem of decoherence in toric code, the double state is likely to overestimate the threshold for the actual transition, and it will be interesting to find a description of the aforementioned transition in p + ip SC directly in terms of the separability properties of the mixed state.
One important subtlety we would like to point out is that we assumed periodic boundary condition in our discussion of the SPT and chiral states.If instead one considers open boundaries such that the boundaries do not break the symmetry responsible for nontrivial SPT/chiral topological order, then the pure (nondecohered) state is always LRE, e.g., due to propagating edge modes or topolgical order at the boundary.In the presence of decoherence, our naive expectation is that the resulting mixed state is not sym-SRE, even if the decoherence breaks the symmetry from exact down to average.It will be interesting to study this aspect in the future.
2. Symmetry broken states: The first example we discussed, primarily to illustrate the distinction between SRE and sym-SRE states, was the Gibbs state of the transverse-field Ising model in any dimension (Sec.III).We discussed an explicit decomposition of this state at a non-zero temperature as a convex sum of pure states which we argued are SRE at any non-zero temperature.This conclusion is consistent with numerical results on Renyi negativity [26], and mean-field arguments [25,27].On the other hand, if one imposes the Ising symmetry on the pure states into which the Gibbs state is being decomposed, we adopted an argument from Ref. [21] to show that these pure states must be long-range entangled for T < T c .This implies that the Gibbs state is sym-LRE for T ≤ T c .In contrast, for T > T c we provided an explicit sym-SRE decomposition of the Gibbs state.The basic idea of the argument is to write e −βH as ϕ e −βH/2 |ϕ⟩⟨ϕ|e −βH/2 where {ϕ} are chosen as complete set of states in z(x) basis if one wants to expand the Gibbs state as a sum of (sym-)SRE pure states.
There are several open questions along this direction.Firstly, the argument we provided for the aforementioned pure states being (sym-)SRE is not mathematically rigorous.To explicitly show that a state is SRE, one needs to construct finite-depth circuit that prepares it starting from a product state.We only provided arguments in the continuum limit that the pure states under consideration have short-range correlations.It will be worthwhile to study the entanglement structure of the pure states we claimed to be SRE using numerical methods (e.g.quantum Monte Carlo), or using a detailed field-theoretic analysis.Secondly, as we discussed, the transverse field Ising model in d ≥ 2 must exhibit a separability transition from a sym-LRE to sym-SRE as a function of temperature.It will be interesting to study the symmetryresolved negativity [99] to quantify the nature of longrange entanglement across this transition.Finally, our arguments only apply to Gibbs states that break a discrete symmetry spontaneously, and it will be interesting to consider generalization to systems with spontaneously broken continuous symmetries that host Goldstone modes at a non-zero temperature.
3. Experimental and numerical implications: It is interesting to contemplate experimental implications of a symmetry-enforced separability transition.One perspective is that symmetry resolved versions of mixedstate entanglement measures such as entanglement negativity or entanglement of formation, that are specifically designed to quantify the lack of separability, would likely experience a singularity across such a phase transition.For example, for the Gibbs state ρ of the transverse-field Ising model (Sec.III), one can in principle prepare the states ρ ± = P ± ρ where P ± are the projectors onto the even and odd sectors of the Ising symmetry.This can be done, e.g., by entangling an ancilla qubit with the system qubits sequentially using CNOT gates, and by measuring the ancilla qubit at the end.As discussed in Sec.III, the resulting mixed state (i.e.ρ + or ρ − , depending on the outcome of the measurement on the ancilla) will show long-range mixed-state entanglement for T < T c , in contrast to the original density matrix ρ, which will be short-ranged entangled for any T > 0. The long-range entanglement of ρ ± can in principle be quantified experimentally using the Renyi negativity [100].
One may also imagine a very patient, gedanken experimentalist who has access to local unitary gates with a finite fidelity, so that they have the ability to prepare only an ensemble of SRE pure states (i.e.pure states preparable with a constant depth unitary).If so, then a separability transition from an SRE to an LRE mixed state is equivalent to the transition from success to failure in preparing the ensemble corresponding to the mixed state.One may similarly characterize a transition from a sym-SRE to a sym-LRE by putting symmetry constraints on the local gates that form the circuit.
Perhaps a more practical implication of our results is that it may allow efficient classical simulation of a class of mixed states.For example, in the context of Gibbs state of the quantum Ising model, we argued that it admits a convex decomposition in terms of SRE pure states at any non-zero temperature if one does not impose any symmetry constraint on the pure states.Since SRE states are typically easier to study, such a representation facilitates the task of simulating the corresponding mixed state.In contrast, if one tries to prepare the Gibbs state of the quantum Ising model starting with a product state (assisted with ancillae), then long-range correlations below T c imply that one necessarily requires a deep quantum circuit [101].We note that the decompositions we study generically involve an exponentially large number of pure states, which may lead to another difficulty in preparation.We can imagine at least two distinct ways to address this.Firstly, if a mixed state is SRE, and does not contain any classical long-range correlations, then it is reasonable to expect that it can be purified into an SRE pure state (using ancilla degrees of freedom).This equivalence between SRE mixed state and SRE purification is discussed for Gibbs states in Ref. [101], and suggested more generally in Ref. [2] although we are not aware of an explicit proof/construction in the non-Gibbs case.If one can indeed find an SRE purification, then an SRE mixed state can be prepared by a finite depth unitary circuit acting on system and ancilla degrees of freedom.An alternative route that is more generally available is to sample from the ensemble of SRE states that enter a given decomposition (assuming that the density matrix is SRE) using Monte Carlo algorithms, instead of preparing each and every SRE state that enters the SRE decomposition.While the sampling task may still be generally hard even in classical mechanics [102,103], the decomposition we have provided clearly simplifies the "quantum hardness" of simulating a mixed state, analogous to the METTS algorithm [58].
On a different note, one way to prepare an ensemble of pure states that may show a mixed-state separability phase transition is via a judicious combination of unitaries and measurements [21,36,37,[94][95][96][97]104].For example, Refs.[36,37] provide a construction of mixed states that are closely related to the mixed states discussed in Sec.II, and which have also been implemented in a recent experiment (Ref.[105]).It will be interesting to design experiments that probe the phase diagram in Fig. 1 using similar ideas, although we suspect it may be comparatively more challenging as it requires measuring non-local observables supplemented with appropriate decoding scheme [36].
note that one can also pair up two majorana fermions to get regular fermions through α j = (ξ 2j−1 + iξ 2j )/2, and the density matrix take the form ρ ∝ e −2β j α † j αj mentioned in Sec.V B.
To proceed, we note that Eq.(B1) will map a Gaussion state to a Gaussian state.A Gaussian fermionic state ρ, whether pure or mixed, is fully specified by the covariance matrix M , defined as M jk = −i tr(ρ(γ j γ k − δ jk )).(B4) To see how one can use Eq.(B4) to deduce the resulting decohered mixed state, let us first explicitly write down the relation between a general density matrix ρ and its covariance matrix M .Let us write the density matrix as ρ ∝ e −Hρ .Since ρ is Gaussian, H ρ can be written as a sum of Majorana bilinears: where K is a 2N × 2N antisymmetric matrix and we denote the number of Majorana modes as 2N .To see how the matrix M is related to the matrix K, we begin by block diagonalizing K: Using −i tr(ρξ 2j−1 ξ 2j ) = tanh((K d ) j,j ) and the relation between ξ j and γ j , the covariance matrix can then be where tanh β = (1 − 2p) 2 and ξ j = (O T ) jk γ k are the majorana operators that block diagonalize the original Hamiltonian H.

Double-state formalism for fermions
In this section, we derive the double-state formalism for fermions.As a warmup, we will first derive how the bosonic density matrices and channels are mapped to pure states and operators, respectively, under We will later derive the similar correspondence for fermions using Grassmann algebra.
For a general bosonic density matrix ρ = p,q ρ p q |p⟩⟨q|, one can compute |ρ⟩ using Eq.(B10): |ρ⟩ = p,q ρ p q |p⟩⟨q| r |r⟩ ⊗ |r⟩ = p,q,r ρ p q δ q r |p⟩ ⊗ |r⟩ = p,q ρ p q |p, q⟩, (B11) which is tantamount to flipping the bra vector of ρ in the original Hilbert space H to the ket vector in H.This intuition can be visualized by expressing Eq.(B10) using tensor network [see Fig. 5(a)].In this sense, one can regard the maximally entangled states as a tool to transform the bra(ket) to ket(bra) spaces, and such a trick is called Choi-Jamio lkowski (C-J) isomorphism.We note that C-J isomorphism was originally used to map quantum channels E (superoperators) to quantum states σ (density matrices) [44,45]:

(B13)
Therefore, one finds On the other hand, the corresponding operator σ E for the channel E is evaluated using Eq.(B12), and one finds where |K α ⟩ = K α ⊗ I|Φ⟩.It is then obvious that N E and σ E are different.For example, if E[•] = K(•)K † only consists of a single Kraus operator K, σ E is necessarily a projector while N E is not in general.
Now we turn to the C-J map for fermions, i.e., the analog of Eq.(B13) for fermions.Note that the derivation of N E in Eq.(B13) required the insertion of a complete basis.For fermions, this can be achieved by using Grassmann algebra.To build intuition, we consider the mixed state density matrix ρ with a single fermionic mode with creation/annihilation operators c † /c (which act on the Hilbert space H in our notation), i.e., ρ = ρ(c, c † ).The maximally entangled state in the double Hilbert space for fermions can then be defined as Interestingly, treating C-J map as a general way to transform the bra(ket) space to ket(bra) space leads to other useful applications for fermionic problems.For example, the fermionic transpose can be naturally defined using C-J map, and we find that this definition is consistent with the fermionic time-reversal, which was proposed in Ref. [108] to resolve the issue that the conventional definition of the fermionic-transpose fails to capture the entanglement negativity due to the formation of the edge Majorana fermions.Specifically, the fermionic transpose can be defined as follows: where {|m⟩} is an arbitrary complete basis for fermions.
One may show that this definition makes sense by expressing Eq.(B23) in terms of tensor network (see Fig.

FIG. 1 .
FIG.1.Cluster states under decoherence in (a) 1d, (b) 2d, and (c) 3d.The first column depicts the Hamiltonian of cluster states.The second column divides the decohered mixed state as a function of error rates into several regimes that have qualitatively different behaviors.The white regions (region (iv)) in the three phase diagrams denote phases where the mixed state is 'sym-SRE' ('trivial'), i.e., it is expressible as a convex sum of symmetric, short-ranged entangled pure states.In contrast, the colored regions or lines (regions (i), (ii), (iii)) denote phases where such a decomposition is not possible ('sym-LRE').There can be phase transitions from one kind of sym-LRE phase to a different kind of sym-LRE phase as depicted by different colors.The phase diagram is obtained by calculating objects of the form [⟨O⟩ 2 ] = Q P (Q) (⟨O⟩Q)2 where O corresponds to an appropriate observable that characterizes symmetry-enforced long-range entanglement, and P (Q) is the probability for obtaining the symmetry charge q.The pc ≈ 0.109 in the second row corresponds to the ferromagnetic to paramagnetic phase transition in the 2d random-bond Ising model along the Nishimori line, while pc ≈ 0.029 in third row corresponds to the critical point in the 3d random plaquette gauge model along the Nishimori line.The third column shows the phase diagram obtained by expressing ρ as a convex sum of symmetric states, where each symmetric state |ψm⟩ = ρ 1/2 |m⟩ with |m⟩ the product state in Pauli-X basis.See main text for more details.
carry no charge under the symmetry U b as |ψ m ⟩ is an eigenstate of U b .Therefore, there are two disjoint possibilities for the operators O L and O R : they are either both charged under the symmetry U b or neither of them are charged under U b .If neither of them are charged under U b , x b ⟩ m is a product state.To check whether |ψ m ⟩ remains SRE for any non-infinite β a and β b , let us consider the 'partition function with respect to |ψ m ⟩' Z m (β a , β b ) = ⟨ψ m |ψ m ⟩ (12) as a function of β.As β a , β b are increased from zero, if the state |ψ m ⟩ becomes long-range entangled, one expects that it will lead to a non-analytic behavior of Z m (β) as a function of β a , β b .The calculation for Z m (β) = ⟨x a , x b |ρ|x a , x b ⟩ is quite similar to the one for tr(ρ Qa,Q b ) detailed in Appendix A, and one finds that Z m (β) is proportional to the product of two partition functions for the 1d classical Ising model at inverse temperatures β a , β b .Therefore, we expect that |ψ m (β)⟩ remains an SRE state as long as both β a , β b < ∞, which confirms our expectation that ρ is sym-SRE for non-infinite β a , β b (i.e.p a , p b > 0).

Fig. 1 (
Fig.1(b).The Hamiltonian has both a zero-form symmetry Z (0) 2 , and a one-form symmetry Z (1) 2 on each cube c in the dual lattice.Now, let's focus on the physical observables that characterize each sector.These are the membrane operators M S = f ∈S (−h (1 ′ ) f ) with S a contractible surface on the original lattice (by contractible surface we mean an openmembrane whose boundary ∂S is non-zero and is a closed loop) and M S = e∈ S (−h (1 ′ ) ) precisely corresponds to the disorder averaged Wilson loop of the 3D random plaquette gauge model (RPGM) along the Nishimori line[69].It follows that [⟨M S ⟩ 2 ] ∼ e −κ|∂S| ('perimeter-law') when p f < p c ≈ 0.029 while [⟨M S ⟩ 2 ] ∼ e −κ|S| ('area-law') when p f > p c .One can also define the average membrane order parameter [⟨M S ⟩ 2 ] for M S , and the results are analogous with the same critical error rate p c .Therefore, using the qualitative behaviors of [⟨M S ⟩ 2 ] and [⟨M S ⟩ 2 ], one can divide the decohered state ρ as a function of p f and p e into four regimes, see Fig.1(c): (i) p f , p e < p c : both [⟨M S ⟩ 2 ] and [⟨M S ⟩ 2 ] satisfy perimeter law.(ii) p f < p c , p e > p c : [⟨M S ⟩ 2 ] satisfies perimeter-law while [⟨M S ⟩ 2 ] satisfies area-law.(iii) p f > p c , p e < p c : [⟨M S ⟩ 2 ] satisfies area-law while [⟨M S ⟩ 2 ] satisfies perimeter-law.(iv) p f , p e > p c : Both [⟨M S ⟩ 2 ] and [⟨M S ⟩ 2 ] satisfy arealaw.

′ ) 2
symmetry, we choose the basis {|m⟩ = |x e , x f ⟩}.Similar to the previous cases, we consider the 'partition function' Z m = ⟨ψ m |ψ m ⟩ whose singularties are expected to indicate the presence of a phase transition.The evaluation of Z m = ⟨x f , x e |ρ|x f , x e ⟩ is quite similar to that for tr ρ Q (1 ′ ) ,Q (1) and one finds that Z m ∼ Z 3D gauge,x f Z 3D gauge,xe .One may also compute the expectation values of two membrane operators and find ⟨ψ m |M S |ψ m ⟩ = ( f ∈S x f )⟨W ∂S ⟩ 3D gauge,x f and ⟨ψ m |M S |ψ m ⟩ = ( e∈ S x e )⟨W ∂ S ⟩ 3D gauge,xe .Using these one may then define an average membrane order parameters [⟨M S ⟩ 2 ] = m P m ⟨ψ m |M S |ψ m ⟩ 2 and [⟨M S ⟩ 2 ] = m P m ⟨ψ m |M S |ψ m ⟩ 2 .Using same arguments as those following Eq.(

FIG. 2 .
FIG. 2. Modular commutator JABC /J0,ABC on a L × L torus as a function of L corresponding to several different pure states |ψm⟩ that enter the convex decomposition of the p + ip SC subjected to decoherence with Kraus operators linear in Majorana fermions (Eq.(27)), as well as the modular commutator of the decohered mixed state itself.We choose error rate p = 0.04, and the following initial states |m⟩ in Eq.(27): |m⟩ = |0, ..., 0⟩ (uniform), |0, 1, 0, 1, .., 0, 1⟩ (staggered), and |m⟩ = a random bit string in the occupational number basis.The inset shows the geometry of regions A, B, C used to define the modular commutator.We use anti-periodic boundary conditions along both directions so that the ground state is unique.
argue that ρ is sym-SRE for β < β c , we now decompose ρ as ρ = m |ψ m ⟩⟨ψ m | where |ψ m ⟩ = e −βH/2 |x m ⟩/ √ Z.The corresponding Z m = ⟨ψ m |ψ m ⟩ can again be expressed in the continuum limit as an imaginary-time path integral Z m ∼ O k |ψ⟩| 2 /⟨m|ψ⟩.Roughly speaking, the type-II strange correlator is devised to capture the case that ρ can be written as an incoherent sum of pure state |ψ p ⟩, where ⟨m|O j O k |ψ p ⟩ is non-trivial but can be either positive or negative depending on |ψ p ⟩.